<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 WASP-80b broadband analysis
3a. Gaussian process hyperparameter estimation I
Hannu Parviainen, Instituto de Astrofísica de Canarias<br>
This notebook works as an appendix to Parviainen et al., Ground based transmission spectroscopy of WASP-80b (2017). The paper covers two analyses
Step1:  Compute residuals
The GP hyperparameters are fit to the residuals from the ckwn white-noise analysis.
Step2:  Plot the light curves
Step3:  Fit the Hyperparameters and plot the GP mean with the data
Step4:  Create a Pandas dataframe and save the hyperparameters
<ASSISTANT_TASK:>
Python Code:
%pylab inline
%run __init__.py
from exotk.utils.misc import fold
from src.extcore import *
Explanation: WASP-80b broadband analysis
3a. Gaussian process hyperparameter estimation I
Hannu Parviainen, Instituto de Astrofísica de Canarias<br>
This notebook works as an appendix to Parviainen et al., Ground based transmission spectroscopy of WASP-80b (2017). The paper covers two analyses: a broadband analysis using three previously published datasets, and a transmission spectroscopy analysis using two GTC-observed spectroscopic time series, and this notebook covers a part of the broadband analysis.
Last (significant) revision: 11.08.2017
Here we estimate the GP hyperparameters (HPs) for the light curves in T13 and M14 datasets. The GP for these two datasets uses a simple kernel with time as the only input parameter.
End of explanation
lpf = LPFTM()
pv0 = pd.read_hdf(RFILE_EXT, 'ckwn/fc').median().values
fluxes_m = lpf.compute_transit(pv0)
residuals = [fo-fm for fo,fm in zip(lpf.fluxes, fluxes_m)]
gps = [GPTime(time, res) for time,res in zip(lpf.times, residuals)]
hps = []
Explanation: Compute residuals
The GP hyperparameters are fit to the residuals from the ckwn white-noise analysis.
End of explanation
phases = list(map(lambda t: fold(t, P, TC, 0.5)-0.5, lpf.times))
fig,axs = subplots(4,3, figsize=(14,14),sharey=True, sharex=True)
for iax,ilc in enumerate(lpf.lcorder):
    a = axs.flat[iax]
    a.plot(phases[ilc], lpf.fluxes[ilc],'.', alpha=0.5)
    a.plot(phases[ilc], fluxes_m[ilc],'k')
    a.plot(phases[ilc], lpf.fluxes[ilc]-fluxes_m[ilc]+0.95,'.', alpha=0.5)
    a.text(0.5, 0.95, lpf.passbands[ilc], ha='center', va='top', size=12, transform=a.transAxes)
setp(axs, ylim=(0.94,1.01), xlim=(-0.035,0.035))
fig.tight_layout()
axs.flat[-1].set_visible(False)
Explanation: Plot the light curves
End of explanation
hps = []
for gp in tqdm(gps, desc='Optimising GP hyperparameters'):
    gp.fit()
    hps.append(gp.hp)
fig,axs = subplots(4,3, figsize=(14,10),sharey=True, sharex=True)
for iax,ilc in enumerate(lpf.lcorder):
    axs.flat[iax].plot(phases[ilc], gps[ilc].flux, '.', alpha=0.5)
    gps[ilc].compute(hps[ilc])
    pr = gps[ilc].predict()
    axs.flat[iax].plot(phases[ilc], pr, 'k')
setp(axs, ylim=(-0.015,.015), xlim=(-0.04,0.04))
fig.tight_layout()
axs.flat[-1].set_visible(False)
fig,axs = subplots(4,3, figsize=(14,10),sharey=True, sharex=True)
for iax,ilc in enumerate(lpf.lcorder):
    axs.flat[iax].plot(phases[ilc], lpf.fluxes[ilc], '.', alpha=0.5)
    gps[ilc].compute(hps[ilc])
    pr = gps[ilc].predict()
    axs.flat[iax].plot(phases[ilc], fluxes_m[ilc]+pr, 'k')
setp(axs, ylim=(0.955,1.015), xlim=(-0.04,0.04))
fig.tight_layout()
axs.flat[-1].set_visible(False)
Explanation: Fit the Hyperparameters and plot the GP mean with the data
End of explanation
with pd.HDFStore(DFILE_EXT) as f:
    ntr = [k[3:] for k in f.keys() if 'lc/triaud' in k]
    nma = [k[3:] for k in f.keys() if 'lc/mancini' in k]
df = pd.DataFrame(hps, columns=gp.names, index=lpf.passbands)
df['lc_name'] = ntr+nma
df
df.ix[:3].to_hdf(RFILE_EXT, 'gphp/triaud2013')
df.ix[3:].to_hdf(RFILE_EXT, 'gphp/mancini2014')
Explanation: Create a Pandas dataframe and save the hyperparameters
End of explanation
<END_TASK>