In [1]:
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
sns.set()
In [2]:
# Create some data
np.random.seed(42)
N = 100
t = 10 * np.random.rand(N)
dy = 0.1 + 0.1 * np.random.rand(N)
y = np.random.normal(np.sin(t) + 0.5 * np.sin(2 * t), dy)
In [3]:
from bombscargle.standard import lomb_scargle, lomb_scargle_huber
from astroML.time_series import lomb_scargle as lomb_scargle_astroml
In [4]:
# Plot the data
plt.errorbar(t, y, dy, fmt='.k', ecolor='gray');
In [5]:
period = np.linspace(0.5, 20, 1000)
omega = 2 * np.pi / period
# Standard Lomb-Scargle result
PLS = lomb_scargle(t, y, dy, omega, Nterms=2,
compute_offset=False)
In [6]:
plt.plot(period, PLS)
plt.ylim(0, 1);
In [7]:
# Add some outliers
i = np.random.random(len(y))
mask = (i < 0.2)
y[mask] += 2 * np.random.randn(mask.sum())
In [8]:
# Plot the data
plt.errorbar(t, y, dy, fmt='.k', ecolor='gray');
In [9]:
period = np.linspace(0.5, 20, 1000)
omega = 2 * np.pi / period
# Standard Lomb-Scargle result
PLS = lomb_scargle(t, y, dy, omega, Nterms=1,
compute_offset=False)
# Huber Lomb-Scargle result
PLS_huber = lomb_scargle_huber(t, y, dy, omega, Nterms=1, c=5,
compute_offset=False)
PLS_astroML = lomb_scargle_astroml(t, y, dy, omega, generalized=False)
In [10]:
plt.plot(period, PLS, label='standard')
plt.plot(period, PLS_huber, label='huber')
#plt.plot(period, PLS_astroML, label='astroML')
plt.legend(frameon=False);
In [11]:
from astroML.datasets import fetch_LINEAR_sample
data = fetch_LINEAR_sample()
In [ ]:
t, y, dy = data.get_light_curve(10022663).T
plt.errorbar(t, y, dy, fmt='o');
In [ ]:
period = np.linspace(0.615, 0.618, 1000)
omega = 2 * np.pi / period
P1 = lomb_scargle(t, y, dy, omega, compute_offset=False)
P2 = lomb_scargle_huber(t, y, dy, omega, compute_offset=False)
P3 = lomb_scargle_astroml(t, y, dy, omega, generalized=False)
In [ ]:
plt.plot(period * 24, P1, label='lomb scargle')
plt.plot(period * 24, P2, label='Huber')
#plt.plot(period * 24, P3, label='astroML version')
plt.legend();
In [ ]:
t, y, dy = data[1004849].T
plt.errorbar(t, y, dy, fmt='o');
In [ ]:
period = np.linspace(0.2, 0.7, 10000)
omega = 2 * np.pi / period
P = lomb_scargle_astroml(t, y, dy, omega)
plt.plot(period, P);
In [ ]:
period = np.linspace(0.457, 0.46, 1000)
omega = 2 * np.pi / period
P = lomb_scargle_astroml(t, y, dy, omega)
plt.plot(period, P);
In [ ]:
period_best = period[np.argmax(P)]
phase = t % period_best
plt.errorbar(phase, y, dy, fmt='.', ecolor='#AABBFF')
plt.gca().invert_yaxis();
In [ ]:
P1 = lomb_scargle(t, y, dy, omega, compute_offset=False)
P2 = lomb_scargle_huber(t, y, dy, omega, compute_offset=False)
plt.plot(period * 24, P1, label='lomb scargle')
plt.plot(period * 24, P2, label='Huber')
plt.legend();
In [ ]: