<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 Brian Larsen - 28 April 2016 balarsen@lanl.gov
Overview
This is meant to be a proof of concept example of the SWx to Mission instrument prediction planned for AFTAC funding.
Abstract
There is currently no accepted method within SNDD for the time dependent quantification of environmental background in this mission instruments, this work aims to provide this capability. This work provides a preliminary proof of concept of the techniques in order to both show utility and provide a baseline for improvement. The methods incorporated here are Bayesian MCMC curve fitting and Bayesian posterior sampling. The method is to utilize data of SWx and mission instruments to find (linear) correlation coefficients using all the data. Then utilizing the relations discovered in the correlation step the predicted response in the mission instruments are then created. 
Goals
Provide a quantifiable prediction of count rate in the mission instruments
Need to study which data are required to provide the best prediction
First step is beginning with data moving toward prediction, second step is a first principle assessment of background based on SAGE
Provide an automatic method to determine a quantifiable deviation from the prediction in order to determine Events of Interest (EOI)
Provide an update and future capability into the method for update and reassessment of accuracy and addition of additional data 
Provide mechanisms for the inclusion of the model and prediction rules into DIORAMA without DIORAMA having to run MCMC.
The Steps are
Perform a correlation analysis between the two data sets figuring out which variables make the best prediction
Use the results of this correlation to predict the response for the Mission instrument from the SWx
Utilize the predicted response to determine quantifiable count rate prediction limits
Develop a technique to quantify excursions of the observed data from the predicted limits to flag as EOI with an interest score
Step1:  Create simulated data that can be used in this proof of concept
Step2:  Now based on the results above we can use this as a prediction of the Y-data from the X-data into the future
Now that we have time dependent data use the results from the above correlation to predict what the Y-instrument would have seen
Step3:  Or do this smarter by sampling the posterior
<ASSISTANT_TASK:>
Python Code:
!date
# Standard Library Python Modules
# Common Python Modules
import matplotlib.pyplot as plt
import numpy as np
import spacepy.plot as spp
import spacepy.toolbox as tb
import pandas as pd
import pymc # this is the MCMC tool
# put plots into this document
%matplotlib inline
Explanation: Brian Larsen - 28 April 2016 balarsen@lanl.gov
Overview
This is meant to be a proof of concept example of the SWx to Mission instrument prediction planned for AFTAC funding.
Abstract
There is currently no accepted method within SNDD for the time dependent quantification of environmental background in this mission instruments, this work aims to provide this capability. This work provides a preliminary proof of concept of the techniques in order to both show utility and provide a baseline for improvement. The methods incorporated here are Bayesian MCMC curve fitting and Bayesian posterior sampling. The method is to utilize data of SWx and mission instruments to find (linear) correlation coefficients using all the data. Then utilizing the relations discovered in the correlation step the predicted response in the mission instruments are then created. 
Goals
Provide a quantifiable prediction of count rate in the mission instruments
Need to study which data are required to provide the best prediction
First step is beginning with data moving toward prediction, second step is a first principle assessment of background based on SAGE
Provide an automatic method to determine a quantifiable deviation from the prediction in order to determine Events of Interest (EOI)
Provide an update and future capability into the method for update and reassessment of accuracy and addition of additional data 
Provide mechanisms for the inclusion of the model and prediction rules into DIORAMA without DIORAMA having to run MCMC.
The Steps are
Perform a correlation analysis between the two data sets figuring out which variables make the best prediction
Use the results of this correlation to predict the response for the Mission instrument from the SWx
Utilize the predicted response to determine quantifiable count rate prediction limits
Develop a technique to quantify excursions of the observed data from the predicted limits to flag as EOI with an interest score
End of explanation
# observed data
from scipy.signal import savgol_filter
# make a time dependent x
t_x = np.random.uniform(0, 1, 500)
t_x = savgol_filter(t_x, 95, 2)[95//2:-95//2]
t_x -= t_x.min()
t_x /= (t_x.max() - t_x.min())
plt.plot(t_x)
plt.xlabel('time')
plt.ylabel('SWx data value')
a = 6
b = 2
sigma = 2.0
y_obs = a*t_x + b + np.random.normal(0, sigma, len(t_x))
data = pd.DataFrame(np.array([t_x, y_obs]).T, columns=['x', 'y'])
x = t_x
data.plot(x='x', y='y', kind='scatter', s=50)
plt.xlabel('SWx inst')
plt.ylabel('Mission Inst')
# define priors
a = pymc.Normal('slope', mu=0, tau=1.0/10**2)
b = pymc.Normal('intercept', mu=0, tau=1.0/10**2)
tau = pymc.Gamma("tau", alpha=0.1, beta=0.1)
# define likelihood
@pymc.deterministic
def mu(a=a, b=b, x=x):
    return a*x + b
y = pymc.Normal('y', mu=mu, tau=tau, value=y_obs, observed=True)
# inference
m = pymc.Model([a, b, tau, x, y])
mc = pymc.MCMC(m)
# run 6 chains
for i in range(6):
    mc.sample(iter=90000, burn=10000)
# plot up the data and overplot the possible fit lines
data.plot(x='x', y='y', kind='scatter', s=50)
xx = np.linspace(data.x.min(), data.x.max(), 10)
for ii in range(0, len(mc.trace('slope', chain=None)[:]), 
                len(mc.trace('slope', chain=None)[:])//400):
    yy = (xx*mc.trace('slope', chain=None)[:][ii] + 
          mc.trace('intercept', chain=None)[:][ii])
    plt.plot(xx,yy, c='r')
    
pymc.Matplot.plot(mc)
pymc.Matplot.summary_plot(mc)
Explanation: Create simulated data that can be used in this proof of concept
End of explanation
int_vals = mc.stats()['intercept']['95% HPD interval']
slope_vals = mc.stats()['slope']['95% HPD interval']
print(int_vals, slope_vals)
y_inst =  np.tile(t_x, (2,1)).T * slope_vals + int_vals
plt.plot(y_inst, c='r')
plt.xlabel('time')
plt.ylabel('Mission inst value')
plt.title('Major upper limit of spread')
Explanation: Now based on the results above we can use this as a prediction of the Y-data from the X-data into the future
Now that we have time dependent data use the results from the above correlation to predict what the Y-instrument would have seen
End of explanation
pred = []
for v in t_x:
    pred.append(np.percentile(v * mc.trace('slope', 
                                           chain=None)[:] + 
                              mc.trace('intercept', 
                                       chain=None)[:], 
                              [2.5, 97.5]))
plt.plot(pred, c='r')
plt.xlabel('time')
plt.ylabel('Predicted Mission inst value')
Explanation: Or do this smarter by sampling the posterior
End of explanation
<END_TASK>