<jupyter_start><jupyter_text>Time Series with Python© <NAME>, 2016  Classification- generating data- feature engineering- rolling windows<jupyter_code>import pandas as pd
import numpy as np
%matplotlib inline
import matplotlib.pyplot as plt
from pylab import rcParams
rcParams['figure.figsize'] = 11, 7
import random
np.random.seed(41)
random.seed(41)<jupyter_output><empty_output><jupyter_text>Motivation:![Spinbot](./spinbot_timeseries.png)  Generate data- smooth curve or some sort<jupyter_code>def gen_curve(l=400, repeats=1):
    x = np.linspace(0, 2*np.pi*repeats, l*repeats)
    lab = np.array([False]*l*repeats)
    coeffs = np.random.uniform(0, np.pi, size=(4,4))
    y = np.zeros(l*repeats)
    for i, c in enumerate(coeffs):
        y += c[0] + c[1]*np.sin(c[2]*x/(0.5*np.pi) + c[3])
    return y
plt.plot(gen_curve())
plt.plot(gen_curve())
plt.plot(gen_curve())<jupyter_output><empty_output><jupyter_text>Generate Anomaly- higher frequency- localized<jupyter_code>def gen_anomaly(l=400, duration=40, scale=1.0):
    start = np.random.randint(0, l-duration)
    y = np.zeros(l)
    y[start:start+duration] = np.hanning(duration)
    x = np.linspace(0, 2*np.pi, l)
    freq = np.random.uniform(300, 800)
    z = np.sin(freq*x/(2*np.pi))
    return scale*y*z
plt.plot(gen_anomaly())
plt.plot(gen_anomaly())
plt.plot(gen_anomaly())<jupyter_output><empty_output><jupyter_text>Combine the two<jupyter_code>plt.plot(gen_curve() + gen_anomaly())<jupyter_output><empty_output><jupyter_text>Generate dataset- 10000 curves- all with the same length- 50% with anomaly- 50% without anomaly<jupyter_code>sequences = []
y = []
for i in xrange(5000):
    y0 = gen_curve()
    y1 = gen_curve() + gen_anomaly()
    sequences.append(y0)
    sequences.append(y1)
    y.extend([False, True])

sequences = pd.DataFrame(sequences)
y = np.array(y)
sequences.head()<jupyter_output><empty_output><jupyter_text>Inspect first 4 curves<jupyter_code>sequences.iloc[:4].transpose().plot()<jupyter_output><empty_output><jupyter_text>Question:how do we build features to detect this particular type of anomaly?  Start simple: rescale<jupyter_code>from sklearn.preprocessing import MinMaxScaler
mms = MinMaxScaler()
sequences_scaled = pd.DataFrame(mms.fit_transform(sequences.transpose()).transpose())
sequences_scaled.iloc[:4].transpose().plot()<jupyter_output><empty_output><jupyter_text>First idea: use standard deviation<jupyter_code>X = pd.DataFrame(sequences_scaled.std(axis=1), columns = ['std'])
from sklearn.linear_model import LogisticRegression
from sklearn.cross_validation import cross_val_score

model = LogisticRegression()
cross_val_score(model, X, y)<jupyter_output><empty_output><jupyter_text>=> better than random, but not so great actually  Histogram<jupyter_code>X[y==False]['std'].hist(bins=50, range=(0, 0.5), alpha =0.3)
X[y==True]['std'].hist(bins=50, range=(0, 0.5), alpha =0.3)
plt.legend(['Normal', 'Anomalous'])<jupyter_output><empty_output><jupyter_text>Next idea: standard deviation of differences<jupyter_code>X['std_diff'] = sequences_scaled.diff(axis=1).std(axis=1)
cross_val_score(model, X, y)
X[y==False]['std_diff'].hist(bins=50, range=(-0.1, 0.1), alpha =0.3)
X[y==True]['std_diff'].hist(bins=50, range=(-0.1, 0.1), alpha =0.3)
plt.legend(['Normal', 'Anomalous'])<jupyter_output><empty_output><jupyter_text>Next idea: use a different model<jupyter_code>from sklearn.tree import DecisionTreeClassifier
model = DecisionTreeClassifier()
cross_val_score(model, X, y)<jupyter_output><empty_output><jupyter_text>That was easy...  What if there's noise?<jupyter_code>def gen_noise(l=400, scale=None):
    if scale is not None:
        noise = np.random.normal(scale=scale, size=l)
    else:
        noise = np.random.normal(scale=np.abs(np.random.normal(scale=0.1)), size=l)
    return noise
sequences = []
y = []
for i in xrange(5000):
    y0 = gen_curve() + gen_noise()
    y1 = gen_curve() + gen_anomaly() + gen_noise()
    sequences.append(y0)
    sequences.append(y1)
    y.extend([False, True])

sequences = pd.DataFrame(sequences)
y = np.array(y)
sequences_scaled = pd.DataFrame(mms.fit_transform(sequences.transpose()).transpose())
sequences_scaled.iloc[:4].transpose().plot()
X = pd.DataFrame(sequences_scaled.std(axis=1), columns = ['std'])
X['std_diff'] = sequences_scaled.diff(axis=1).std(axis=1)
cross_val_score(model, X, y)<jupyter_output><empty_output><jupyter_text>Exercise:- can you think of any other feature to add?- try with some other models and see if you can improve the score  Fourier features to the rescueIn this example, the anomaly has a higher frequency than the base signal. Could we isolate it based on the frequency spectrum?<jupyter_code>from scipy.signal import welch
fff, Pxx = welch(sequences_scaled, nperseg=64, axis=1)
X = pd.DataFrame(Pxx)
X.head()<jupyter_output><empty_output><jupyter_text>First few curves<jupyter_code>X.iloc[:4].transpose().plot(logy=True)
plt.legend(['Normal', 'Anomalous', 'Normal', 'Anomalous'])
cross_val_score(model, X, y)
model.fit(X, y)<jupyter_output><empty_output><jupyter_text>Anomaly detectionLet's see if we can use the model we just trained to detect anomalies in a long sequence.<jupyter_code># generate sequence
long_series = gen_curve(l=400, repeats=250) + gen_noise(l=400*250)

# pick 20 indices at random
anom_indexes = np.array(sorted(np.random.choice(xrange(len(long_series-400)), size=20)))

# insert anomaly at each index
for i in anom_indexes:
    anomaly = gen_anomaly()
    l = len(anomaly)
    long_series[i:i+l] += anomaly<jupyter_output><empty_output><jupyter_text>Visualize one such anomaly...<jupyter_code>plt.plot(long_series[anom_indexes[0]: anom_indexes[0]+400])<jupyter_output><empty_output><jupyter_text>Rolling windows- choose windows of size 400- overlap by 200<jupyter_code>def window_stack(a, stepsize=200, width=400):
    return np.vstack( a[i:1+i-width or None:stepsize] for i in range(0,width)).transpose()

wins = window_stack(long_series)
wins.shape
plt.plot(wins[:4].transpose())
wins_scaled = mms.fit_transform(wins.transpose()).transpose()
plt.plot(wins_scaled[:4].transpose())
fff, Pxx = welch(wins_scaled, nperseg=64, axis=1)
X2 = pd.DataFrame(Pxx)
y_pred = model.predict(X2)
sum(y_pred)<jupyter_output><empty_output><jupyter_text>Did we find them?<jupyter_code>plt.plot(wins_scaled[y_pred == True][:4].transpose())<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.14647179003614647, 'fraction_numerical': 0.025616847399025616, 'mean_word_length': 6.685990338164252, 'pattern_counts': {'":': 0, '<': 74, '<?xml version=': 0, '>': 75, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '18', 'dir': 'jupyter-structured-clean-dedup', 'id': '266298', 'license': 'Apache-2.0', 'n_tokens_mistral': 2516, 'n_tokens_neox': 2352, 'n_words': 453, 'path': 'myapps/jupyter/TimeSeries/Time series classification.ipynb', 'repo_name': 'alonsoir/pipeline', 'size': '6374'}