Scala versions: 2.13 2.12 2.11 2.10

Makes it easy to create, manipulate and sample probability distributions.

## Installation

Include this in your sbt config:

``````"org.jliszka" %% "probability-monad" % "1.0.3"
``````

## Examples

Here's how you would code up the following problem: You are given either a fair coin or a biased coin with equal probability. If you flip it 5 times and it comes up heads each time, what is the probability you have the fair coin?

```case class Trial(haveFairCoin: Boolean, flips: List[Coin])

def bayesianCoin(nflips: Int): Distribution[Trial] = {
for {
haveFairCoin <- tf()
c = if (haveFairCoin) coin else biasedCoin(0.9)
flips <- c.repeat(nflips)
} yield Trial(haveFairCoin, flips)
}

bayesianCoin(5).given(_.flips.forall(_ == H)).pr(_.haveFairCoin)```

Or: You repeatedly roll a 6-sided die and keep a running sum. What is the probability the sum reaches exactly 30?

```def dieSum(rolls: Int): Distribution[List[Int]] = {
markov(rolls, List(0))(runningSum => for {
d <- die
} yield (d + runningSum.head) :: runningSum)
}

dieSum(30).pr(_ contains 30)```

Or: Each family has children until it has a boy, and then stops. What is the expected fraction of girls in the population?

```sealed abstract class Child
case object Boy extends Child
case object Girl extends Child

def family = {
discreteUniform(List(Boy, Girl)).until(_ contains Boy)
}

def population(families: Int) = {
for {
children <- family.repeat(families).map(_.flatten)
val girls = children.count(_ == Girl)
} yield 1.0 * girls / children.length
}

population(4).ev```

## How it works

A `Distribution[T]` represents a random variable that, when sampled, produces values of type `T` according to a particular probability distribution. For example, `Distribution.uniform` is a `Distribution[Double]` that produces `Double` values between 0.0 and 1.0, uniformly distributed. `Distribution.coin` is a `Distribution[Coin]` that produces the values `H` and `T` with equal probability, and `Distribution.biasedCoin(0.3)` is a `Distribution[Coin]` that produces the value `H` 30% of the time and the value `T` 70% of the time.

You can think of a `Distribution[T]` as a collection like any other scala collection that you can `map`, `flatMap` and `filter` over. The presence of these methods allow you to use scala's for-comprehensions to manipulate distributions. For example, here's how you would create a distribution that represents the sum of 2 die rolls:

```val dice = for {
d1 <- die
d2 <- die
} yield d1 + d2```

Here, `die` is a `Distribution[Int]`, and `d1` and `d2` are both `Int`s. The type of `dice` is `Distribution[Int]`. You can see that for-comprehensions are an easy way to define new a distribution in terms of individual samples from other distributions.

You can visualize a distribution with `hist`:

``````scala> dice.hist
2  2.61% ##
3  5.48% #####
4  8.70% ########
5 10.53% ##########
6 14.21% ##############
7 16.90% ################
8 13.90% #############
9 11.43% ###########
10  8.35% ########
11  5.17% #####
12  2.72% ##
``````

If you want more control over the display of continuous distributions, use `bucketedHist`:

``````scala> normal.map(_ * 2 + 1).bucketedHist(20)   // 20 buckets, min & max determined automatically
-7.0  0.02%
-6.0  0.03%
-5.0  0.22%
-4.0  1.01% #
-3.0  2.69% ##
-2.0  6.43% ######
-1.0 12.19% ############
0.0 17.07% #################
1.0 19.74% ###################
2.0 17.55% #################
3.0 12.17% ############
4.0  6.55% ######
5.0  2.92% ##
6.0  1.10% #
7.0  0.23%
8.0  0.06%
9.0  0.01%
10.0  0.01%

scala> cauchy.bucketedHist(-10, 10, 20)   // min=-10, max=10, #buckets=20
-10.0  0.20%
-9.0  0.38%
-8.0  0.44%
-7.0  0.55%
-6.0  0.82%
-5.0  1.23% #
-4.0  1.85% #
-3.0  2.92% ##
-2.0  6.78% ######
-1.0 16.78% ################
0.0 30.04% ##############################
1.0 16.64% ################
2.0  6.22% ######
3.0  3.06% ###
4.0  1.76% #
5.0  1.26% #
6.0  0.84%
7.0  0.67%
8.0  0.48%
9.0  0.42%
10.0  0.14%
``````

This probability monad is based on sampling, so the values and plots produced will be inexact and will vary between runs.

``````scala> normal.stdev
res9: Double = 1.0044818262040809

scala> normal.stdev
res10: Double = 1.0071194147525722
``````

# Code and examples

Distribution.scala contains code for creating and manipulating probability distributions. Built-in distributions include:

• uniform discrete (including die and fair coin)
• weighted discrete (biased coin, uses the alias method)
• bernoulli
• geometric
• binomial
• negative binomial
• poisson
• zipf
• uniform continuous
• normal
• cauchy
• chi2
• pareto
• exponential
• lognormal
• student's t-distribution
• gamma
• beta

Methods for manipulating distributions include:

• adding (convolution), subracting (cross-correlation), multiplying and dividing distributions
• joint distributions (a flatMap)
• marginal distributions (a filter)
• creating Markov chains (an iterated flatMap)
• finding the probability of arbitrary predicates, conditional probabililty
• finding expected values (mean), standard deviation, variance, skewness and kurtosis
• sampling, histogram

Examples.scala contains some example uses, and possibly a RISK simulator.

To try out some examples, do

``````\$ ./sbt console

scala> runBayesianCoin(5)
``````

Contributions welcome!