2019년 12월 4일 수요일

LiPO 배터리 충전 에러 문제

LiPO 베터리 충전시 BATTERY INVALID ERROR 란 에러가 발생하는 경우가 있다. 베터리는 여러개의 셀로 연결되어 있은 데, 이 경우, 특정 셀이 충전되기 어려울 만큼 과전압이거나 저전압일 때 발생한다. 

다음은 이런 경우 해결할 수 있는 방법을 보여준다. 

2019년 12월 3일 화요일

NVIDIA JETSON NANO 배터리 전압 문제 솔류션

이 글은 NVIDIA JETSON NANO 배터리 파워 문제 솔류션과 관련된 내용이다.

엔비디아 젯슨 나노는 매우 저렴한 가격(10만원)으로 모바일 가능한 작은 크기의 딥러닝 솔류션을 지원한다. 전력 소모량이 적어 5V, 3-4A에서도 동작한다. 연결할 센서 등 장치가 많은 전류를 소모하면 4A이상 필요하다. 하지만, 나노는 정확한 전압 및 전류량에 동작하도록 설계되는 바람에 배터리 연결 시 제대로 부팅조차 되지 않는 등 문제가 있어, 많은 유저가 시행착오를 하고 있다(엔비디아 위키 페이지에도 정확한 솔류션은 나와있지 않다).

배터리를 통해 나노를 동작시키고, 나노에 전력소모량이 많은 센서나 장치를 연결하려면, 최소 5V, 4A가 필요하다. 이 전압과 전류를 처리하기 위해 다음과 같이 DC-DC 전압강하(STEP DOWN) 컨버터(CONVERTER)가 필요하다.
DC-DC XL4016 STEP DOWN CONVERTER

DC-DC 회로는 큰 입력 전압(40-10V)를 다른 전압으로 변환해 주는 회로이다. 이 회로는 POWER REGULATOR가 포함되어 있다. 전압은 맞춰줄 수 있지만, 나노에 필요한 전류는 나노에서 소모하는 양에 따라 달라지므로, 시행착오가 필요하다.

이 회로는 입력 전압, 출력 전압 핀이 두개 양쪽에 있다. 오른쪽 하단에 작은 출력전압 나사 조절기가 있고, 이 나사를 조정해 실제 전압을 맞춰줘야 한다. 전압은 눈으로 확인할 수 없으므로, 전압과 전류를 측정할 수 있는 디지털 미터기가 있으면 편리하다. 본인의 경우, 전압은 5.63V, 전류는 0.024A 수준으로 설정하였다.
전압 조정
전류 조정

이를 이용해, 5V보다 약간 높은 전압(5.5-5.7V)를 맞춘 후, 배터리를 입력단자에 연결하고, 출력단자를 나노 POWER 입력에 연결한다. 만약 과전압이나 저전압을 나노에 제공하면, 부팅이 안될 것이다. 다시 조정하는 과정을 거쳐 적절한 전압을 찾는다.
배터리, 전압 변환 회로 및 나노 연결

이렇게 회로 출력전압 나사를 조정해 나노의 동작전압을 찾아낸다. 다음은 이렇게 출력전압을 조정한 후 나노에 RGBD 카메라를 연결해 동작시킨 모습이다. 배터리만 연결한 상태에서 3시간 이상 나노와 센서가 동작되는 것을 확인하였다.
나노 및 센서 정상 동작 확인

나노는 가성비가 높은 임베디드 보드이다. 하지만, 이와 같이 까다로운 요구조건이 있어 처음 사용시 고생을 한다. 참고로, 다음에 이외 해결 방법에 대한 글을 링크한다.

2019년 12월 2일 월요일

아나콘다 완전 삭제 방법

아나콘다는 손쉬운 패키지 설치 및 설정으로 많이 활용되고 있다. 하지만, 버전이 엉키면, 패키지를 빌드하고 설치할 때 의존성 문제가 발생하여 불편한 경우가 발생한다.

이 경우, 아나콘다를 완전히 삭제하고 재설치해야 한다. 다음은 이를 위한 방법이다.

conda install anaconda-clean
anaconda-clean --yes
rm -rf ~/anaconda3             # removes the entire anaconda directory
rm -rf ~/.anaconda_backup

상세 내용은 여기를 참고한다. 

Keras 를 이용한 대중교통 사용 데이터 이상 탐지

이 글은 Keras 를 이용한 대중교통 사용 데이터 이상 탐지 방법을 간략히 소개한다. 좀 더 상세한 내용은 레퍼런스를 참고한다.

딥러닝 도구는 Google CoLAB을 사용하겠다. 데이터는 다음 링크에서 다운로드 한다.
이상 데이터 패턴은 다음과 같다. 

데이터를 다운로드 받고, CoLAB에 업로드 한다.
from google.colab import files
uploaded = files.upload()
for fn in uploaded.keys():
  print('User uploaded file "{name}" with length {length} bytes'.format(
      name=fn, length=len(uploaded[fn])))

소스 코드는 다음과 같다.

1. pands를 이용해 데이터 프레임을 생성한다.
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
import tqdm

from sklearn.preprocessing import StandardScaler, MinMaxScaler
from sklearn.metrics import mean_squared_log_error, mean_squared_error
import re
import datetime

from keras.models import *
from keras.layers import *
from keras.layers.core import Lambda
from keras import backend as K

### READ DATA AND CREATE FEATURES FOR year, month, day, hour ###
df = pd.read_csv('nyc_taxi.csv')

date = pd.DataFrame(np.vstack(df.timestamp.apply(lambda x: re.sub(r':|-' ,' ', x).split()))[:,:4].astype(int))

df = pd.concat([df, date], axis=1)
df.columns = ['timestamp','value','yr','mt','d','H']
df.head()
print(df.shape)
df.head()

### PLOT SAMPLE OF DATA ###
#df.iloc[4000:4000+7*48,:].plot(y='value',x='timestamp',figsize=(8,6))
df.iloc[8400:8400+7*48,:].plot(y='value',x='timestamp',figsize=(8,6))
plt.xticks(rotation=70)



2. LSTM을 위한 학습용 데이터를 생성한다.
### WEEKLY AUTOCORR PLOT (10 WEEKS DEPTH) ###
timeLags = np.arange(1,10*48*7)
autoCorr = [df.value.autocorr(lag=dt) for dt in timeLags]

plt.figure(figsize=(19,8))
plt.plot(1.0/(48*7)*timeLags, autoCorr);
plt.xlabel('time lag [weeks]'); plt.ylabel('correlation coeff', fontsize=12);

### CREATE WEEKDAY FEATURE AND COMPUTE THE MEAN FOR WEEKDAYS AT EVERY HOURS ###
weekday = df[['yr', 'mt', 'd']].apply(lambda x: datetime.datetime(x['yr'], x['mt'], x['d']).weekday(),axis=1).values
# weekend = weekday.copy()
# weekend[np.logical_and(weekend != 5, weekend != 6)] = 0
# weekend[weekend != 0] = 1

df['weekday'] = weekday
df['weekday_hour'] = df.weekday.astype(str) +' '+ df.H.astype(str)
df['m_weekday'] = df.weekday_hour.replace(df[:5000].groupby('weekday_hour')['value'].mean().to_dict())

### CREATE GENERATOR FOR LSTM ###
sequence_length = 48

def gen_index(id_df, seq_length, seq_cols):
    data_matrix =  id_df[seq_cols]
    num_elements = data_matrix.shape[0]
    for start, stop in zip(range(0, num_elements-seq_length, 1), range(seq_length, num_elements, 1)):
        yield data_matrix[stop-sequence_length:stop].values.reshape((-1,len(seq_cols)))

### CREATE AND STANDARDIZE DATA FOR LSTM ### 
cnt, mean = [], []
for sequence in gen_index(df, sequence_length, ['value']):
    cnt.append(sequence)

for sequence in gen_index(df, sequence_length, ['m_weekday']):
    mean.append(sequence)

cnt, mean = np.log(cnt), np.log(mean)
cnt = cnt - mean
cnt.shape

### CREATE AND STANDARDIZE LABEL FOR LSTM ###
init = df.m_weekday[sequence_length:].apply(np.log).values
label = df.value[sequence_length:].apply(np.log).values - init
label.shape

### DEFINE QUANTILE LOSS ###
def q_loss(q,y,f):
    e = (y-f)
    return K.mean(K.maximum(q*e, (q-1)*e), axis=-1)

### TRAIN TEST SPLIT ###
X_train, X_test = cnt[:5000], cnt[5000:]
y_train, y_test = label[:5000], label[5000:]

3. LSTM 모델을 생성한다. ### CREATE MODEL ###
losses = [lambda y,f: q_loss(0.1,y,f), lambda y,f: q_loss(0.5,y,f), lambda y,f: q_loss(0.9,y,f)]

inputs = Input(shape=(X_train.shape[1], X_train.shape[2]))
lstm = Bidirectional(LSTM(64, return_sequences=True, dropout=0.3))(inputs, training = True)
lstm = Bidirectional(LSTM(16, return_sequences=False, dropout=0.3))(lstm, training = True)
dense = Dense(50)(lstm)
out10 = Dense(1)(dense)
out50 = Dense(1)(dense)
out90 = Dense(1)(dense)

model = Model(inputs, [out10,out50,out90])
model.compile(loss=losses, optimizer='adam', loss_weights = [0.3,0.3,0.3])



4. 데이터를 학습한다.
history = model.fit(X_train, [y_train,y_train,y_train], epochs=50, batch_size=128, verbose=2, shuffle=True)


5. 각 주기별(10, 50, 90) 데이터 예측치 확인. keras backend function을 이용해 값 예측.
### QUANTILEs BOOTSTRAPPING ###
pred_10, pred_50, pred_90 = [], [], []
NN = K.function([model.layers[0].input, K.learning_phase()], 
                [model.layers[-3].output,model.layers[-2].output,model.layers[-1].output])

for i in tqdm.tqdm(range(0,100)):
    predd = NN([X_test, 0.5])
    pred_10.append(predd[0])
    pred_50.append(predd[1])
    pred_90.append(predd[2])

pred_10 = np.asarray(pred_10)[:,:,0] 
pred_50 = np.asarray(pred_50)[:,:,0]
pred_90 = np.asarray(pred_90)[:,:,0]

### REVERSE TRANSFORM PREDICTIONS ###
pred_90_m = np.exp(np.quantile(pred_90,0.9,axis=0) + init[5000:])
pred_50_m = np.exp(pred_50.mean(axis=0) + init[5000:])
pred_10_m = np.exp(np.quantile(pred_10,0.1,axis=0) + init[5000:])

6. 이상 데이터 확인 ### EVALUATION METRIC ###
mean_squared_log_error(np.exp(y_test + init[5000:]), pred_50_m)

### PLOT QUANTILE PREDICTIONS ###
plt.figure(figsize=(16,8))
plt.plot(pred_90_m, color='cyan')
plt.plot(pred_50_m, color='blue')
plt.plot(pred_10_m, color='green')

### CROSSOVER CHECK ###
plt.scatter(np.where(np.logical_or(pred_50_m>pred_90_m, pred_50_m<pred_10_m))[0], 
            pred_50_m[np.logical_or(pred_50_m>pred_90_m, pred_50_m<pred_10_m)], c='red', s=50)

### PLOT UNCERTAINTY INTERVAL LENGHT WITH REAL DATA ###
plt.figure(figsize=(16,8))
plt.plot(np.exp(y_test + init[5000:]), color='red', alpha=0.4)
plt.scatter(range(len(pred_10_m)), pred_90_m - pred_10_m)

그럼 다음과 같은 결과를 얻을 수 있다. 아래에서 청색 지점은 이상이 상대적으로 많은 곳이다.
레퍼런스

LSTM을 이용한 주가 이상 탐지

이 글은 LSTM을 이용한 주가 이상 탐지 방법을 간단히 요약한다.

예측 데이터는 The S&P 500 지수이다. 이 데이터는 https://www.kaggle.com/pdquant/sp500-daily-19862018 에서 다운로드 받는다.

데이터를 우선 다음과 같이 다운로드한다.
!gdown --id 10vdMg_RazoIatwrT7azKFX4P02OebU76 --output spx.csv

소스코드는 다음과 같다.

df = pd.read_csv('spx.csv', parse_dates=['date'], index_col='date')
train_size = int(len(df) * 0.95)
test_size = len(df) - train_size
train, test = df.iloc[0:train_size], df.iloc[train_size:len(df)]
print(train.shape, test.shape)

from sklearn.preprocessing import StandardScaler
scaler = StandardScaler()
scaler = scaler.fit(train[['close']])

train['close'] = scaler.transform(train[['close']])
test['close'] = scaler.transform(test[['close']])

def create_dataset(X, y, time_steps=1):
    Xs, ys = [], []
    for i in range(len(X) - time_steps):
        v = X.iloc[i:(i + time_steps)].values
        Xs.append(v)
        ys.append(y.iloc[i + time_steps])
    return np.array(Xs), np.array(ys)

TIME_STEPS = 30

# reshape to [samples, time_steps, n_features]
X_train, y_train = create_dataset(
  train[['close']],
  train.close,
  TIME_STEPS
)

X_test, y_test = create_dataset(
  test[['close']],
  test.close,
  TIME_STEPS
)

print(X_train.shape)

model = keras.Sequential()
model.add(keras.layers.LSTM(
    units=64,
    input_shape=(X_train.shape[1], X_train.shape[2])
))

model.add(keras.layers.Dropout(rate=0.2))
model.add(keras.layers.RepeatVector(n=X_train.shape[1]))
model.add(keras.layers.LSTM(units=64, return_sequences=True))
model.add(keras.layers.Dropout(rate=0.2))
model.add(
  keras.layers.TimeDistributed(
    keras.layers.Dense(units=X_train.shape[2])
  )
)

model.compile(loss='mae', optimizer='adam')

history = model.fit(
    X_train, y_train,
    epochs=10,
    batch_size=32,
    validation_split=0.1,
    shuffle=False
)

X_train_pred = model.predict(X_train)
train_mae_loss = np.mean(np.abs(X_train_pred - X_train), axis=1)

THRESHOLD = 0.65
X_test_pred = model.predict(X_test)
test_mae_loss = np.mean(np.abs(X_test_pred - X_test), axis=1)

test_score_df = pd.DataFrame(index=test[TIME_STEPS:].index)
test_score_df['loss'] = test_mae_loss
test_score_df['threshold'] = THRESHOLD
test_score_df['anomaly'] = test_score_df.loss > test_score_df.threshold
test_score_df['close'] = test[TIME_STEPS:].close

anomalies = test_score_df[test_score_df.anomaly == True]

다음은 주식 예측에서 이상탐지된 결과이다.

레퍼런스

2019년 11월 30일 토요일

어린이도 손쉽게 만드는 인공지능 딥러닝 도구 - 구글 Teachable machine

더이상 인공지능 딥러닝은 사용하기 어려운 기술이 아니다. 좀 더 정확히 말하면, 단순한 데이터 입력-학습-사용 정도의 딥러닝은 기술이 아닌 어린이도 할 수 있는 세상이 되어가고 있다.

구글에서 만든 teachable machine은 누구나 딥러닝을 쉽게 할 수 있도록 한다.
teachable machine site

사용 방법은 매우 간단하다.

1. 해당 사이트를 방문하고, 이미지 인식, 시그널 분류 등 준비된 템플릿을 선택한다.
2. 이제 데이터를 입력하고, 학습하고, 실행한 결과를 본다. 데이터는 웹캠으로 입력 가능하다.

본인은 건설 공학에서 고전적인 문제인 균열 탐지 부분을 간단히 데이터 만들어 넣어 보았다.
딥러닝 학습 중인 모습

학습 후 균열 검출

정확도가 상당히 잘 나온다. 보통, 40~50개 정도 이미지를 준비하면 분류 모델 만드는 데 큰 문제가 없었다. 클래스를 추가하는 것도 간단하다. 다음과 같이 smile 이미지를 준비해 학습했다.
Smile 이미지 학습
인식 결과

학습 후 인식하면 인식이 잘 되는 것을 알 수 있다.

이렇게 생성된 딥러닝 모델은 keras, tensorflow js 등의 모델로 다운로드 받아 다양한 기기에서 딥러닝을 사용할 수 있다. 애써 만든 학습용 데이터는 다운로드해 다시 재활용할 수 있다.
Export model

저장된 모델을 이용해 다양한 방식으로 배포할 수 있다.

이를 이용하면, 스마트폰 등 다양한 기기에서 딥러닝 실행 결과를 손쉽게 얻을 수 있다.
ㅎ 이제 학습용 데이터만 준비 잘 하면 되겠다.^^
이미지 딥러닝 분류 영상

조만간 이를 이용한 다양한 사례가 쏟아져 나올 듯 보인다. 템플릿이 겨우 3개이며, 복잡한 딥러닝을 수행할 수는 없다. 하지만, 구글은 템플릿을 계속 개발하고 있으며, 다양한 딥러닝 모델 학습 기능을 추가할 계획이다.

딥러닝은 이제 누가 손쉽게 모델을 빨리 만들 수 있고, 사용할 수 있도록 하는 지가 기술이 될 듯 하다.

2019년 11월 28일 목요일

LSTM 딥러닝 모델 기반 이상 데이터 탐지

이 글은 LSTM 기반 이상신호 탐지용 딥러닝 모델 개발 방법을 간략히 설명한다. 이상탐지는 장비 모니터링, 사기 탐지, 해킹 등 다양한 분야에서 응용될 수 있다. 이 글은 LSTM(Long Short Term Memory) 딥러닝 모델을 이용한 이상탐지 방법을 간략히 설명하고 구현하는 방법을 설명한다.

머리말
시계열에서 이상탐지는 예상되는 미래값이 실제 얻은 값과 편차를 가지게 될 때 이상으로 판단한다. 학습 단계는 다음과 같다.

1. 시계열 데이터의 과거값들을 이용해, 한단계 바로 뒤의 값들을 예측하도록 학습 데이터를 만든다.
2. 다변량 가우스 분포를 이용해 오차 벡터를 계산한다.
3. 만약, 예측된 값과 실제값의 오차가 다변량 가우스 분포의 끝 부분에 위치해 있다면, 일반적으로 예측될 수 있는 값의 범위를 넘어가므로, 이상이라 판단한다.

다음 그림은 이를 보여준다.

Mahalanobis 거리
분포내 예측과 오차에 대한 희소한 양을 표시하기 위해 마할라노비스 거리란 함수를 사용한다. 수식은 다음과 같다.

실제 코딩해 보기
간단한 예제를 코딩해 보도록 한다. 아래와 같은 패키지가 우선 설치되어 있어야 한다.
  • ReNom 2.6.2
  • matplotlib 2.2.2
  • pandas 0.23.1
  • numpy 1.14.5
  • scikit-learn 0.19.1
  • scipy 1.1.0
제대로 설치되어 있다면, 다음 코드를 실행해 본다.
%matplotlib inline
import matplotlib.pyplot as plt
import pandas as pd 
import numpy as np 
from sklearn.model_selection import train_test_split
from copy import deepcopy
from sklearn.preprocessing import StandardScaler

import renom as rm
from renom.optimizer import Adam
from renom.cuda import set_cuda_active
set_cuda_active(False)

데이터셋은 다은 링크에서 ECG dataset, qtdb/sel102 ECG dataset 을 사용한다.
http://www.cs.ucr.edu/~eamonn/discords/

ECG data를 정규화하기 위해, 시간을 0 - 5000 부분만 플롯팅해 그 구조를 확인해 본다.
df = pd.read_csv('data/qtdbsel102.txt', header=None, delimiter='\t')
ecg = df.iloc[:,2].values
ecg = ecg.reshape(len(ecg), -1)
print('length of ECG data : ', len(ecg))

# standardize
scaler = StandardScaler()
std_ecg = scaler.fit_transform(ecg)

plt.style.use('ggplot')
plt.figure(figsize=(15,5))
plt.xlabel('time')
plt.ylabel('ECG\'s value')
plt.plot(np.arange(5000), std_ecg[:5000], color='b')
plt.ylim(-3, 3)
x = np.arange(4200,4400)
y1 = [-3]*len(x)
y2 = [3]*len(x)
plt.fill_between(x, y1, y2, facecolor='g', alpha=.3)
plt.show()

이 데이터는 주기적이다. 대략 4250 시간 근처에 주기가 붕괴되어 변화가 있다.

1단계 - LSTM 학습
1단계의 목적은 LSTM으로 정상 데이터를 학습하는 것이다. 그러므로, 5000이후 데이터를 학습용 데이터로 사용한다.
normal_cycle = std_ecg[5000:]

plt.figure(figsize=(10,5))
plt.title("training data")
plt.xlabel('time')
plt.ylabel('ECG\'s value')
plt.plot(np.arange(5000,8000), normal_cycle[:3000], color='b')# stop plot at 8000 times for friendly visual
plt.show()


이제, 다음 코드로 d 길이만큼 서브시퀀스를 생성하고, l 차원의 라벨을 생성한다.
 create data of the "look_back" length from time-series, "ts"
# and the next "pred_length" values as labels
def create_subseq(ts, look_back, pred_length):
    sub_seq, next_values = [], []
    for i in range(len(ts)-look_back-pred_length):  
        sub_seq.append(ts[i:i+look_back])
        next_values.append(ts[i+look_back:i+look_back+pred_length].T[0])
    return sub_seq, next_values

이 예제에서는 d = 10, l = 3이다. 아래 코드를 실행해, sub_seq와 예측해야할 라벨이 있는 next_values 시퀀스를 만든 후, 훈련 및 검증 데이터로 분리한다.
look_back = 10
pred_length = 3

sub_seq, next_values = create_subseq(normal_cycle, look_back, pred_length)


X_train, X_test, y_train, y_test = train_test_split(
    sub_seq, next_values, test_size=0.2)

X_train = np.array(X_train)
X_test = np.array(X_test)
y_train = np.array(y_train)
y_test = np.array(y_test)


train_size = X_train.shape[0]
test_size = X_test.shape[0]
print('train size:{}, test size:{}'.format(train_size, test_size))

LSTM 모델을 정의한다.
# model definition
model = rm.Sequential([
    rm.Lstm(35),
    rm.Relu(),
    rm.Lstm(35),
    rm.Relu(),
    rm.Dense(pred_length)
    ])

LSTM에서 35개 벡터를 출력하고, RELU 함수를 거친 후, 다시 LSTM 35개 벡터 출력, RELU 함수 처리, 3개 벡터값을 예측하기 위해 Dense 레이어를 추가했다.

배치 크기와 최대 세대를 설정하고, 최적화 함수는 ADAM으로 한다.
# params
batch_size = 100
max_epoch = 2000
period = 10 # early stopping checking period
optimizer = Adam()

이제, 다음과 같이 학습시킨다. 참고로, 아래 코드에서 np.random.permutation()은 학습을 잘 시켜주기 위해, 순열을 섞어주는 함수이다(링크). 큰 순서는 다음과 같다.
1. model.train()에서 학습한다.
2. 모든 데이터에 대해, 바로 학습된 모델에 학습데이터를 입력하여, 출력된 z와 실제 라벨값인 batch_y와 차이를 MSE(Mean Square Error)로 계산하여, train_loss 값을 누적한다.
3. 모든 데이터에 대해, 학습된 모델을 이용해 예측 에러가 얼마인지 계산한다.
앞의 과정을 테스트 손실이 특정 조건이 될때까지 반복해 학습시킨다.
# Train Loop
epoch = 0
loss_prev = np.inf

learning_curve, test_curve = [], []

while(epoch < max_epoch):
    epoch += 1

    perm = np.random.permutation(train_size)
    train_loss = 0

    for i in range(train_size // batch_size):
        batch_x = X_train[perm[i*batch_size:(i+1)*batch_size]]
        batch_y = y_train[perm[i*batch_size:(i+1)*batch_size]]

        # Forward propagation
        l = 0
        z = 0
        with model.train():
            for t in range(look_back):
                z = model(batch_x[:,t])
                l = rm.mse(z, batch_y)
            model.truncate()
        l.grad().update(optimizer)
        train_loss += l.as_ndarray()

    train_loss /= (train_size // batch_size)
    learning_curve.append(train_loss)

    # test
    l = 0
    z = 0
    for t in range(look_back):
        z = model(X_test[:,t])
        l = rm.mse(z, y_test)
    model.truncate()
    test_loss = l.as_ndarray()
    test_curve.append(test_loss)

    # check early stopping
    if epoch % period == 0:
        print('epoch:{} train loss:{} test loss:{}'.format(epoch, train_loss, test_loss))
        if test_loss > loss_prev*0.99:
            print('Stop learning')
            break
        else:
            loss_prev = deepcopy(test_loss)

plt.figure(figsize=(10,5))
plt.plot(learning_curve, color='b', label='learning curve')
plt.plot(test_curve, color='orange', label='test curve')
plt.xlabel('epoch')
plt.ylabel('loss')
plt.legend(fontsize=20)
plt.show()

실행 결과는 다음과 같다.

2단계 - 가우스 분포 피팅
이상 탐지를 위해 가우스 분포 계산을 한다.
# computing errors
for t in range(look_back):
    pred = model(X_test[:,t])
model.truncate()
errors = y_test - pred

mean = sum(errors)/len(errors)

cov = 0
for e in errors:
    cov += np.dot((e-mean).reshape(len(e), 1), (e-mean).reshape(1, len(e)))
cov /= len(errors)

print('mean : ', mean)
print('cov : ', cov)

mean :  [-0.00471252  0.00561184  0.01125641]
cov :  [[0.00093565 0.00088413 0.00097755]
 [0.00088413 0.00208558 0.0025572 ]
 [0.00097755 0.0025572  0.00498106]]

3단계 - 이상탐지
이 예측 모델이 미지의 데이터에 대해서도 제대로 동작하는 지 검증한다. 우선, 서브 시퀀스와 라벨을 생성하고 에러 벡터를 계산한다.
# calculate Mahalanobis distance
def Mahala_distantce(x,mean,cov):
    d = np.dot(x-mean,np.linalg.inv(cov))
    d = np.dot(d, (x-mean).T)
    return d

# anomaly detection
sub_seq, next_values = create_subseq(std_ecg[:5000], look_back, pred_length)
sub_seq = np.array(sub_seq)
next_values = np.array(next_values)

for t in range(look_back):
    pred = model(sub_seq[:,t])
model.truncate()
errors = next_values - pred

다음으로 Mahalanobis 거리를 출력해 얼마나 이상한 값인지를 판단한다.
m_dist = [0]*look_back 
for e in errors:
    m_dist.append(Mahala_distantce(e,mean,cov))

fig, axes = plt.subplots(nrows=2, figsize=(15,10))

axes[0].plot(std_ecg[:5000],color='b',label='original data')
axes[0].set_xlabel('time')
axes[0].set_ylabel('ECG\'s value' )
axes[0].set_ylim(-3, 3)
x = np.arange(4200,4400)
y1 = [-3]*len(x)
y2 = [3]*len(x)
axes[0].fill_between(x, y1, y2, facecolor='g', alpha=.3)

axes[1].plot(m_dist, color='r',label='Mahalanobis Distance')
axes[1].set_xlabel('time')
axes[1].set_ylabel('Mahalanobis Distance')
axes[1].set_ylim(0, 1000)
y1 = [0]*len(x)
y2 = [1000]*len(x)
axes[1].fill_between(x, y1, y2, facecolor='g', alpha=.3)

plt.legend(fontsize=15)
plt.show()


앞의 그림을 보면 4250 근처에 Mahalanobis 거리가 가장 크다. 이를 통해, 주기적 데이터가 언제 붕괴가 일어나는 지를 확인할 수 있다.

이런 방식으로 어떤 데이터든지 Anomaly Detection이 가능하다.


레퍼런스
  • ReNom, 2018, LSTM based anomaly detection
  • Malhotra, Pankaj, et al. "Long short term memory networks for anomaly detection in time series." Proceedings. Presses universitaires de Louvain, 2015.
  • E. Keogh, J. Lin and A. Fu (2005). HOT SAX: Efficiently Finding the Most Unusual Time Series Subsequence. In The Fifth IEEE International Conference on Data Mining.