docs/recipes/Synthethic_generation_of_light_curves.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Synthethic generation of light curves\n",
    "\n",
    "This notebooks illustrates the synthethic generation of light curves with mutis.\n",
    "\n",
    "There are three methods implemented with mutis, with some variations:\n",
    "- Generating signals by **sampling**. This signals have same statisical distribution.\n",
    "- Gnerating signals by randomization of the Fourier transform from the **PSD** (Power Spectral Distribution). This signals have the same variability, and same mean and std.\n",
    "- Generating signals by integration of an **stochastic** process (Orstein-Uhlenbeck process). Signals have similar shape."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# I'm doing a lot of changes to MUTIS while writting this, better reload automatically.\n",
    "%load_ext autoreload\n",
    "%autoreload\n",
    "\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "\n",
    "import matplotlib as mplt\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "from astropy.time import Time\n",
    "\n",
    "import  mutis"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Load test data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "data = {}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "data['3mm'] = pd.read_csv('data/mm-I.dat', comment='!')\n",
    "data['3mm']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "data['gamma'] = pd.read_csv('data/gamma-I.dat', comment='!')\n",
    "data['gamma']"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3mm"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import mutis\n",
    "sig3mm = mutis.Signal(data['3mm']['jyear'], data['3mm']['I'], data['3mm']['dI'])\n",
    "sig3mm.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### `lc_gen_samp`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We first check the generation with simple sampling of the signals. We see that the shape is not similar."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sig3mm.check_gen('lc_gen_samp')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### `lc_gen_psd_*`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sig3mm.check_gen('lc_gen_psd_fft')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We see that the signals have similar shape, and almost identical variability. The statistical distribution is not exactly the same.\n",
    "\n",
    "We now check using the non-uniform fourier transform. The generation with the FFT was not technically correct since the signal was not evenly sampled in time."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sig3mm.check_gen('lc_gen_psd_nft')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Still another method can be used to reconstruct the signal; the Lomb-Scargle method to compute the PSD."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sig3mm.check_gen('lc_gen_psd_lombscargle')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**From these results, we see that the shape for 3mm is better reproduced with `lc_gen_psd_fft`.**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### `lc_gen_OU`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To use the stochastic OU method, first we need to find suitable parameters:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sig3mm.OU_fit()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sig3mm.check_gen('lc_gen_ou', fgen_params={'mu':1.8, 'sigma':1.01, 'theta':1.96})"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sig3mm.OU_mu, sig3mm.OU_sigma, sig3mm.OU_theta = 1.8, 1.01, 1.96"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We see that these value produce synthethic light curves that have similar shape, distribution and variability."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Gamma"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sigGamma = mutis.Signal(data['gamma']['jyear'][np.isfinite(data['gamma']['CFlux'])], 1e6*data['gamma']['CFlux'][np.isfinite(data['gamma']['CFlux'])], 1e6*data['gamma']['CFluxErr'][np.isfinite(data['gamma']['CFlux'])])\n",
    "sigGamma.plot()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### `lc_gen_psd_*`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sigGamma.check_gen('lc_gen_psd_fft')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here we see that `lc_gen_psd_fft` does not generate signals with similar shapes, specially during flares."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sigGamma.check_gen('lc_gen_psd_nft')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### `lc_gen_OU`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sigGamma.OU_fit()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sigGamma.check_gen('lc_gen_ou',\n",
    "                   fgen_params={'mu':0.21, 'sigma':4.9, 'theta':7.6, \n",
    "                               })#'scale':np.std(sigGamma.values), 'loc':np.mean(sigGamma.values)})"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Set and generate synthethics"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now that we have checked which generation methods are the best, we set them and generate the synthethic light curves."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%time\n",
    "sig3mm.fgen = 'lc_gen_psd_fft'\n",
    "sig3mm.gen_synth(400)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%time\n",
    "sigGamma.fgen = 'lc_gen_ou'\n",
    "sigGamma.OU_mu, sigGamma.OU_sigma, sigGamma.OU_theta = 0.2, 4.9, 7\n",
    "sigGamma.gen_synth(400)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Correlation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "corr3mmGamma = mutis.Correlation(sig3mm, sigGamma, 'welsh')\n",
    "corr3mmGamma.plot_signals()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "corr3mmGamma.gen_times(ftimes='uniform', tmin=-200/365, tmax=+200/365, n=50, nbinsmin=11)\n",
    "corr3mmGamma.plot_times()\n",
    "plt.xlim([-200/365,+200/365])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "corr3mmGamma.samples = 400"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%time\n",
    "corr3mmGamma.gen_corr(uncert=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "corr3mmGamma.plot_corr(uncert=False)\n",
    "plt.xlim([-200/365,+200/365])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}