¿Cuándo deben usarse log1p y expm1?

30

Tengo una pregunta simple que es realmente difícil para Google (además del canónico Lo que todo informático debe saber sobre el papel de aritmética de punto flotante ).

¿Cuándo deberían usarse funciones como log1po en expm1lugar de logy exp? ¿Cuándo no deberían ser utilizados? ¿Cómo difieren las diferentes implementaciones de esas funciones en términos de su uso?

Tim
fuente
2
Bienvenido a Scicomp.SE! Esa es una pregunta muy razonable, pero sería más fácil de responder si explicara un poco a qué log1p se refiere (especialmente cómo se implementa, por lo que no tenemos que adivinar).
Christian Clason
44
Para argumentos de valor real, log1p y expm1 deben usarse cuando es pequeño, por ejemplo, cuando en precisión de coma flotante. Ver, por ejemplo, docs.scipy.org/doc/numpy/reference/generated/numpy.expm1.html y docs.scipy.org/doc/numpy/reference/generated/numpy.log1p.html . (x)(x)x1+x=1
GoHokies
@ChristianClason gracias, me estoy refiriendo principalmente a C ++ std o R, pero cuando me preguntan, empiezo a pensar que aprender sobre las diferencias en las implementaciones también sería muy interesante.
Tim
Aquí hay un ejemplo: scicomp.stackexchange.com/questions/8371/…
Juan M. Bello-Rivas
1
@ user2186862 "cuando es pequeño" es correcto, pero no solo "cuando en precisión de coma flotante" (que ocurre para en la aritmética de doble precisión habitual). Las páginas de documentación que vinculó muestran que ya son útiles para , por ejemplo. x1+x=1x1016x1010
Federico Poloni

Respuestas:

25

Todos sabemos que implica que para , tenemos . Esto significa que si tenemos que evaluar en coma flotante , para puede producir una cancelación catastrófica.

exp(x)=n=0xnn!=1+x+12x2+
|x|1exp(x)1+xexp(x)1|x|1

Esto se puede demostrar fácilmente en python:

>>> from math import (exp, expm1)

>>> x = 1e-8
>>> exp(x) - 1
9.99999993922529e-09
>>> expm1(x)
1.0000000050000001e-08

>>> x = 1e-22
>>> exp(x) - 1
0.0
>>> expm1(x)
1e-22

Los valores exactos son

exp(108)1=0.000000010000000050000000166666667083333334166666668exp(1022)1=0.000000000000000000000100000000000000000000005000000

En general, una implementación "precisa" de expy expm1debe ser correcta a no más de 1ULP (es decir, una unidad del último lugar). Sin embargo, dado que lograr esta precisión da como resultado un código "lento", a veces hay disponible una implementación rápida y menos precisa. Por ejemplo en CUDA tenemos expfy expm1f, donde fsignifica rápido. De acuerdo con la guía de programación CUDA C, aplicación. D el expftiene un error de 2ULP.

Si no le importan los errores en el orden de pocas ULPS, generalmente las implementaciones diferentes de la función exponencial son equivalentes, pero tenga en cuenta que los errores pueden estar ocultos en alguna parte ... (¿Recuerda el error Pentium FDIV ?)

Por lo tanto, está bastante claro que expm1debería usarse para calcular para pequeña . Usarlo para general no es dañino, ya que se espera que sea preciso en todo su rango:exp(x)1xxexpm1

>>> exp(200)-1 == exp(200) == expm1(200)
True

(En el ejemplo anterior está muy por debajo de 1ULP de , por lo que las tres expresiones devuelven exactamente el mismo número de coma flotante).1exp(200)

Una discusión similar es válida para las funciones inversas logy log1pdesde para .log(1+x)x|x|1

Stefano M
fuente
1
Esta respuesta ya estaba contenida en los comentarios a la pregunta de OP. Sin embargo, me pareció útil dar una cuenta más larga (aunque básica) solo por claridad, con la esperanza de que sea útil para algunos lectores.
Stefano M
OK, pero entonces uno puede simplemente concluir "para que siempre pueda usar expm1 en lugar de exp" ...
Tim
1
@tim su conclusión es incorrecta: siempre puede usar en expm1(x)lugar de exp(x)-1. Por supuesto que exp(x) == exp(x) - 1no se cumple en general.
Stefano M
OK, eso está claro. ¿Y hay algún criterio claro para ? x1
Tim
1
@Tim no hay un umbral claro, y la respuesta depende de la precisión de la implementación de coma flotante (y el problema que se está resolviendo). Si bien expm1(x)debe ser preciso a 1ULP en todo el rango , pierde progresivamente la precisión de algunos ULP cuando a un desglose completo cuando , donde es máquina-epsilon. 0x1exp(x) - 1x1x<ϵϵ
Stefano M
1

Para ampliar la diferencia entre logy log1ppodría ayudar recordar el gráfico si el logaritmo:

Logaritmo

Si sus datos contienen ceros, entonces probablemente no quiera usarlos, logya que no están definidos en cero. Y a medida que acerca a , el valor de acerca a . Entonces, si sus valores de están cerca de , entonces el valor de es potencialmente un número negativo grande. Por ejemplo y y así sucesivamente. Esto puede ser útil, pero también puede distorsionar sus datos hacia números negativos grandes, especialmente si su conjunto de datos también contiene números mucho más grandes que cero.x0ln(x)x0ln(x)ln(1e)=1ln(1e10)=10

Por otro lado, cuando acerca a , el valor de acerca a desde la dirección positiva. Por ejemplo y . Por lo tanto, produce solo valores positivos y elimina el "peligro" de los grandes números negativos. Esto generalmente asegura una distribución más homogénea cuando un conjunto de datos contiene números cercanos a cero.x0ln(x+1)0ln(1+1e)0.31ln(1+1e10)0.000045log1p

En resumen, si el conjunto de datos es mayor que , generalmente está bien. Pero, si el conjunto de datos tiene números entre y , generalmente es mejor.10 1log01log1p

sfmiller940
fuente