0

I am trying to calculate the value of the Gaussian function using the following definition

from math import e

def function(x):
    return 100 / (e**((x-2000)**2))**1/16

But I am getting an overflow with values of 20000, 30000, 1000000, etc. Do you know a workaround Thanks!

meowgoesthedog
  • 14,670
  • 4
  • 27
  • 40
polmonroig
  • 937
  • 4
  • 12
  • 22
  • 2
    Possible duplicate of [OverflowError: (34, 'Result too large')](https://stackoverflow.com/questions/20201706/overflowerror-34-result-too-large) – Dinko Pehar Jan 22 '19 at 12:56
  • `e**(500)` is `1.4035922178528e+217`. How could you possibly need `e**(18000)`, let alone squaring the exponent and using even higher exponents? – roganjosh Jan 22 '19 at 12:57
  • Just think about it - `x = 20000 => (x-2000)^2 = 324000000` and then you calculate `e^324000000`. No wonder you get an overflow... – Tobias Brösamle Jan 22 '19 at 12:57
  • https://stackoverflow.com/questions/20201706/overflowerror-34-result-too-large – naivepredictor Jan 22 '19 at 12:59

2 Answers2

2

  1. Use math.exp:

    100 / math.exp((x-2000)**2))**1/16

    This is more performant and numerically accurate than e**.

  2. Note that ...**1/16 probably doesn't do what you want. ** has higher operator precedence than /, so the above is equivalent to (...**1)/16. Wrap the fraction in parentheses:

    100 / math.exp((x-2000)**2)**(1/16.0)

    Writing 16 in floating point format (16.0) to prevent integer division.

  3. Use exponent rules to incorporate the 1/16:

    100 / math.exp((x-2000)**2/16.0)

    This reduces large arguments, which helps to improve precision.

  4. Flip the sign of the argument and change the division to a multiplication:

    100 * math.exp(-(x-2000)**2/16.0)

    This will underflow instead of overflow in extreme cases. Underflows are much more well-behaved in this situation – you'll just get zero instead of an OverflowError.

meowgoesthedog
  • 14,670
  • 4
  • 27
  • 40
  • I'd recommend changing 2/16 to 0.125 if the complier isn't doing it for you. – duffymo Jan 22 '19 at 13:14
  • @duffymo `**` has [higher precedence](https://docs.python.org/3/reference/expressions.html#operator-precedence) than `/`; also would it really help (performance or accuracy -wise) to do this? – meowgoesthedog Jan 22 '19 at 13:17
  • 1
    No one will notice any performance issues. I worry about integer versus floating point division. When I enter 2/16 into my Python REPL it returns zero instead of 1/8 = 0.125. Looks like I have Python 2.7.15 running. – duffymo Jan 22 '19 at 13:18
  • @duffymo you are right, thanks for the catch; I assumed that `**` would perform a conversion to floating point. – meowgoesthedog Jan 22 '19 at 13:20
  • I'm glad to help. Edit your answer. – duffymo Jan 22 '19 at 13:21
  • @duffymo your edit produces different results because `**` has higher precedence than `/` – i.e. `(x - 2000)**2/16 == ((x - 2000)**2)/16 != (x-2000)**(2/16)`. Also changing the denominator to `16.0` works too (and the value of the standard deviation is clearer when written in this way). – meowgoesthedog Jan 22 '19 at 13:39
  • Minor nitpicking: instead of "performant", say "faster" or "uses less memory" or whatever specific behavior you have in mind. – Robert Dodier Jan 22 '19 at 18:19
1

My advice is to rework the formula by taking the logarithm of it (i.e. convert the formula with the product of terms into the sum of the logarithm of each term) and then take the exponential of the result.

Robert Dodier
  • 16,905
  • 2
  • 31
  • 48