numero - acos(1) devuelve NaN para algunos valores, no otros

math sign

Question

Tengo una lista de valores de latitud y longitud, y estoy tratando de encontrar la distancia entre ellos. Utilizando un método estándar de gran círculo, necesito encontrar:

acos(sin(lat1)*sin(lat2) + cos(lat1)*cos(lat2) * cos(long2-long1))

Y multiplica esto por el radio de la tierra, en las unidades que estoy usando. Esto es válido siempre que los valores que tomamos el acos de estén en el rango [-1,1]. Si están ligeramente fuera de este rango, devolverá NaN , incluso si la diferencia se debe al redondeo.

El problema que tengo es que a veces, cuando dos valores lat / long son idénticos, esto me da un error NaN . No siempre, incluso para el mismo par de números, pero siempre los mismos en una lista. Por ejemplo, tengo una persona parada en una carretera en el desierto:

Time  |lat     |long
1:00PM|35.08646|-117.5023
1:01PM|35.08646|-117.5023
1:02PM|35.08646|-117.5023
1:03PM|35.08646|-117.5023
1:04PM|35.08646|-117.5023

Cuando calculo la distancia entre los puntos consecutivos, el tercer valor, por ejemplo, siempre será NaN , aunque los otros no lo sean. Esto parece ser un error extraño con el redondeo R.

Answer 1

No puedo decir exactamente sin ver tus datos (prueba dput ), pero esto es probablemente una consecuencia de la Pregunta frecuente 7.31 .

(x1 <- 1)
## [1] 1
(x2 <- 1+1e-16)
## [1] 1
(x3 <- 1+1e-8)
## [1] 1
acos(x1)
## [1] 0
acos(x2)
## [1] 0
acos(x3)
## [1] NaN

Es decir, incluso si sus valores son tan similares que sus representaciones impresas son las mismas, aún pueden diferir: algunas estarán dentro de .Machine$double.eps y otras no ...

Una forma de asegurarse de que los valores de entrada estén limitados por [-1,1] es usar pmax y pmin : acos(pmin(pmax(x,-1.0),1.0))

Answer 2

Una solución simple es usar pmin (), como este:

acos(pmin(sin(lat1)*sin(lat2) + cos(lat1)*cos(lat2) * cos(long2-long1),1))

Ahora se asegura de que la pérdida de precisión conduzca a un valor no mayor que exactamente 1.

Sin embargo, esto no explica lo que está sucediendo.

(Editar: Matthew Lundberg señaló que necesito usar pmin para que funcione con entradas vectorizadas. Esto soluciona el problema de hacer que funcione, pero todavía no estoy seguro de por qué está redondeando incorrectamente).