<jupyter_start><jupyter_text>K2c02 YSO Summary Statistics *Clean*  Calculate the summary statistics for all 1678 YSO K2 C2 lightcurves, but in a restricted range of time to avoid instrumental artifactsMichael Gully-Santiago  January-February 2016There are 1678 light curve files for sources that were proposed for by Kraus, Hillenbrand, Covey, or (?) to be YSOs or YSO candidates in K2 Cycle 2.  See my other notebook for how I downloaded them.  Once we have the files it's a matter of what to do with them.  Here we read in each file and construct the basic summary statistics using the pandas `describe()` function.  We then save that series and some metadata into a data frame for each lightcurve.  **The result is a $1658 \times 12 $ dataset of summary statistics for each K2 C2 lightcurve.**  We had to drop 20 lightcurves that had NaN's for some reason.Finally, we construct some plots of the new data just to explore it.<jupyter_code>%matplotlib inline
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
%config InlineBackend.figure_format = 'retina' 
from astropy.utils.console import ProgressBar
import warnings<jupyter_output><empty_output><jupyter_text>Make a list of all the lightcurve filenames.There's a wget script in the `data` directory that will automatically fetch each file. ```python!ls ../data/hl* > file_list.csv```<jupyter_code>file_list_raw = pd.read_csv('file_list.csv', names=['fname'])
file_list = file_list_raw.fname.values<jupyter_output><empty_output><jupyter_text>Compute the descriptive statistics for all 1678 of the YSO lightcurves from K2 Cycle 02.<jupyter_code>columns = ['count', 'mean', 'std', 'min', '25%', '50%', '75%', 'max']
agg_data = pd.DataFrame(columns=columns)

n_files = len(file_list)<jupyter_output><empty_output><jupyter_text>Clean the signal by taking only a subset of the lightcurve. There's a big discontinuous jump in many of the light curves, close to the middle.  I think it's due to the spacecraft momentum dumps.  We should just keep the first half for now.<jupyter_code>file = file_list[5]
raw_lc = pd.read_csv(file, index_col=False)
subset_i = (raw_lc['BJD - 2454833'] > 2063.0) & (raw_lc['BJD - 2454833'] < 2102.0)
plt.plot(raw_lc['BJD - 2454833'][subset_i], raw_lc[' Corrected Flux'][subset_i])<jupyter_output><empty_output><jupyter_text>This takes about 1-2 minutes:<jupyter_code>with ProgressBar(n_files, ipython_widget=True) as bar:
    for i in range(n_files):
        bar.update()
        file = file_list[i]
        raw_lc = pd.read_csv(file, index_col=False)
        #subset_i = (raw_lc['BJD - 2454833'] > 2063.0) & (raw_lc['BJD - 2454833'] < 2102.0)
        agg_data.loc[i] = raw_lc[' Corrected Flux'].describe()
cln_data = pd.DataFrame(columns=columns)
with ProgressBar(n_files, ipython_widget=True) as bar:
    for i in range(n_files):
        bar.update()
        file = file_list[i]
        raw_lc = pd.read_csv(file, index_col=False)
        subset_i = (raw_lc['BJD - 2454833'] > 2063.0) & (raw_lc['BJD - 2454833'] < 2102.0)
        cln_data.loc[i] = raw_lc[' Corrected Flux'][subset_i].describe()<jupyter_output><jupyter_text>Voilà:<jupyter_code>agg_data['fname'] = file_list_raw.fname
agg_data['EPIC_ID'] = agg_data['fname'].str[33:33+9]
agg_data.head(3)
cln_data['fname'] = file_list_raw.fname
cln_data['EPIC_ID'] = cln_data['fname'].str[33:33+9]
cln_data.head(3)<jupyter_output><empty_output><jupyter_text>Take the log of the std and IQR.<jupyter_code>agg_data['logstd'] = np.log10(agg_data['std'])
agg_data['logiqr'] = np.log10(agg_data['75%']-agg_data['25%'])
cln_data['logstd'] = np.log10(cln_data['std'])
cln_data['logiqr'] = np.log10(cln_data['75%']-cln_data['25%'])<jupyter_output><empty_output><jupyter_text>Save the cleaned data.<jupyter_code>new_col_order = ['EPIC_ID','count', 'mean', 'std', 'min', '25%', '50%', '75%', 'max', 'logstd', 'logiqr','fname']
cln_out = cln_data[new_col_order]
#cln_out.head()
cln_out.to_csv('../analysis/K2C02_YSO_cln_data.csv', index=False)<jupyter_output><empty_output><jupyter_text>Plot the standard deviation $\sigma$ of the light curve versus its interquartile range ($Q_3-Q_1$). HACK TO Visualize the cleaned up data!<jupyter_code>#agg_data = cln_data
tp = sns.jointplot(x="logiqr", y="logstd", alpha=0.2, kind='scatter',data=agg_data, xlim=(-4, 0), ylim=(-4, 0));
import bokeh
#tp = sns.jointplot(x="logiqr", y="logstd", alpha=0.2, kind='scatter',data=agg_data, xlim=(-4, 2));
sns.jointplot(x="logiqr", y="logstd", alpha=0.2, kind='scatter',data=cln_data, xlim=(-4, 0), ylim=(-4, 0));<jupyter_output>//anaconda/lib/python3.4/site-packages/matplotlib/collections.py:590: FutureWarning: elementwise comparison failed; returning scalar instead, but in the future will perform elementwise comparison
  if self._edgecolors == str('face'):<jupyter_text>The scatter of standard deviation above the interquartile range reflects sources with extra flicker noise beyond just secular sinusoidal variation.<jupyter_code>ad = agg_data.dropna()
agg_data.shape, ad.shape<jupyter_output><empty_output><jupyter_text>I had to drop 20 lightcurves with `NaN` values.  Dunno what the deal is there.<jupyter_code>#sns.distplot(ad['logiqr'].values, label='IQR')<jupyter_output><empty_output><jupyter_text>Search for sources in a specific range of variability.<jupyter_code>gi = (ad.logiqr < - 1.1) & (ad.logiqr > -1.2) & (ad.logstd <-1.4)
tp = sns.jointplot(x="logiqr", y="logstd", data=ad, xlim=(-2.0, -1.0), ylim=(-2.0, -1.0), alpha=0.5)
tp.ax_joint.plot(ad.logiqr[gi], ad.logstd[gi], '.', alpha=1.0)
print("There are {} sources that meet our above selection criterion".format(gi.sum()))
ad.EPIC_ID[gi].head()
sns.set_context('paper', font_scale=1.3)
sns.set_style('ticks')<jupyter_output><empty_output><jupyter_text>Spot-check the light curves by eye-- look for sinusoidal variations.<jupyter_code>ii = 226

plt.figure(figsize=(10, 4))
file = ad.fname[ii]
raw_lc = pd.read_csv(file, index_col=False)

plt.plot(raw_lc['BJD - 2454833'], raw_lc[' Corrected Flux'])
plt.ylim(0.90, 1.1)
plt.title("EPIC {} K2 C02 Light Curve".format(ad.EPIC_ID[ii]))
plt.xlabel('BJD - 2454833')
plt.ylabel('Flux');<jupyter_output><empty_output><jupyter_text>Good candidates:  226, 271, 313, 403, 474 Good, but messy candidates:421, 440  Save the dataframe.<jupyter_code>new_col_order = ['EPIC_ID','count', 'mean', 'std', 'min', '25%', '50%', '75%', 'max', 'logstd', 'logiqr','fname']
ad = ad[new_col_order]
ad.head(2)
!mkdir ../analysis
ad.to_csv('../analysis/K2C02_YSO_agg_data.csv', index=False)<jupyter_output><empty_output>
StarcoderdataJupyterStructuredCleanDedup