Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

docs: update documentation for LogProbability functions #817

Merged
merged 3 commits into from
Mar 22, 2022
Merged
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
17 changes: 7 additions & 10 deletions src/main/scala/com/fulcrumgenomics/util/NumericTypes.scala
Original file line number Diff line number Diff line change
Expand Up @@ -146,24 +146,24 @@ object NumericTypes {
}

/**
* Precise computation of log(1-exp(x)) for use in the `aOrNotB` method.
* Precise computation of log(1-exp(-x)) for use in the `aOrNotB` method.
* See Equation (7) in https://cran.r-project.org/web/packages/Rmpfr/vignettes/log1mexp-note.pdf
wmchad marked this conversation as resolved.
Show resolved Hide resolved
*/
private def log1mexp(value: Double): Double = {
if (value < 0) throw new IllegalArgumentException("value was less than zero: " + value)
if (value <= 0) throw new IllegalArgumentException("value was less than or equal to zero: " + value)
else if (value <= LnTwo) log(-expm1(-value))
else log1p(-exp(-value)) // value > LnTwo
}

/** Computes the probability of a or b, where a and b are independent events: Pr(A or B) = Pr(A) + Pr(B). */
/** Computes the probability of a or b, where a and b are mutually exclusive events: Pr(A or B) = Pr(A) + Pr(B). */
def or(a: LogProbability, b: LogProbability): LogProbability = {
if (a.isNegInfinity) b
else if (b.isNegInfinity) a
else if (b < a) or(b, a)
else a + log1pexp(b - a) // for precision we use log1pexp, which is equivalent to log(1+exp(x))
}

/** Computes the probability of any of the given independent events occurring: Pr(AB..N) = Pr(A)+Pr(B)+...+Pr(N). */
/** Computes the probability of any of the given mutually exclusive events occurring: Pr(A,B,..N) = Pr(A)+Pr(B)+...+Pr(N). */
def or(values: Array[LogProbability]): LogProbability = {
if (values.forall(_.isNegInfinity)) Double.NegativeInfinity
else {
Expand All @@ -178,16 +178,13 @@ object NumericTypes {
}
}

/** Computes the probability of a and b, where a and b are independent events: Pr(AB) = Pr(A)*Pr(B). */
/** Computes the probability of a and b, where a and b are independent events: Pr(A,B) = Pr(A)*Pr(B). */
def and(a: LogProbability, b: LogProbability): LogProbability = a + b

/** Computes the probability of the given independent events co-occurring: Pr(AB..N) = Pr(A)*Pr(B)*...*Pr(N). */
/** Computes the probability of the given independent events co-occurring: Pr(A,B,..N) = Pr(A)*Pr(B)*...*Pr(N). */
def and(values: Array[Double]): LogProbability = values.sum

/**
* Computes the probability Pr(A OR not B) = Pr(A) - Pr(B). While this could be computed using
* or(a, not(b)), this form is more efficient and more precise.
*/
/** Computes the probability Pr(A AND not B) = Pr(A) - Pr(B), where B is a subset of A. */
def aOrNotB(a: LogProbability, b: LogProbability): LogProbability = {
if (b.isNegInfinity) a
else if (a == b) Double.NegativeInfinity
Expand Down