class MoreMath::TDistribution

def inverse_probability(p)

for the probability +p+.
Returns the inverse cumulative probability (t-value) of the TDistribution
def inverse_probability(p)
  case
  when p <= 0
    -1 / 0.0
  when p >= 1
    1 / 0.0
  else 
    begin
      bisect = NewtonBisection.new { |x| probability(x) - p }
      range = bisect.bracket(-10..10)
      bisect.solve(range, 1_000_000)
    rescue
      0 / 0.0
    end
  end
end