Handling Edge Cases / Extreme Cases

Hi all,

I am using rstan and in one of my functions, I have the following code, where B is a 4x4 matrix.
The idea is that I am summing over 16 different combinations of values.

n=B[1,1];
ll = log(sum(exp(B-n))/16)+n;
return(ll);

What would happen if at some point, the entries of the matrix have very different values, i.e., ranging from 100 to 4000? Then, the exponential of a very high number would potentially give infinity, right?

Depending on the initial values, this does not necessarily result in an error of `` Log probability evaluates to log(0), i.e. negative infinity.ā€˜ā€™.

Hence, would I realize that there is a problem? Or is this even a problem?

Iā€™m happy about ideas about how to deal with this.

You can use the log_sum_exp function - see Matrix Operations for the function signature and Floating Point Arithmetic for more explanation. If the entries in B span that many orders of magnitude, then the result will mostly depend on the largest value(s).

2 Likes

If a > b and a - b > 500, then log_sum_exp(a, b) = a in double-precision floating point arithmetic, because it gets computed as follows by pulling out the maximum value and then subtracting before evaluating exp:

log_sum_exp(a, b) 
  = a + log_sum_exp(a - a, b - a)
  = a + log(exp(0) + exp(b - a))
  = a + log(1 + exp(b - a))

The problem is that if b - a < 500, then exp(b - a) underflows to 0, and the result reduces to a.

1 Like