Example usage for org.apache.commons.math3.special Gamma Gamma

List of usage examples for org.apache.commons.math3.special Gamma Gamma

Introduction

In this page you can find the example usage for org.apache.commons.math3.special Gamma Gamma.

Prototype

private Gamma() 

Source Link

Document

Default constructor.

Usage

From source file:org.magicdgs.popgenlib.diversity.NucleotideDiversity.java

/**
 * Gets the denominator of Watterson's θ
 * (formula 3 in <a href="http://www.genetics.org/content/123/3/585">Tajima (1989)</a>), which
 * is the {@code numberOfSamples}<SUP>th</SUP> -1 Harmonic Number.
 *
 * This method returns the exact value for up to 49 samples. Otherwise, it uses the
 * EulerMascheroni constant (&gamma;) and an approximation for the digamma distribution
 * (&psi;) for computing the n<SUP>th</SUP> Harmonic Number:
 *
 * <p>&gamma; + &psi;({@code numberOfSamples})
 *
 * @param numberOfSamples number of samples to compute the denominator.
 *
 * @throws IllegalArgumentException if the number of samples is lower than 2.
 * @see Gamma#digamma(double)//from w  w  w . j a va  2 s. c  o m
 */
public static double wattersonsDenominatorApproximation(final int numberOfSamples) {
    // Defined as sum(1/j) for j=1 to j=n-1; where n is the number of samples
    // This is actually the (n-1)th Harmonic number (https://en.wikipedia.org/wiki/Harmonic_number)
    // This could be more efficiently computed using the formula implying the digamma distribution,
    // which is implemented in an approximated way in commons-math3 (good enough for our purposes).
    // In common-math3 they only use the fast algorithm if n is >= 49
    // thus, here the limit for switching to the fast algorithm is 49 too,
    // because it will use the same number of iterations if not
    if (numberOfSamples < 49) {
        return IntStream.range(1, numberOfSamples).mapToDouble(i -> 1d / i).sum();
    }
    // The formula of the nth Harmonic number using the digamma function is defined as:
    // gamma + psi(n + 1); where gamma is the Euler-Mascheroni constant and psi the digamma function
    // because here we want the number numberOfSamples-1 Harmonic number, we can use directly
    // gamma + psi(numberOfSamples) if 49 or more samples were provided (efficient computation)
    return Gamma.GAMMA + Gamma.digamma(numberOfSamples);
}

From source file:reflex.module.ReflexGamma.java

public ReflexValue gamma(List<ReflexValue> params) {
    if (params.size() == 0) {
        return new ReflexValue(Gamma.GAMMA);
    } else {/*  w  w  w . j a v  a  2  s .  c o m*/
        throw new ReflexException(-1, "Where is gamma?");
    }
}