List of usage examples for org.apache.commons.math3.special Gamma Gamma
private Gamma()
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 (γ) and an approximation for the digamma distribution * (ψ) for computing the n<SUP>th</SUP> Harmonic Number: * * <p>γ + ψ({@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?"); } }