(Attempt at) the Kumaraswamy distribution in Rainier.

Kumaraswamy.scala

// Doesn't work yet!
object Kumaraswamy {
  def apply(a: Real, b: Real): Continuous = new Continuous {

    def realLogDensity(x: Real): Real =
      If(x <= 0,
         Real.zero.log,
         If(x >= 1, Real.zero.log, (a - 1) * x.log + (b - 1) * (1 - x.pow(a))))

    def param: RandomVariable[Real] = {
      for {
        u <- Uniform(0, 1).param
      } yield (1 - (1 - u).pow(1.0 / a)).pow(1.0 / b)
    }

    val generator: Generator[Double] = Generator.from { (r, n) =>
      val u = r.standardUniform
      math.pow(1 - math.pow(1 - u, 1 / n.toDouble(b)), 1 / n.toDouble(a))
    }
  }
}