Skip to content

Commit

Permalink
address review comments
Browse files Browse the repository at this point in the history
  • Loading branch information
akiezun committed Jun 10, 2015
1 parent 466faf9 commit afa5001
Show file tree
Hide file tree
Showing 22 changed files with 862 additions and 581 deletions.
71 changes: 45 additions & 26 deletions src/main/java/org/broadinstitute/hellbender/utils/MathUtils.java
Original file line number Diff line number Diff line change
Expand Up @@ -80,18 +80,20 @@ private static final class JacobianLogTable {
public static final double MAX_TOLERANCE = 8.0;

public static double get(final double difference) {
if (cache == null)
if (cache == null) {
initialize();
}
final int index = fastRound(difference * INV_STEP);
return cache[index];
}

private static synchronized void initialize() {
private static void initialize() {
if (cache == null) {
final int tableSize = (int) (MAX_TOLERANCE / TABLE_STEP) + 1;
cache = new double[tableSize];
for (int k = 0; k < cache.length; k++)
for (int k = 0; k < cache.length; k++) {
cache[k] = Math.log10(1.0 + Math.pow(10.0, -((double) k) * TABLE_STEP));
}
}
}

Expand All @@ -117,8 +119,9 @@ public static double approximateLog10SumLog10(final double[] vals, final int end
double approxSum = vals[maxElementIndex];

for (int i = 0; i < endIndex; i++) {
if (i == maxElementIndex || vals[i] == Double.NEGATIVE_INFINITY)
if (i == maxElementIndex || vals[i] == Double.NEGATIVE_INFINITY) {
continue;
}

final double diff = approxSum - vals[i];
if (diff < JacobianLogTable.MAX_TOLERANCE) {
Expand All @@ -134,20 +137,20 @@ public static double approximateLog10SumLog10(final double a, final double b, fi
return approximateLog10SumLog10(a, approximateLog10SumLog10(b, c));
}

public static double approximateLog10SumLog10(double small, double big) {
public static double approximateLog10SumLog10(final double small, final double big) {

This comment has been minimized.

Copy link
@vruano

vruano Jun 15, 2015

Contributor

Please add javadoc indicating that small does not necessarily need to be equal or smaller than big, or change names to something smarter.

This comment has been minimized.

Copy link
@akiezun

akiezun Jun 17, 2015

Author Contributor

ok

// make sure small is really the smaller value
if (small > big) {
final double t = big;
big = small;
small = t;
return approximateLog10SumLog10(big, small);
}

if (small == Double.NEGATIVE_INFINITY || big == Double.NEGATIVE_INFINITY)
if (small == Double.NEGATIVE_INFINITY || big == Double.NEGATIVE_INFINITY) {

This comment has been minimized.

Copy link
@vruano

vruano Jun 15, 2015

Contributor

Suspect that || big == Double.NEGATIVE_INFINITY) could simply be removed here ... if you are adding up 0 + X the result is X. So return big would always work.

This comment has been minimized.

Copy link
@akiezun

akiezun Jun 17, 2015

Author Contributor

done. added tests too.

return big;
}

final double diff = big - small;
if (diff >= JacobianLogTable.MAX_TOLERANCE)
if (diff >= JacobianLogTable.MAX_TOLERANCE) {
return big;
}

// OK, so |y-x| < tol: we use the following identity then:
// we need to compute log10(10^x + 10^y)
Expand Down Expand Up @@ -233,6 +236,22 @@ public static double log10BinomialCoefficient(final int n, final int k) {
return log10Factorial(n) - log10Factorial(k) - log10Factorial(n - k);
}

/**
* Computes a binomial probability. This is computed using the formula
* <p/>
* B(k; n; p) = [ n! / ( k! (n - k)! ) ] (p^k)( (1-p)^k )
* <p/>
* where n is the number of trials, k is the number of successes, and p is the probability of success
*
* @param n number of Bernoulli trials
* @param k number of successes
* @param p probability of success
* @return the binomial probability of the specified configuration. Computes values down to about 1e-237.
*/
public static double binomialProbability(final int n, final int k, final double p) {
return Math.pow(10, log10BinomialProbability(n, k, Math.log10(p)));
}

/**
* binomial Probability(int, int, double) with log10 applied to result
*/
Expand All @@ -253,21 +272,21 @@ public static double log10sumLog10(final double[] log10values) {


public static double log10sumLog10(final double[] log10p, final int start, final int finish) {

if (start >= finish)
if (start >= finish) {
return Double.NEGATIVE_INFINITY;
final int maxElementIndex = MathUtils.maxElementIndex(log10p, start, finish);
}
final int maxElementIndex = maxElementIndex(log10p, start, finish);
final double maxValue = log10p[maxElementIndex];
if(maxValue == Double.NEGATIVE_INFINITY)
if(maxValue == Double.NEGATIVE_INFINITY) {
return maxValue;
}
double sum = 1.0;
for (int i = start; i < finish; i++) {
double curVal = log10p[i];
double scaled_val = curVal - maxValue;
if (i == maxElementIndex || curVal == Double.NEGATIVE_INFINITY) {
continue;
}
else {
} else {
sum += Math.pow(10.0, scaled_val);
}
}
Expand All @@ -277,6 +296,16 @@ public static double log10sumLog10(final double[] log10p, final int start, final
return maxValue + (sum != 1.0 ? Math.log10(sum) : 0.0);
}

/**
* normalizes the log10-based array. ASSUMES THAT ALL ARRAY ENTRIES ARE <= 0 (<= 1 IN REAL-SPACE).
*
* @param array the array to be normalized
* @return a newly allocated array corresponding the normalized values in array
*/
public static double[] normalizeFromLog10(final double[] array) {
return normalizeFromLog10(array, false);
}

/**
* normalizes the log10-based array. ASSUMES THAT ALL ARRAY ENTRIES ARE <= 0 (<= 1 IN REAL-SPACE).
*
Expand Down Expand Up @@ -334,16 +363,6 @@ public static double[] normalizeFromLog10(final double[] array, final boolean ta
return normalized;
}

/**
* normalizes the log10-based array. ASSUMES THAT ALL ARRAY ENTRIES ARE <= 0 (<= 1 IN REAL-SPACE).
*
* @param array the array to be normalized
* @return a newly allocated array corresponding the normalized values in array
*/
public static double[] normalizeFromLog10(final double[] array) {
return normalizeFromLog10(array, false);
}

/**
* normalizes the real-space probability array.
*
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -13,7 +13,7 @@
*/
public final class RandomDNA {

private Random random;
private final Random random;

/**
* Constructs a new random DNA generator.
Expand All @@ -23,7 +23,6 @@ public final class RandomDNA {
* described in {@link Random} documentation.
* </p>
*/
@SuppressWarnings("unused")
public RandomDNA() {
random = new Random();

This comment has been minimized.

Copy link
@vruano

vruano Jun 15, 2015

Contributor

Changed to this(new Random())

This comment has been minimized.

Copy link
@akiezun

akiezun Jun 17, 2015

Author Contributor

done

}
Expand All @@ -36,8 +35,9 @@ public RandomDNA() {
* @throws IllegalArgumentException if {@code rnd} is {@code null}.
*/
public RandomDNA(final Random rnd) {
if (rnd == null)
if (rnd == null) {

This comment has been minimized.

Copy link
@vruano

vruano Jun 15, 2015

Contributor

Write instead random = Utils.nonNull(rnd,"the random number generator cannot be null")

This comment has been minimized.

Copy link
@akiezun

akiezun Jun 17, 2015

Author Contributor

done

throw new IllegalArgumentException("the random number generator cannot be null");
}
random = rnd;
}

Expand Down
Loading

0 comments on commit afa5001

Please sign in to comment.