Proofs about random programs

Stream: general

[David Ledvinka]

Is there any work on proving facts about the distributions of random programs? Speaking loosely because I am out of my depth here (as a mathematican who just read FPIL) I envision some sort of monadic approach is possible where I can write a program where I can swap the monad between something that "assumes" true random distributions (Giry monad?) that I can prove things about and a monad that lets me actually run the program (IO based), without changing the logic of the program (basically showing that if the input randomness was true random then what I proved about the program is true).

[Shreyas Srinivas]

You can run prove things inside the Id monad

[David Ledvinka]

Can you elaborate? I can see how this allows you to prove things about the deterministic aspects of the program but not how to prove things about the program as a distribution.

[David Ledvinka]

I realized it wasn't totally clear in my initial post but I want to specifically prove things about the distribution of outputs of the program if you assumed the random number generation was true random. An example application would be a program that runs a Monte Carlo or Markov Chain Monte Carlo simulation where you can prove convergence results (in particular that it actually estimates what you want to estimate), confidence intervals etc...

[(deleted)]

1. With a monad it is possible to abstract away the retrieval of random values. Then, we can write a handler to supply the random values to the program to get the output.
2. Now let's assume that our program terminates. Things get more interesting when there is a possibility of non-termination. We now face the first hurdle: we don't know how many random values the program requests. To overcome this hurdle, the handler needs to return an error when not enough random values or too many random values are being supplied, and return the actual program output otherwise.
3. With this, we now have the set of valid lists of random values. Then, by running the program on every list in the set, we also have a multiset of possible outputs. With this multiset, we can start talking about the behavior of programs with random number generation.

[(deleted)]

I'm watching a movie right now. But I can provide code to illustrate what I'm saying.

[Shreyas Srinivas]

To add to the above, you can employ two tricks. Abstract the random variable away and prove things about the resulting function by assuming a distribution on this parameter. The second trick is to use Yao’s minimax principle to prove something about the distribution of inputs instead of distribution of algorithms that arises from randomly chosen parameters.

[Geoffrey Irving]

Here is an example of proof that applies to randomized computations: in the query oracle model, sorting takes at least log2 n! queries in the average case, including for randomized algorithms.

https://github.com/girving/debate/blob/0dab8ad4b58b7a0e82f1738e86135d444462b953/Comp/Sort.lean#L176

[David Ledvinka]

Inspired by this paper and this code I was able to create a rough draft of the approach I was envisioning

import Mathlib

open ENNReal

class RandomBool (m : Type → Type) : Type where
  randomBool : m Bool

instance : RandomBool IO where
  randomBool := IO.rand 0 1

def Distr (α : Type) : Type := α → ℝ≥0∞

-- assume all types are finite so I can sum over `α`
axiom finite_types (α : Type) : Fintype α

noncomputable instance {α : Type} : Fintype α := finite_types _

-- Giry Monad on finite types
noncomputable instance : Monad Distr where
  pure x := fun y => open Classical in if x = y then 1 else 0
  bind {α _} μ f := open Classical in fun y ↦ ∑ x : α, μ x * f x y

noncomputable instance : RandomBool Distr where
  randomBool
  | true => 1 / 2
  | false => 1 / 2

def random_Bool {m} [Monad m] [RandomBool m] : m Bool := do
  return ← RandomBool.randomBool

#eval random_Bool (m := IO) -- returns a random Bool

def random_notBool {m} [Monad m] [RandomBool m] : m Bool := do
  return not (← RandomBool.randomBool)

#eval random_notBool (m := IO) -- returns a random not Bool

/- prove that under the assumption of uniformly distributed Booleans,
random_Bool has the same distribution as random_notBool -/
example : random_Bool (m := Distr) = random_notBool (m := Distr) := by
  simp [Distr, random_Bool, random_notBool, RandomBool.randomBool,
    Bind.bind, Pure.pure]
  ext x
  cases x <;> simp

Obviously it contains some nonsense (I assume all types are fintype so I can take the sum). The idea though is that by replacing the monad m in the same program I can go from a computable program that can actually generate "random" outputs to a noncomputable program that I can reason about (for example proving that random_Bool and random_notBool have the same distribution). One would want to extend this so that you can take a sequence of "independent" random inputs and ideally be able to use Measures in mathlib in place of Distr so that you can use the full strength of mathlib to prove stuff (unfortunately I am not sure if this is easy or doable). Does this approach seem reasonable though?

[(deleted)]

import Mathlib

inductive GetRandom (x : Type)
| Done (v : x)
| Read (continuation : Bool → GetRandom x)

def bind {x} (a : GetRandom x) (f : x → GetRandom x) :=
  match a with
  | GetRandom.Done v => f v
  | GetRandom.Read continuation => GetRandom.Read (fun x => bind (continuation x) f)

-- up to this point it's standard monadic stuff

inductive Outcome {x : Type}
| TooLittle
| TooMuch
| JustRight (a : x)

def supply {x} (l : List Bool) (a : GetRandom x) :=
  match l, a with
  | [], GetRandom.Done k => Outcome.JustRight k
  | [], GetRandom.Read _ => Outcome.TooLittle
  | (_ :: _), GetRandom.Done _ => Outcome.TooMuch
  | (a :: tail), GetRandom.Read k => supply tail (k a)

def possibleOutputs {x} (a : GetRandom x) : Multiset x :=
  match a with
  | GetRandom.Done v => {v}
  | GetRandom.Read continuation =>
    possibleOutputs (continuation true) + possibleOutputs (continuation false)

[(deleted)]

Important caveat: if you have to generate an integer in a range, you have to come up with a different monad. You can't just generate bits and then use rejection sampling. This is because when you use rejection sampling, the program isn't guaranteed to terminate.

[(deleted)]

I feel my approach is fundamentally different from yours. Because my mind is shaped by the industry, and I prefer approaches that just work rather than elegant approaches.

[Geoffrey Irving]

The finite time rejection sampling thing is why for my purposes above the computations are allowed to sample from arbitrary finitely supported distributions. This doesn’t precisely model what real computers do, but is a more convenient setting for a lot of theoretical computer science results.

[Devon Tuma]

You could also use the new free monad stuff in mathlib to write the programs, and embed separately in PMF monad to discuss distributions and IO monad to actually run it. You don't have a direct correspondence between the semantics and the actual run then though. PMF or SPMF is essentially just David's Distr definition

[David Ledvinka]

Oh thanks! I wasn't aware of the PMF monad. This should let me remove the fintype hack I used!

[Devon Tuma]

It uses tsum for infinite cases which can be a bitter tougher to work with. It reduces to finite sums though if you do have the instance around

[David Ledvinka]

Which is potentially a good thing though because it lets you work with infinite (discrete) distributions. I was trying to get this working with Measure but it doesnt work because you need a measurable space instance on your types so you can't make it a monad but I suppose this works because you can give every type the discrete topology and sigma algebra.

[Devon Tuma]

It does also have functions to turn it into measures. But yeah they aren't that useful except for occasionally pulling back lemmas from measure theory, since you need those instances

[(deleted)]

Though I feel using complicated machinery isn't worth it for terminating programs.

[Thomas Zhu]

We have the coin flips measure on (Nat -> Bool) at docs#MeasureTheory.Measure.infinitePiNat and a random program can also be seen as a deterministic program that takes a seed (Nat -> Bool) as additional input. The pushforward of the coin flips measure by this program will be the theoretical behavior of the program and perhaps you may also be able to compose with some IO (Nat -> Bool) to be able to run it?