<jupyter_start><jupyter_text>Simulating Microlensing graphs$I_{obs} \left( t \right) = I_{s0} + I_{A} + N \left( \sigma^{2} \right) $$I_{s0}$ is baseline, $I_{A}$ is lens magnification<jupyter_code>import numpy as np
from matplotlib import pyplot as plt<jupyter_output><empty_output><jupyter_text>Baseline $I_{s0} = m_{A} \mbox{sin} \left( \frac{2 \pi t}{T} + \phi \right) $where$m_{A}$ is amplitude$T$ is period$\phi$ is phase<jupyter_code>def intensity_baseline(period, amplitude, phase, t):
  return amplitude * np.sin((2 * np.pi * t / period) + phase)
is0_ex1 = [intensity_baseline(5, 0.00005, 0, x) for x in range(100)]
plt.plot(is0_ex1)<jupyter_output><empty_output><jupyter_text>Lens Magnification $A \left( t \right) = \frac{\left[ u(t) \right]^{2} + 2}{u(t) + \sqrt{[u(t)]^2 + 4}} $$u(t) = u_0 + \left| \frac{t-t_0}{T_E} \right|$where$u_0$ is minimum impact parameter$T_E$ and $t_0$ determins the radius and shift in time<jupyter_code>def lens_magnification(u0, t0, tE, t):
  ut = u0 + np.abs((t - t0) / tE)
  return (ut ** 2 + 2) / (ut * np.sqrt(ut**2 + 4))
iA_ex1 = [lens_magnification(10, 50, 2, x) for x in range(100)]
plt.plot(iA_ex1)<jupyter_output><empty_output><jupyter_text>Putting all together<jupyter_code>iobs_ex1 = np.array(iA_ex1) + np.array(is0_ex1) + np.random.normal(scale=0.00003, size=100)
plt.plot(iobs_ex1)<jupyter_output><empty_output><jupyter_text>Perhaps its a good idea to normalize the graph so that it will have mean=0 and stdev=1<jupyter_code>iobs_ex1_normalized = iobs_ex1 / np.linalg.norm(iobs_ex1)
plt.plot(iobs_ex1_normalized)
import numpy as np
import matplotlib.pyplot as plt

def intensity_baseline(period, amplitude, phase, t):
  return amplitude * np.sin((2 * np.pi * t / period) + phase)


def lens_magnification(u0, t0, tE, t):
  ut = u0 + np.abs((t - t0) / tE)
  return (ut ** 2 + 2) / (ut * np.sqrt(ut**2 + 4))


def simulate_microlensing(baseline_period, baseline_amplitude, baseline_phase,
                          lens_min_impact, lens_shift, lens_radius, noise, t_range=range(100)):
  i_s0 = np.array([intensity_baseline(baseline_period, baseline_amplitude, baseline_phase, t) for t in t_range])
  i_A = np.array([lens_magnification(lens_min_impact, lens_shift, lens_radius, t) for t in t_range])
  n = np.random.normal(scale=noise, size=len(t_range))
  return i_s0 + i_A + n
plt.plot(simulate_microlensing(baseline_period=2, baseline_amplitude=0.00001, baseline_phase=0,
                               lens_min_impact=10, lens_shift=75, lens_radius=1, noise=0.00003, t_range=range(100)))<jupyter_output><empty_output><jupyter_text>Fitting sampled data into ARIMA Installing pyflux (library for time series analysis). Will be using ARIMA class for fitting the sampled data into ARIMA model and make prediction.<jupyter_code>pip install pyflux 
import pyflux as pf<jupyter_output><empty_output><jupyter_text>Creating the model<jupyter_code>data = simulate_microlensing(baseline_period=2, baseline_amplitude=0.00001, baseline_phase=0,
                               lens_min_impact=10, lens_shift=75, lens_radius=1, noise=0.00003, t_range=range(100))

model = pf.ARIMA(data=data, ar=10, integ=1, ma=10, target='microlense', family=pf.Normal())<jupyter_output><empty_output><jupyter_text>Fitting the data and getting the summary<jupyter_code>fit = model.fit("MLE")
fit.summary()<jupyter_output>Normal ARIMA(10,1,10)                                                                                     
======================================================= ==================================================
Dependent Variable: Differenced Series                  Method: MLE                                       
Start Date: 10                                          Log Likelihood: 772.7007                          
End Date: 99                                            AIC: -1501.4013                                   
Number of observations: 89                              BIC: -1446.6513                                   
==========================================================================================================
Latent Variable                          Estimate   Std Error  z        P>|z|    95% C.I.                 
======================================== ========== ========== ======== ======== =========================
Constant                             [...]<jupyter_text>Plotting the fit-<jupyter_code>model.plot_fit(figsize=(15,10))
model.plot_predict(h=20,past_values=20,figsize=(15,5))
np.array(model.predict(h=70)['Differenced Series'])
model = pf.ARIMA(data=data, ar=4, integ=1, ma=4, target='microlense', family=pf.Normal())
model.fit("MLE")
model.plot_predict_is(h=70)<jupyter_output><empty_output><jupyter_text>data_with_peak = simulate_microlensing(baseline_period=2, baseline_amplitude=0.00001, baseline_phase=0,                               lens_min_impact=10, lens_shift=75, lens_radius=1, noise=0.00003, t_range=range(100))<jupyter_code>model = pf.ARIMA(data=data_with_peak, ar=4, integ=1, ma=4, target='microlense', family=pf.t())
model.fit("MLE")
model.plot_predict_is(h=70)
plt.plot(data_with_peak)<jupyter_output><empty_output><jupyter_text>Generate Train and Test Data<jupyter_code>def simulate_noevent(baseline_period, baseline_amplitude, baseline_phase, noise, t_range=range(100)):
    i_s0 = np.array([intensity_baseline(baseline_period, baseline_amplitude, baseline_phase, t) for t in t_range])
    n = np.random.normal(scale=noise, size=len(t_range))
    return i_s0 + n
from tqdm import tqdm

def generate_dataset(num_dataset):
    X_out = np.zeros((num_dataset, 100, 2))
    y_out = np.zeros((num_dataset, 100))

    for x in tqdm(range(num_dataset)):
        seq=None
        ys=None
        if np.random.uniform() > 0.5:
            peak_t = int(np.random.uniform(low=50, high=80))
            seq = simulate_microlensing(baseline_period=2, 
                                        baseline_amplitude=0.00001, 
                                        baseline_phase=0,
                                        lens_min_impact=10, 
                                        lens_shift=peak_t, 
                                        lens_radius=1, 
                                        noise=0.00003, 
                                        t_range=range(100))
            seq_avg = np.mean(seq[:30])
            seq_std = np.std(seq[:30])
            seq = (seq - seq_avg) / seq_std
            ys = np.zeros(100)
            ys[peak_t-10:peak_t] = 1
        else:
            seq = simulate_noevent(baseline_period=2, 
                                baseline_amplitude=0.00001,
                                baseline_phase=0,
                                noise=0.00003,
                                t_range=range(100))
            seq_avg = np.mean(seq[:30])
            seq_std = np.std(seq[:30])
            seq = (seq - seq_avg) / seq_std
            ys = np.zeros(100)
        model=pf.ARIMA(data=seq,ar=4, integ=1, ma=4, family=pf.Normal())
        model.fit('MLE')
        hs = model.predict_is(h=70)
        seq_sim = np.concatenate((seq[:30], hs['Differenced Series']))
        X_out[x,:,0] = seq
        X_out[x,:,1] = seq_sim
        y_out[x] = ys
    
    return X_out, y_out
np.random.seed(420)
X_train, y_train = generate_dataset(100)
np.random.seed(520)
X_test, y_test = generate_dataset(10)<jupyter_output>100%|██████████| 10/10 [00:24<00:00,  2.49s/it]<jupyter_text>Injecting the residuals into LSTM  Building the LSTM model<jupyter_code>from tensorflow.keras import backend as K

def recall_m(y_true, y_pred):
        true_positives = K.sum(K.round(K.clip(y_true * y_pred, 0, 1)))
        possible_positives = K.sum(K.round(K.clip(y_true, 0, 1)))
        recall = true_positives / (possible_positives + K.epsilon())
        return recall

def precision_m(y_true, y_pred):
        true_positives = K.sum(K.round(K.clip(y_true * y_pred, 0, 1)))
        predicted_positives = K.sum(K.round(K.clip(y_pred, 0, 1)))
        precision = true_positives / (predicted_positives + K.epsilon())
        return precision
    
def f1_m(y_true, y_pred):
    precision = precision_m(y_true, y_pred)
    recall = recall_m(y_true, y_pred)
    return 2*((precision*recall)/(precision+recall+K.epsilon()))
%tensorflow_version 2.x
import tensorflow as tf
from tensorflow.keras import layers
# Initialising the RNN
rnn_lstm = tf.keras.Sequential()
# Adding the LSTM layers and some Dropout regularisation

# Adding the first layer
rnn_lstm.add(layers.LSTM(units=300, return_sequences=True, input_shape=(None, 2)))
rnn_lstm.add(layers.Dropout(0.1))

# Output layer
rnn_lstm.add(layers.Dense(units=1, activation='sigmoid'))
# Compiling the RNN
rnn_lstm.compile(optimizer='adam', loss='mean_squared_error', metrics=['acc', f1_m])
rnn_lstm.summary()
# Fitting the RNN to training set
rnn_lstm.fit(X_train, y_train, batch_size=32, epochs=10, validation_data=(X_test, y_test))
rnn_lstm.evaluate(X_test, y_test)
sum(sum(y_test)) / 1000
y_pred=rnn_lstm.predict_proba(X_test)
plt.plot(y_pred[4][30:])
plt.plot(X_test[4,:,0][30:])<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.15585416896133936, 'fraction_numerical': 0.03910122260160811, 'mean_word_length': 3.0864086408640863, 'pattern_counts': {'":': 0, '<': 48, '<?xml version=': 0, '>': 49, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 16, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '12', 'dir': 'jupyter-structured-clean-dedup', 'id': '296595', 'license': 'MIT', 'n_tokens_mistral': 3399, 'n_tokens_neox': 3082, 'n_words': 590, 'path': 'simulate_curve_lstm_gru.ipynb', 'repo_name': 'dibyajit30/ExoplanetHunting', 'size': '9079'}