Python vs Julia autocorrelación

Estoy tratando de hacer autocorrelación con Julia y compararlo con el resultado de Python. ¿Cómo es que dan resultados diferentes?

Código de Julia

using StatsBase

t = range(0, stop=10, length=10)
test_data = sin.(exp.(t.^2))

acf = StatsBase.autocor(test_data)

da

10-element Array{Float64,1}:
  1.0                   
  0.13254954979179642   
 -0.2030283419321465    
  0.00029587850872956104
 -0.06629381497277881   
  0.031309038331589614  
 -0.16633393452504994   
 -0.08482388975165675   
  0.0006905628640697538 
 -0.1443650483145533

Código de Python

from statsmodels.tsa.stattools import acf
import numpy as np

t = np.linspace(0,10,10)
test_data = np.sin(np.exp(t**2))

acf_result = acf(test_data)

da

array([ 1.        ,  0.14589844, -0.10412699,  0.07817509, -0.12916543,
       -0.03469143, -0.129255  , -0.15982435, -0.02067688, -0.14633346])

Esto se debe a que tu test_dataes diferente:

Pitón:

array([ 0.84147098, -0.29102733,  0.96323736,  0.75441021, -0.37291918,
        0.85600145,  0.89676529, -0.34006519, -0.75811102, -0.99910501])

Julia

[0.8414709848078965, -0.2910273263243299, 0.963237364649543, 0.7544102058854344,
 -0.3729191776326039, 0.8560014512776061, 0.9841238290665676, 0.1665709194875013,
 -0.7581110212957692, -0.9991050130774393]

Esto sucede porque estás tomando sincantidades enormes. Por ejemplo, con el último número en t10, exp(10^2)es ~ 2.7 * 10 ^ 43. En esta escala, las imprecisiones de coma flotante son aproximadamente 3 * 10 ^ 9. Entonces, incluso si el bit menos significativo es diferente para Python y Julia, el sinvalor estará muy lejos.

De hecho, podemos inspeccionar los valores binarios subyacentes de la matriz inicial t. Por ejemplo, difieren en el tercer último valor:

Julia

julia> reinterpret(Int, range(0, stop=10, length=10)[end-2])
4620443017702830535

Pitón:

>>> import struct
>>> s = struct.pack('>d', np.linspace(0,10,10)[-3])
>>> struct.unpack('>q', s)[0]
4620443017702830536

De hecho, podemos ver que no están de acuerdo con exactamente una máquina épsilon. Y si usamos Julia tomamos sinel valor obtenido por Python:

julia> sin(exp(reinterpret(Float64, 4620443017702830536)^2))
-0.3400651855865199

Obtenemos el mismo valor que Python.

Solo para expandir un poco la respuesta (agregar como respuesta ya que es demasiado larga para un comentario). En Julia tienes lo siguiente:

julia> t = collect(range(0, stop=10, length=10))
10-element Array{Float64,1}:
  0.0               
  1.1111111111111112
  2.2222222222222223
  3.3333333333333335
  4.444444444444445 
  5.555555555555555 
  6.666666666666667 
  7.777777777777778 
  8.88888888888889  
 10.0               

julia> t .- [10*i / 9 for i in 0:9]
10-element Array{Float64,1}:
 0.0
 0.0
 0.0
 0.0
 0.0
 0.0
 0.0
 0.0
 0.0
 0.0

mientras en Python:

>>> t = np.linspace(0,10,10)
>>> t - [10*i/9 for i in range(10)]
array([0.0000000e+00, 0.0000000e+00, 0.0000000e+00, 0.0000000e+00,
       0.0000000e+00, 0.0000000e+00, 0.0000000e+00, 8.8817842e-16,
       0.0000000e+00, 0.0000000e+00])

y ves que el octavo número en Python es una aproximación inexacta de 70/9, mientras que en Julia en este caso obtienes la secuencia de aproximaciones más cercanas de 10*i/9uso Float64.

Entonces parece que debido a que las secuencias originales difieren, el resto sigue lo que comentó @Jakob Nissen.

Sin embargo, las cosas no son tan simples. Como las expfunciones en Julia y Python difieren un poco en lo que producen. Ver Python:

>>> from math import exp
>>> from mpmath import mp
>>> mp.dps = 1000
>>> float(mp.exp((20/3)**2) - exp((20/3)**2))
-1957.096392544307

mientras en Julia:

julia> setprecision(1000)
1000

julia> Float64(exp(big((20/3)^2)) - exp((20/3)^2))
2138.903607455693

julia> Float64(exp(big((20/3)^2)) - nextfloat(exp((20/3)^2)))
-1957.096392544307

(puede verificar que (20/3)^2es lo mismo Float64tanto en Julia como en Python).

Entonces, en este caso, expPython es un poco más preciso que Julia. Por lo tanto, incluso la fijación t(que es fácil usando una comprensión en Python en lugar de linspace) no hará que el ACF sea igual.

En general, la conclusión es lo que @Jakob Nissen comentó para valores tan grandes que los resultados estarán fuertemente influenciados por las imprecisiones numéricas.