Improve runtime performance via fewer calls

We'll pre-compute and pass around some values which we are highly likely
to need more than once. There's very little overhead in the use once
case. I suppose if it is common to call these and not use the
functions.values ever then this might be more expensive.

I have no idea what kind of side-effects might be expected, though

Some of this code was already done in this style, so I assume this is
readable enough and perhaps welcome.

Author: mwmiller

lib/statistics/distributions/beta.ex

@@ -17,13 +17,15 @@ defmodule Statistics.Distributions.Beta do
   """
   @spec pdf(number, number) :: fun
   def pdf(a, b) do
+    bab = Functions.beta(a, b)
+
     fn x ->
       cond do
         x <= 0.0 ->
           0.0
 
         true ->
-          Math.pow(x, a - 1) * Math.pow(1 - x, b - 1) / Functions.beta(a, b)
+          Math.pow(x, a - 1) * Math.pow(1 - x, b - 1) / bab
       end
     end
   end
@@ -91,15 +93,17 @@ defmodule Statistics.Distributions.Beta do
 
   """
   @spec rand(number, number) :: number
-  def rand(a, b) do
+  def rand(a, b), do: rand(pdf(a, b))
+
+  defp rand(rpdf) do
     # beta only exists between 0 and 1
     x = Math.rand()
 
-    if pdf(a, b).(x) > Math.rand() do
+    if rpdf.(x) > Math.rand() do
       x
     else
       # keep trying
-      rand(a, b)
+      rand(rpdf)
     end
   end
 end

lib/statistics/distributions/binomial.ex

@@ -66,18 +66,18 @@ defmodule Statistics.Distributions.Binomial do
   @spec ppf(non_neg_integer, number) :: fun
   def ppf(n, p) do
     fn x ->
-      ppf_tande(x, n, p)
+      ppf_tande(x, n, p, cdf(n, p), 0)
     end
   end
 
   # trial-and-error method which refines guesses
   # to arbitrary number of decimal places
-  defp ppf_tande(x, n, p, g \\ 0) do
-    g_cdf = cdf(n, p).(g)
+  defp ppf_tande(x, n, p, npcdf, g) do
+    g_cdf = npcdf.(g)
 
     cond do
       x > g_cdf ->
-        ppf_tande(x, n, p, g + 1)
+        ppf_tande(x, n, p, npcdf, g + 1)
 
       x <= g_cdf ->
         g
@@ -96,14 +96,16 @@ defmodule Statistics.Distributions.Binomial do
 
   """
   @spec rand(non_neg_integer, number) :: non_neg_integer
-  def rand(n, p) do
+  def rand(n, p), do: rand(n, p, pmf(n, p))
+
+  defp rand(n, p, rpmf) do
     x = Math.rand() * n
 
-    if pmf(n, p).(x) > Math.rand() do
+    if rpmf.(x) > Math.rand() do
       Float.round(x)
     else
       # keep trying
-      rand(n, p)
+      rand(n, p, rpmf)
     end
   end
 end

lib/statistics/distributions/chisq.ex

@@ -19,10 +19,10 @@ defmodule Statistics.Distributions.Chisq do
   """
   @spec pdf(non_neg_integer) :: fun
   def pdf(df) do
-    fn x ->
-      1 / (Math.pow(2, df / 2) * Functions.gamma(df / 2)) * Math.pow(x, df / 2 - 1) *
-        Math.exp(-1 * x / 2)
-    end
+    hdf = df / 2
+    g = Math.pow(2, hdf) * Functions.gamma(hdf)
+
+    fn x -> 1 / g * Math.pow(x, hdf - 1) * Math.exp(-1 * x / 2) end
   end
 
   @doc """
@@ -36,9 +36,11 @@ defmodule Statistics.Distributions.Chisq do
   """
   @spec cdf(non_neg_integer) :: fun
   def cdf(df) do
+    hdf = df / 2.0
+    g = Functions.gamma(hdf)
+
     fn x ->
-      g = Functions.gamma(df / 2.0)
-      b = Functions.gammainc(df / 2.0, x / 2.0)
+      b = Functions.gammainc(hdf, x / 2.0)
       b / g
     end
   end
@@ -55,28 +57,28 @@ defmodule Statistics.Distributions.Chisq do
   @spec ppf(non_neg_integer) :: fun
   def ppf(df) do
     fn x ->
-      ppf_tande(x, df)
+      ppf_tande(x, cdf(df))
     end
   end
 
   # trial-and-error method which refines guesses
   # to arbitrary number of decimal places
-  defp ppf_tande(x, df, precision \\ 14) do
-    ppf_tande(x, df, 0, precision + 2, 0)
+  defp ppf_tande(x, tcdf, precision \\ 14) do
+    ppf_tande(x, tcdf, 0, precision + 2, 0)
   end
 
   defp ppf_tande(_, _, g, precision, precision) do
     g
   end
 
-  defp ppf_tande(x, df, g, precision, p) do
+  defp ppf_tande(x, tcdf, g, precision, p) do
     increment = 100 / Math.pow(10, p)
     guess = g + increment
 
-    if x < cdf(df).(guess) do
-      ppf_tande(x, df, g, precision, p + 1)
+    if x < tcdf.(guess) do
+      ppf_tande(x, tcdf, g, precision, p + 1)
     else
-      ppf_tande(x, df, guess, precision, p)
+      ppf_tande(x, tcdf, guess, precision, p)
     end
   end
 
@@ -92,14 +94,16 @@ defmodule Statistics.Distributions.Chisq do
 
   """
   @spec rand(non_neg_integer) :: number
-  def rand(df) do
+  def rand(df), do: rand(df, cdf(df))
+
+  defp rand(df, rcdf) do
     x = Math.rand() * 100
 
-    if pdf(df).(x) > Math.rand() do
+    if rcdf.(x) > Math.rand() do
       x
     else
       # keep trying
-      rand(df)
+      rand(df, rcdf)
     end
   end
 end

lib/statistics/distributions/f.ex

@@ -20,11 +20,14 @@ defmodule Statistics.Distributions.F do
   """
   @spec pdf(number, number) :: fun
   def pdf(d1, d2) do
+    powa = Math.pow(d2, d2)
+    cfac = Functions.beta(d1 / 2, d2 / 2)
+
     fn x ->
       # create components
-      a = Math.pow(d1 * x, d1) * Math.pow(d2, d2)
+      a = Math.pow(d1 * x, d1) * powa
       b = Math.pow(d1 * x + d2, d1 + d2)
-      c = x * Functions.beta(d1 / 2, d2 / 2)
+      c = x * cfac
       # for the equation
       Math.sqrt(a / b) / c
     end
@@ -47,9 +50,11 @@ defmodule Statistics.Distributions.F do
   # which is the CDF of the Beta distribution
   @spec cdf(number, number) :: fun
   def cdf(d1, d2) do
+    bcdf = Beta.cdf(d1 / 2, d2 / 2)
+
     fn x ->
       xx = d1 * x / (d1 * x + d2)
-      Beta.cdf(d1 / 2, d2 / 2).(xx)
+      bcdf.(xx)
     end
   end
 

lib/statistics/distributions/hypergeometric.ex

@@ -17,6 +17,8 @@ defmodule Statistics.Distributions.Hypergeometric do
   """
   @spec pmf(non_neg_integer, non_neg_integer, non_neg_integer) :: fun
   def pmf(pn, pk, n) do
+    combos = Math.combination(pn, n)
+
     fn k ->
       cond do
         n < k ->
@@ -33,7 +35,7 @@ defmodule Statistics.Distributions.Hypergeometric do
 
         true ->
           xk = Math.to_int(k)
-          Math.combination(pk, xk) * Math.combination(pn - pk, n - xk) / Math.combination(pn, n)
+          Math.combination(pk, xk) * Math.combination(pn - pk, n - xk) / combos
       end
     end
   end
@@ -46,10 +48,12 @@ defmodule Statistics.Distributions.Hypergeometric do
   """
   @spec cdf(non_neg_integer, non_neg_integer, non_neg_integer) :: fun
   def cdf(pn, pk, n) do
+    cpmf = pmf(pn, pk, n)
+
     fn k ->
       0..Math.to_int(Math.floor(k))
       |> Enum.to_list()
-      |> Enum.map(fn i -> pmf(pn, pk, n).(i) end)
+      |> Enum.map(fn i -> cpmf.(i) end)
       |> Enum.sum()
     end
   end
@@ -63,19 +67,19 @@ defmodule Statistics.Distributions.Hypergeometric do
   @spec ppf(non_neg_integer, non_neg_integer, non_neg_integer) :: fun
   def ppf(pn, pk, n) do
     fn x ->
-      ppf_tande(x, pn, pk, n)
+      ppf_tande(x, cdf(pn, pk, n), 0)
     end
   end
 
   # trial-and-error method which refines guesses
   # to arbitrary number of decimal places
 
-  defp ppf_tande(x, pn, pk, n, guess \\ 0) do
-    g_cdf = cdf(pn, pk, n).(guess)
+  defp ppf_tande(x, tcdf, guess) do
+    g_cdf = tcdf.(guess)
 
     cond do
       x > g_cdf ->
-        ppf_tande(x, pn, pk, n, guess + 1)
+        ppf_tande(x, tcdf, guess + 1)
 
       x <= g_cdf ->
         guess
@@ -86,14 +90,16 @@ defmodule Statistics.Distributions.Hypergeometric do
   Draw a random number from hypergeometric distribution
   """
   @spec rand(non_neg_integer, non_neg_integer, non_neg_integer) :: non_neg_integer
-  def rand(pn, pk, n) do
+  def rand(pn, pk, n), do: rand(pk, pmf(pn, pk, n))
+
+  defp rand(pk, rpmf) do
     x = Math.floor(Math.rand() * pk)
 
-    if pmf(pn, pk, n).(x) > Math.rand() do
+    if rpmf.(x) > Math.rand() do
       Float.round(x)
     else
       # keep trying
-      rand(pn, pk, n)
+      rand(pk, rpmf)
     end
   end
 end

lib/statistics/distributions/normal.ex

@@ -61,8 +61,10 @@ defmodule Statistics.Distributions.Normal do
 
   @spec cdf(number, number) :: fun
   def cdf(mu, sigma) do
+    denom = sigma * Math.sqrt(2)
+
     fn x ->
-      0.5 * (1.0 + Functions.erf((x - mu) / (sigma * Math.sqrt(2))))
+      0.5 * (1.0 + Functions.erf((x - mu) / denom))
     end
   end
 
@@ -127,27 +129,25 @@ defmodule Statistics.Distributions.Normal do
   end
 
   @spec rand(number, number) :: number
-  def rand(mu, sigma) do
+  def rand(mu, sigma), do: rand(mu, sigma, pdf(0, 1))
+
+  defp rand(mu, sigma, rpdf) do
     # Note: an alternate method exists and may be better
     # Inverse transform sampling - https://en.wikipedia.org/wiki/Inverse_transform_sampling
     # ----
     # Generate a random number between -10,+10
     # (probability of 10 ocurring in a Normal(0,1) distribution is
     # too small to calculate with the precision available to us)
     x = Math.rand() * 20 - 10
-    rmu = 0
-    rsigma = 1
 
     cond do
-      pdf(rmu, rsigma).(x) > Math.rand() ->
-        # get z-score
-        z = (rmu - x) / rsigma
+      rpdf.(x) > Math.rand() ->
         # transpose to specified distribution
-        mu + z * sigma
+        mu - x * sigma
 
       true ->
         # keep trying
-        rand(mu, sigma)
+        rand(mu, sigma, rpdf)
     end
   end
 end

lib/statistics/distributions/poisson.ex

@@ -21,8 +21,10 @@ defmodule Statistics.Distributions.Poisson do
   """
   @spec pmf(number) :: fun
   def pmf(lambda) do
+    nexp = Math.exp(-lambda)
+
     fn k ->
-      Math.pow(lambda, k) / Math.factorial(k) * Math.exp(-lambda)
+      Math.pow(lambda, k) / Math.factorial(k) * nexp
     end
   end
 
@@ -37,12 +39,14 @@ defmodule Statistics.Distributions.Poisson do
   """
   @spec cdf(number) :: fun
   def cdf(lambda) do
+    nexp = Math.exp(-1 * lambda)
+
     fn k ->
       s =
         Enum.map(0..Math.to_int(k), fn x -> Math.pow(lambda, x) / Math.factorial(x) end)
         |> Enum.sum()
 
-      Math.exp(-lambda) * s
+      nexp * s
     end
   end
 
@@ -61,15 +65,17 @@ defmodule Statistics.Distributions.Poisson do
   """
   @spec ppf(number) :: fun
   def ppf(lambda) do
+    lcdf = cdf(lambda)
+
     fn x ->
-      ppf_tande(x, lambda, 0.0)
+      ppf_tande(x, lcdf, 0.0)
     end
   end
 
   # the trusty trial-and-error method
-  defp ppf_tande(x, lambda, guess) do
-    if x > cdf(lambda).(guess) do
-      ppf_tande(x, lambda, guess + 1)
+  defp ppf_tande(x, lcdf, guess) do
+    if x > lcdf.(guess) do
+      ppf_tande(x, lcdf, guess + 1)
     else
       guess
     end
@@ -87,14 +93,16 @@ defmodule Statistics.Distributions.Poisson do
 
   """
   @spec rand(number) :: number
-  def rand(lambda) do
+  def rand(lambda), do: rand(lambda, pmf(lambda))
+
+  defp rand(lambda, lpmf) do
     x = (Math.rand() * 100 + lambda) |> Math.floor()
 
-    if pmf(lambda).(x) > Math.rand() do
+    if lpmf.(x) > Math.rand() do
       x
     else
       # keep trying
-      rand(lambda)
+      rand(lambda, lpmf)
     end
   end
 end