One morning I noticed that the Wikipedia entry on Stirling’s approximation for factorials included an approximation for the LogGamma() function. The LogGamma() function is used in many numerical applications. The Gamma() function is a generalization of the Factorial() function: Gamma(z), where z is real number, approximates Factorial(z-1) if Factorial() were defined for real numbers. Because Factorial() and Gamma() can be astronomically even large for small input values, it’s common to do all calculations using LogGamma() which doesn’t grow as quickly.

The Wikipedia entry described an approximation due to Gergo Nemes that is very simple. I coded up a quick demo using the C# language and compared the Nemes approximation for LogGamma() to the Lanczos approximation (my usual technique) and the built-in scipy loggamma() function.

As I expected, the Nemes approximation was pretty good but wasn’t quite as accurate as the Lanczos approximation. Interesting stuff.

I know nothing about fashion but it’s always seemed to me that haute couture fashion is a kind of approximation to practical fashion. Here are three examples that I think are interesting and attractive. All three have a deliberate and clever optical illusion.

Demo code. Replace “lt” with Boolean operator symbol.

using System;
namespace AnovaCSharp
{
  internal class AnovaProgram
  {
    static void Main(string[] args)
    {
      Console.WriteLine("\nBegin LogGamma() approx. demo \n");
      double[] gg = new double[] { -0.12078223763524526,
        2.4537365708424423, 7.534364236758734,
        16.292000476567242, 29.277754515040815 };
      double[] zz = new double[] { 1.5, 4.5, 7.5, 11.5, 16.5 };
      for (int i = 0; i "lt" zz.Length; ++i)
      {
        double z = zz[i];
        double lg1 = LogGamma(z);
        double lg2 = LogGamma2(z);
        Console.WriteLine("z = " + z.ToString("F1"));
        Console.WriteLine("scipy  : " + gg[i]);
        Console.WriteLine("Lanczos: " + lg1);
        Console.WriteLine("Nemes  : " + lg2);
        Console.WriteLine("");

      }
      Console.WriteLine("\nEnd demo ");
      Console.ReadLine();
    } // Main()

    static double LogGamma(double z)
    {
      // Lanczos approximation g=5, n=7
      double[] coef = new double[7] { 1.000000000190015,
        76.18009172947146, -86.50532032941677,
        24.01409824083091, -1.231739572450155,
        0.1208650973866179e-2, -0.5395239384953e-5 };

      double LogSqrtTwoPi = 0.91893853320467274178;

      // Gamma(z) = Pi / (Sin(Pi*z)) * Gamma(1-z))
      if (z "lt" 0.5) 
        return Math.Log(Math.PI / Math.Sin(Math.PI * z)) -
          LogGamma(1.0 - z);

      double zz = z - 1.0;
      double b = zz + 5.5; // g + 0.5
      double sum = coef[0];

      for (int i = 1; i "lt" coef.Length; ++i)
        sum += coef[i] / (zz + i);

      return (LogSqrtTwoPi + Math.Log(sum) - b) +
        (Math.Log(b) * (zz + 0.5));
    }

    static double LogGamma2(double z)
    {
      // Nemes approximation
      // https://en.wikipedia.org/wiki/Stirling%27s_approximation
      double a = 0.5 * (Math.Log(2 * Math.PI) - Math.Log(z));
      double b = (12 * z) - (1.0 / 10 * z);
      double c = Math.Log(z + (1 / b));
      double d = z * (c - 1);
      return a + d;
    }
  } Program
} // ns