이상치 탐색

이상치 탐색¶

---

주가 데이터 분석¶

• pip install requests
• mkdir -p ~/aiffel/anomaly_detection/kospi

In [4]:

import requests
import os

# 아래 url은 yahoo finance 서버에 우리나라 코스피 데이터를 요청하는 주소입니다.
url = "https://query1.finance.yahoo.com/v7/finance/download/%5EKS11?period1=867715200&period2=1597276800&interval=1d&events=history"

# 데이터를 요청하고 그 결과를 response로 받습니다.
response = requests.get(url)

csv_file = os.getenv('HOME')+'/aiffel/anomaly_detection/kospi/kospi.csv'

# response 의 컨텐츠를 csv로 저장합니다.
with open(csv_file, "w") as fp:  
    fp.write(response.text)

In [5]:

# 저장한 csv를 읽어옵니다.
import pandas as pd
df = pd.read_csv(csv_file)
df.head(2)

Out[5]:

	Date	Open	High	Low	Close	Adj Close	Volume
0	1997-07-01	744.979980	760.520020	744.669983	758.030029	758.030029	37000.0
1	1997-07-02	763.390015	777.289978	763.390015	777.289978	777.289978	52200.0

In [6]:

# 날짜데이터를 Datetime 형식으로 바꿔줍니다.
df.loc[:,'Date'] = pd.to_datetime(df.Date)

# 데이터의 정합성을 확인합니다
df.isna().sum()

Out[6]:

Date           0
Open         142
High         142
Low          142
Close        142
Adj Close    142
Volume       142
dtype: int64

In [7]:

print("삭제 전 데이터 길이(일자수):",len(df))
df = df.dropna(axis=0).reset_index(drop=True)

print("삭제 후 데이터 길이(일자수):",len(df))
df.isna().sum()

삭제 전 데이터 길이(일자수): 5842
삭제 후 데이터 길이(일자수): 5700

Out[7]:

Date         0
Open         0
High         0
Low          0
Close        0
Adj Close    0
Volume       0
dtype: int64

In [8]:

import matplotlib.pyplot as plt
from matplotlib.pylab import rcParams

plt.rcParams["figure.figsize"] = (10,5)
# Line Graph by matplotlib with wide-form DataFrame

plt.plot(df.Date, df.Close, marker='s', color='r')
plt.plot(df.Date, df.High, marker='o', color='g')
plt.plot(df.Date, df.Low, marker='*', color='b')
plt.plot(df.Date, df.Open, marker='+', color='y')

plt.title('KOSPI ', fontsize=20) 
plt.ylabel('Stock', fontsize=14)
plt.xlabel('Date', fontsize=14)
plt.legend(['Close', 'High', 'Low', 'Open'], fontsize=12, loc='best')

plt.show()

In [9]:

df.loc[df.Low > df.High]

Out[9]:

	Date	Open	High	Low	Close	Adj Close	Volume

1. Outlier / Novelties 구분하기¶

In [10]:

# Outlier 제거
# 카카오 주식차트 결과로 대체합니다. 
df.loc[df.Date == '2020-05-06', 'Low'] = 1903 

# 비정상데이터가 제거되었는지 다시 확인해 봅니다.
df.loc[df.Low>df.High]

Out[10]:

	Date	Open	High	Low	Close	Adj Close	Volume

중간 정리¶

• 이상치(Anomailies) = 극단치(Outlier) + 특이치(Novelites)
• 극단치는 제거해야 모형에 좋습니다.
• 특이치는 남겨둬야 모형에 좋습니다.

• pip install statsmodels (통계 모델 사용 라이브러리)
• 이상치를 찾는 방법
  • 정규분포를 따르는지 확인하는 방법 : z-test
  • 시계열 데이터 중에서 정규분포에 가까운 데이터를 뽑아내는 방법 : Time series decomposition

Z-Test¶

In [11]:

%matplotlib inline
from scipy.stats import norm
import numpy as np
import matplotlib.pyplot as plt
x = np.arange(-5, 5, 0.001)
y = norm.pdf(x,0,1) 
# 평균이 0이고, 표준편차가 1인 정규분포를 그립니다.

# build the plot
fig, ax = plt.subplots(figsize=(9,6))
ax.fill_between(x,y,0, alpha=0.3, color='b')
ax.set_xlim([-4,4])
ax.set_title('normal distribution')
plt.show()

In [12]:

# 주식 종가 plot
fig, ax = plt.subplots(figsize=(9,6))
_ = plt.hist(df.Close, 100, density=True, alpha=0.75)

In [13]:

from statsmodels.stats.weightstats import ztest
_, p = ztest(df.Close)
print(p)

0.0

z-test를 통한 정규분포 신뢰도 측정¶

• p가 0.05 이하로 나왔으면 normal distribution과 거리가 멀다는 뜻으로 이런 데이터로는 '정규분포'를 가정한 통계적 추정이 어렵다. 그래서 이 데이터에 대해서는 정규분포를 가정한 신뢰구간 분석은 적용하기 어려울 것입니다.

Time series decomposition¶

시계열 분해에 대한 설명

In [14]:

from statsmodels.tsa.seasonal import seasonal_decompose
result = seasonal_decompose(df.Close, model='additive', two_sided=True, 
                            period=50, extrapolate_trend='freq') # 계절적 성분 50일로 가정
result.plot()
plt.show()

위 함수에 대한 설명¶

• two_sided=True란?
  Trend 성분을 만들기 위해 freq 길이에 해당하는 샘플이 필요합니다. 이걸 rolling window로 만들어냅니다.
  우리의 예제에서는 period=50 이기 때문에
  two_sided=True라면 049번째 데이터로 25번째 Trend값을, 150번째 데이터로 26번째 Trend값을 생성합니다
  two_sided=False라면 049번째 데이터로 50번째 Trend값을, 150번째 데이터로 51번째 Trend값을 생성합니다.

• extrapolate_trend='freq'란?
  위의 Trend 성분을 만들기 위한 rolling window 때문에 필연적으로 trend, resid에는 Nan 값이 발생합니다
  우리의 예제에서는 period=50 이기 때문에
  two_sided=True라면 맨앞에서 25개, 맨뒤에서 25개가 Nan이고, two_sided=False라면 맨앞에서 50개가 Nan입니다.
  이렇게 발생하는 Nan 값을 채워주는 옵션이 extraplate_trend 입니다
  이 옵션을 빼면 Trend와 Resid에 Nan 값들이 생깁니다.

• additive 는 Observed[t] = trend[t] + seasonal[t] + resid[t]
  multicative 는 Observed[t] = trend[t] x seasonal[t] x resid[t]
  로그를 취하게 되면 동일해진다.

In [15]:

#그래프가 너무 작아서 안보인다면

fig, axes = plt.subplots(ncols=1, nrows=4, sharex=True, figsize=(12,9))
result.observed.plot(ax=axes[0], legend=False)
axes[0].set_ylabel('Observed')
result.trend.plot(ax=axes[1], legend=False)
axes[1].set_ylabel('Trend')
result.seasonal.plot(ax=axes[2], legend=False)
axes[2].set_ylabel('Seasonal')
result.resid.plot(ax=axes[3], legend=False)
axes[3].set_ylabel('Residual')
plt.show()

In [16]:

# seasonal 성분은 너무 빼곡하게 보이네요. 다시 확인해보겠습니다.
result.seasonal[:100].plot()
#-8 에서 10 사이를 주기적으로 반복하는게 보이네요.

Out[16]:

<AxesSubplot:>

In [17]:

fig, ax = plt.subplots(figsize=(9,6))
_ = plt.hist(result.resid, 100, density=True, alpha=0.75)

In [18]:

# p가 0.5이상인것을 보아 데이터 분포가 정규분포를 잘 따르고 있다
r = result.resid.values
st, p = ztest(r)
print(st,p)

-0.670120907482477 0.5027807179026166

In [19]:

#3σ  기준 신뢰구간으로 이상치 찾기
# 평균과 표준편차 출력
mu, std = result.resid.mean(), result.resid.std()
print("평균:", mu, "표준편차:", std)

# 3-sigma(표준편차)를 기준으로 이상치 판단
print("이상치 갯수:", len(result.resid[(result.resid>mu+3*std)|(result.resid<mu-3*std)]))

평균: -0.35098081908284745 표준편차: 39.54282047378896
이상치 갯수: 73

In [20]:

df.Date[result.resid[
    (result.resid>mu+3*std)|(result.resid<mu-3*std)].index]

Out[20]:

2476   2007-07-20
2477   2007-07-23
2478   2007-07-24
2479   2007-07-25
2494   2007-08-16
          ...    
5607   2020-04-01
5650   2020-06-05
5651   2020-06-08
5652   2020-06-09
5653   2020-06-10
Name: Date, Length: 73, dtype: datetime64[ns]

신뢰구간 방법의 한계를 극복하기 위한 방법 Multi-variable Anomaly Detection¶

1. Clustering : 클러스터링으로 묶으면 정상인 데이터끼리 이상한 애들끼리 그룹핑되니 이상한 그룹을 찾는다.
2. Forecasting : 시계열 예측모델을 만들어서, 예측 오차가 크게 발생하는 지점은 이상한 상황이다. 일반적으로 Auto-Encoder로 탐색한다.

클러스터링으로 이상치 찾기 : K-means와 DBSCAN¶

In [21]:

# 데이터 전처리
def my_decompose(df, features, freq=50):
    trend = pd.DataFrame()
    seasonal = pd.DataFrame()
    resid = pd.DataFrame()
    
    # 사용할 feature 마다 decompose 를 수행합니다.
    for f in features:
        result = seasonal_decompose(df[f], 
                                    model='additive', period=freq, extrapolate_trend=freq)
        trend[f] = result.trend.values
        seasonal[f] = result.seasonal.values
        resid[f] = result.resid.values
        
    return trend, seasonal, resid

# 각 변수별 트렌드/계절적/잔차
tdf, sdf, rdf = my_decompose(df, features=['Open','High','Low','Close','Volume']) # 5가지 데이터를 모두 사용해서 찾는다
tdf.describe()

Out[21]:

	Open	High	Low	Close	Volume
count	5700.000000	5700.000000	5700.000000	5700.000000	5.700000e+03
mean	1463.676902	1472.666844	1452.416198	1462.928822	4.019559e+05
std	624.889151	625.455025	623.236019	624.480962	1.756523e+05
min	313.251601	317.456402	309.410898	313.594600	2.547328e+04
25%	829.360720	836.791400	819.983175	828.532226	3.084350e+05
50%	1650.650549	1659.088303	1638.829350	1650.863803	3.788582e+05
75%	1994.407723	2000.954524	1983.616602	1992.979502	4.754910e+05
max	2507.482712	2517.003208	2493.430989	2503.484602	1.200116e+06

In [22]:

rdf.describe()

Out[22]:

	Open	High	Low	Close	Volume
count	5700.000000	5700.000000	5700.000000	5700.000000	5.700000e+03
mean	-0.393551	-0.372158	-0.389955	-0.365469	8.841281e+01
std	39.263848	37.736831	40.474172	39.548444	1.039947e+05
min	-388.426223	-370.136754	-409.158640	-415.479447	-7.010681e+05
25%	-21.465941	-21.051663	-21.515153	-21.481658	-5.030125e+04
50%	1.330474	1.174430	1.742668	1.196649	-7.988840e+03
75%	22.539148	22.049242	22.731070	22.612395	3.829684e+04
max	165.608215	148.808476	168.754715	165.264163	1.198752e+06

In [23]:

# Volume이 숫자가 커서 한쪽으로 쏠리는 현상을 방지하기 위해 정규화
# 표준정규화
from sklearn.preprocessing import StandardScaler

scaler = StandardScaler()
scaler.fit(rdf)
print(scaler.mean_)
norm_rdf = scaler.transform(rdf)
norm_rdf

[-0.39355059 -0.3721578  -0.38995467 -0.36546896 88.41280953]

Out[23]:

array([[-1.99084772, -1.6946113 , -1.94666555, -1.67001221,  0.15229384],
       [-1.53038288, -1.27085417, -1.46921446, -1.20937409,  0.40060153],
       [-1.08911619, -1.10066391, -1.17272681, -1.19383211,  0.21720631],
       ...,
       [ 1.42518133,  2.04460423,  1.81422569,  2.12014492,  0.59771771],
       [ 1.77827091,  2.0759502 ,  1.67321373,  2.43501358, -0.12447359],
       [ 2.72673619,  2.59982453,  2.03457503,  2.37882082,  1.2708417 ]])

K-means 와 DBSCAN¶

• k-means 와 DBSCAN 은 대표적인 unsupervised clustering 알고리즘으로 입력된 데이터들을 유사한 몇개의 그룹으로 분류해준다.
• k-means 은 몇 개의 그룹으로 묶는지 미리 지정해 주는 반면, DBSCAN은 지정해 줄 필요가 없다.

K-means¶

In [24]:

from sklearn.cluster import KMeans

# random_state는 K-means가 난수를 사용하기때문에 동일한 분석결과를 얻기 위해서는 난수설정값을 통일시켜주는 과정이 필요해서
kmeans = KMeans(n_clusters=2, random_state=0).fit(norm_rdf)
print(kmeans.labels_) # 분류된 라벨은 이렇게 kemans.labels_ 로 확인합니다.

[1 1 1 ... 0 0 0]

In [25]:

# 라벨은 몇번 그룹인지 뜻합니다. 
# return_counts=True 를 해서 몇개의 샘플이 몇번 그룹에 할당되었는지 확인해요
lbl, cnt = np.unique(kmeans.labels_,return_counts=True) 
print(lbl) # 0번 그룹, 1번 그룹으로 나뉘어졌네요
print(cnt)

[0 1]
[3183 2517]

In [26]:

from sklearn.cluster import KMeans

kmeans = KMeans(n_clusters=15, random_state=0).fit(norm_rdf)
lbl, cnt = np.unique(kmeans.labels_,return_counts=True,)
['group:{}-count:{}'.format(group, count) for group, count in zip(lbl, cnt)]

Out[26]:

['group:0-count:99',
 'group:1-count:244',
 'group:2-count:430',
 'group:3-count:889',
 'group:4-count:104',
 'group:5-count:296',
 'group:6-count:687',
 'group:7-count:716',
 'group:8-count:23',
 'group:9-count:161',
 'group:10-count:336',
 'group:11-count:9',
 'group:12-count:968',
 'group:13-count:676',
 'group:14-count:62']

In [27]:

# 어떤 날들이 분류된 건지 탐색(원본 코스피 데이터에서 탐색) 
df[(kmeans.labels_==8)|(kmeans.labels_==11)|(kmeans.labels_==14)]

Out[27]:

	Date	Open	High	Low	Close	Adj Close	Volume
494	1999-07-08	1010.770020	1022.880005	997.989990	999.109985	999.109985	442000.0
495	1999-07-09	1004.679993	1027.930054	989.479980	1027.930054	1027.930054	448700.0
496	1999-07-12	1037.739990	1052.599976	1001.409973	1004.250000	1004.250000	489200.0
501	1999-07-19	1034.410034	1044.280029	1001.789978	1024.579956	1024.579956	379300.0
1233	2002-07-12	780.489990	793.440002	776.450012	792.929993	792.929993	1805000.0
...	...	...	...	...	...	...	...
5651	2020-06-08	2215.600098	2217.209961	2174.959961	2184.290039	2184.290039	854000.0
5652	2020-06-09	2206.229980	2212.169922	2166.010010	2188.919922	2188.919922	817100.0
5653	2020-06-10	2187.909912	2200.679932	2178.120117	2195.689941	2195.689941	706300.0
5654	2020-06-11	2184.360107	2200.719971	2148.510010	2176.780029	2176.780029	827400.0
5699	2020-08-13	2455.280029	2458.169922	2412.489990	2437.530029	2437.530029	907200.0

94 rows × 7 columns

In [28]:

df.describe()

Out[28]:

	Open	High	Low	Close	Adj Close	Volume
count	5700.000000	5700.000000	5700.000000	5700.000000	5700.000000	5.700000e+03
mean	1463.283352	1472.294686	1452.026243	1462.563353	1462.563353	4.020443e+05
std	627.378909	627.827355	625.856088	627.015237	627.015237	2.089416e+05
min	283.410004	291.010010	277.369995	280.000000	280.000000	3.016000e+03
25%	828.497513	837.232498	815.889999	826.982514	826.982514	2.839000e+05
50%	1653.059998	1659.780029	1635.830017	1647.359985	1647.359985	3.677000e+05
75%	1999.852478	2006.935028	1988.872498	1997.312500	1997.312500	4.809250e+05
max	2590.409912	2607.100098	2587.550049	2598.189941	2598.189941	2.379300e+06

In [29]:

# 2004-04-14 주변 정황 ??
df.iloc[1660:1670]

Out[29]:

	Date	Open	High	Low	Close	Adj Close	Volume
1660	2004-04-08	913.520020	916.859985	908.070007	916.859985	916.859985	431800.0
1661	2004-04-09	910.479980	910.479980	901.239990	905.440002	905.440002	352100.0
1662	2004-04-12	905.440002	926.070007	905.440002	918.859985	918.859985	366200.0
1663	2004-04-13	925.190002	927.669983	916.200012	917.630005	917.630005	468900.0
1664	2004-04-14	909.950012	919.380005	907.359985	916.309998	916.309998	407200.0
1665	2004-04-16	903.840027	908.719971	897.729980	898.880005	898.880005	391700.0
1666	2004-04-19	900.479980	902.099976	891.270020	902.099976	902.099976	476100.0
1667	2004-04-20	909.909973	919.130005	902.469971	918.900024	918.900024	456800.0
1668	2004-04-21	910.679993	931.210022	910.679993	929.950012	929.950012	554000.0
1669	2004-04-22	931.090027	933.840027	921.969971	924.010010	924.010010	519700.0

In [30]:

# 2. 각 그룹은 어떤 특징을 갖고 있는지
# 각 그룹의 중심부는 어떤 값을 가지고 있는지 확인
pd.DataFrame(kmeans.cluster_centers_, columns=['Open','High','Low','Close','Volume'])

Out[30]:

	Open	High	Low	Close	Volume
0	-0.365323	-0.393770	-0.329102	-0.369401	-2.629301
1	0.926568	0.966440	0.888317	0.908881	1.602310
2	-0.996339	-0.994707	-1.009813	-0.990381	0.333163
3	0.388424	0.401057	0.387116	0.383502	0.275210
4	-2.842317	-2.804001	-3.013889	-2.871066	0.162692
5	1.618193	1.677801	1.618695	1.622420	0.274198
6	0.194224	0.180621	0.210986	0.199175	-0.695088
7	-0.616413	-0.645679	-0.583279	-0.598685	-0.635260
8	0.420870	0.488848	0.535200	0.563688	6.947442
9	-0.139704	-0.123442	-0.205953	-0.172649	2.250147
10	-1.557925	-1.612050	-1.535077	-1.557732	-0.557971
11	-6.991717	-6.698342	-7.489608	-7.129549	-0.470756
12	-0.194926	-0.200217	-0.189814	-0.192799	0.127429
13	0.920997	0.930375	0.917167	0.917292	-0.206653
14	3.029655	3.031578	2.990053	2.961237	0.071519

count 값이 이상했던 8, 11, 14 그룹을 보면 8은 거래량이 유독 많았고, 11은 다른 값들이 유독 작으며, 14는 값들이 유독 큰것을 볼수 있다.

In [31]:

fig = plt.figure(figsize=(15,9))
ax = fig.add_subplot(111)
df.Close.plot(ax=ax, label='Observed', legend=True, color='b')
tdf.Close.plot(ax=ax, label='Trend', legend=True, color='r')
rdf.Close.plot(ax=ax,label='Resid', legend=True, color='y')
plt.show()

파란색은 실제 코스피지수고, 빨간색은 트렌드, 그리고 잔차가 노란색으로 표시됩니다.

COVID19 때문에 2020년 3월에 발생한 폭락장으로 인해 아래로 깊게 파인 잔차가 보입니다.

기대했던 추세보다 너무 많이 하락한 상태를 나타냅니다.

즉, 예상하지 못한 이상치가 맞다는 것을 확인할 수 있습니다.

DBSCAN 으로 이상치 찾기¶

DBSCAN에 대한 설명

In [32]:

from sklearn.cluster import DBSCAN
clustering = DBSCAN(eps=0.7, min_samples=2).fit(norm_rdf)
clustering

Out[32]:

DBSCAN(eps=0.7, min_samples=2)

In [33]:

# 분류된 라벨들은 이렇게 확인할 수 있어요
print(clustering.labels_)

[0 0 0 ... 0 0 0]

In [34]:

lbl, cnt = np.unique(clustering.labels_,return_counts=True)
['group:{}-count:{}'.format(group, count) for group, count in zip(lbl, cnt)]

Out[34]:

['group:-1-count:66',
 'group:0-count:5590',
 'group:1-count:2',
 'group:2-count:5',
 'group:3-count:6',
 'group:4-count:2',
 'group:5-count:6',
 'group:6-count:2',
 'group:7-count:8',
 'group:8-count:2',
 'group:9-count:2',
 'group:10-count:3',
 'group:11-count:2',
 'group:12-count:2',
 'group:13-count:2']

중간 정리¶

• 이상치 탐색의 기본 로직

  • Time series decomposition 을 이용해 Trend/Seasonal/Residual 값으로 분리한다
  • Residual 값의 분포를 이용해서 이상치를 탐지해낸다

• 통계적 방법과 Unsupervised Clustering 기법의 장점

  • 적은 샘플 수로도 분석이 가능하다 (샘플 수가 50개 이상이면 적용 가능)
  • 하드웨어 제약이 거의 없다

• 통계적 방법과 Unsupervised Clustering 기법의 단점

  • Time series decomposition에 의존해야 한다는 점
  • 분석자의 주관이 필요하다는 단점이 있었습니다. (통계기법:몇 배수의 표준편차?, k-means:몇 개 그룹으로 클러스터링?, Time Series decompose: freq는 몇으로?)

• 딥러닝 이상치 탐색의 장점

  • Time series decompostion 없이 분석 수행이 가능하다
  • Trend, Seasonal 데이터를 포함하고 있기 때문에 Trend, Seasonal의 변화도 이상치로써 탐색이 가능하다.
• 딥러닝 이상치 탐색의 단점
  • 분석자의 주관이 어느 정도 필요하긴 하다 (Threshold, window)
  • 학습에 활용할 만큼 충분한 데이터가 확보되어야 한다.

Auto-Encoder를 이용한 이상치 탐색¶

In [35]:

# 필요한 라이브러리를 불러옵니다

import tensorflow as tf
from tensorflow.keras.preprocessing.sequence import TimeseriesGenerator
from tensorflow.keras.models import Sequential
from tensorflow.keras.layers import Dense, Dropout, LSTM, RepeatVector, TimeDistributed
from tensorflow.keras.losses import Huber
from tensorflow.keras.callbacks import ModelCheckpoint, EarlyStopping

In [36]:

# 모델 재현성을 위해 난수 시드설정을 해줍니다
tf.random.set_seed(777)
np.random.seed(777)

LSTM 을 이용한 오토인코더 모델 만들기¶

In [37]:

from sklearn.preprocessing import StandardScaler

# 데이터 전처리 - 하이퍼 파라미터
window_size = 10
batch_size = 32
features = ['Open','High','Low','Close','Volume']
n_features = len(features)
TRAIN_SIZE = int(len(df)*0.7)

# 데이터 전처리
# 표준정규분포화합니다.

scaler = StandardScaler()
scaler = scaler.fit(df.loc[:TRAIN_SIZE,features].values)
scaled = scaler.transform(df[features].values)

keras TimeseriesGenerator 를 이용해서 간편하게 데이터 셋 만들기¶

• data : 입력데이터(x)를 넣어줍니다. 우리는 'Open','High','Low','Close','Volume' 5가지 인풋을 사용합니다.

• targets : 출력데이터를 넣어줍니다. 우리는 5가지 인풋 그대로 예측하기 때문에 data와 동일한 걸 넣어줍니다.

• length : 몇 개의 time_step을 참고할지 입력합니다.

• stride : time_step 사이즈를 결정합니다.

In [38]:

# keras TimeseriesGenerator 를 이용해서 간편하게 데이터 셋을 만듭니다

train_gen = TimeseriesGenerator(
    data = scaled,
    targets = scaled,
    length = window_size,
    stride=1,
    sampling_rate=1,
    batch_size= batch_size,
    shuffle=False,
    start_index=0,
    end_index=None,
)

valid_gen = TimeseriesGenerator(
    data = scaled,
    targets = scaled,
    length = window_size,
    stride=1,
    sampling_rate=1,
    batch_size=batch_size,
    shuffle=False,
    start_index=TRAIN_SIZE,
    end_index=None,
)

print(train_gen[0][0].shape)
print(train_gen[0][1].shape)

(32, 10, 5)
(32, 5)

모델 만들기¶

• 2개 층의 LSTM으로 인코더를 만듭니다
• RepeatVector는 input을 window_size 만큼 복사해줍니다.

In [39]:

model = Sequential([
    # >> 인코더 시작
    LSTM(64, activation='relu', return_sequences=True, 
         input_shape=(window_size, n_features)),
    LSTM(16, activation='relu', return_sequences=False),
    ## << 인코더 끝
    ## >> Bottleneck
    RepeatVector(window_size),
    ## << Bottleneck
    ## >> 디코더 시작
    LSTM(16, activation='relu', return_sequences=True),
    LSTM(64, activation='relu', return_sequences=False),
    Dense(n_features)
    ## << 디코더 끝
])

model.summary()

WARNING:tensorflow:Layer lstm will not use cuDNN kernel since it doesn't meet the cuDNN kernel criteria. It will use generic GPU kernel as fallback when running on GPU
WARNING:tensorflow:Layer lstm_1 will not use cuDNN kernel since it doesn't meet the cuDNN kernel criteria. It will use generic GPU kernel as fallback when running on GPU
WARNING:tensorflow:Layer lstm_2 will not use cuDNN kernel since it doesn't meet the cuDNN kernel criteria. It will use generic GPU kernel as fallback when running on GPU
WARNING:tensorflow:Layer lstm_3 will not use cuDNN kernel since it doesn't meet the cuDNN kernel criteria. It will use generic GPU kernel as fallback when running on GPU
Model: "sequential"
_________________________________________________________________
Layer (type)                 Output Shape              Param #   
=================================================================
lstm (LSTM)                  (None, 10, 64)            17920     
_________________________________________________________________
lstm_1 (LSTM)                (None, 16)                5184      
_________________________________________________________________
repeat_vector (RepeatVector) (None, 10, 16)            0         
_________________________________________________________________
lstm_2 (LSTM)                (None, 10, 16)            2112      
_________________________________________________________________
lstm_3 (LSTM)                (None, 64)                20736     
_________________________________________________________________
dense (Dense)                (None, 5)                 325       
=================================================================
Total params: 46,277
Trainable params: 46,277
Non-trainable params: 0
_________________________________________________________________

In [40]:

# 체크포인트
# 학습을 진행하며 validation 결과가 가장 좋은 모델을 저장해둠
import os

checkpoint_path = os.getenv('HOME')+'/aiffel/anomaly_detection/kospi/mymodel.ckpt'
checkpoint = ModelCheckpoint(checkpoint_path, 
                             save_weights_only=True, 
                             save_best_only=True, 
                             monitor='val_loss',
                             verbose=1)

# 얼리스탑
# 학습을 진행하며 validation 결과가 나빠지면 스톱. patience 횟수만큼은 참고 지켜본다

early_stop = EarlyStopping(monitor='val_loss', patience=5)

model.compile(loss='mae', optimizer='adam',metrics=["mae"])

In [41]:

hist = model.fit(train_gen, 
          validation_data=valid_gen,
          steps_per_epoch=len(train_gen),
          validation_steps=len(valid_gen),
          epochs=50, 
          callbacks=[checkpoint, early_stop])

Epoch 1/50
177/178 [============================>.] - ETA: 0s - loss: 0.5756 - mae: 0.5756
Epoch 00001: val_loss improved from inf to 0.23469, saving model to /home/ssac24/aiffel/anomaly_detection/kospi/mymodel.ckpt
178/178 [==============================] - 7s 37ms/step - loss: 0.5736 - mae: 0.5736 - val_loss: 0.2347 - val_mae: 0.2347
Epoch 2/50
177/178 [============================>.] - ETA: 0s - loss: 0.2315 - mae: 0.2315
Epoch 00002: val_loss improved from 0.23469 to 0.16202, saving model to /home/ssac24/aiffel/anomaly_detection/kospi/mymodel.ckpt
178/178 [==============================] - 6s 34ms/step - loss: 0.2312 - mae: 0.2312 - val_loss: 0.1620 - val_mae: 0.1620
Epoch 3/50
177/178 [============================>.] - ETA: 0s - loss: 0.2192 - mae: 0.2192
Epoch 00003: val_loss did not improve from 0.16202
178/178 [==============================] - 6s 34ms/step - loss: 0.2185 - mae: 0.2185 - val_loss: 0.2664 - val_mae: 0.2664
Epoch 4/50
177/178 [============================>.] - ETA: 0s - loss: 0.2057 - mae: 0.2057
Epoch 00004: val_loss improved from 0.16202 to 0.15212, saving model to /home/ssac24/aiffel/anomaly_detection/kospi/mymodel.ckpt
178/178 [==============================] - 6s 34ms/step - loss: 0.2055 - mae: 0.2055 - val_loss: 0.1521 - val_mae: 0.1521
Epoch 5/50
177/178 [============================>.] - ETA: 0s - loss: 0.1966 - mae: 0.1966
Epoch 00005: val_loss improved from 0.15212 to 0.13299, saving model to /home/ssac24/aiffel/anomaly_detection/kospi/mymodel.ckpt
178/178 [==============================] - 6s 34ms/step - loss: 0.1963 - mae: 0.1963 - val_loss: 0.1330 - val_mae: 0.1330
Epoch 6/50
177/178 [============================>.] - ETA: 0s - loss: 0.1756 - mae: 0.1756
Epoch 00006: val_loss improved from 0.13299 to 0.13220, saving model to /home/ssac24/aiffel/anomaly_detection/kospi/mymodel.ckpt
178/178 [==============================] - 6s 33ms/step - loss: 0.1753 - mae: 0.1753 - val_loss: 0.1322 - val_mae: 0.1322
Epoch 7/50
177/178 [============================>.] - ETA: 0s - loss: 0.1731 - mae: 0.1731
Epoch 00007: val_loss did not improve from 0.13220
178/178 [==============================] - 6s 33ms/step - loss: 0.1733 - mae: 0.1733 - val_loss: 0.2592 - val_mae: 0.2592
Epoch 8/50
177/178 [============================>.] - ETA: 0s - loss: 0.1704 - mae: 0.1704
Epoch 00008: val_loss improved from 0.13220 to 0.12425, saving model to /home/ssac24/aiffel/anomaly_detection/kospi/mymodel.ckpt
178/178 [==============================] - 6s 34ms/step - loss: 0.1699 - mae: 0.1699 - val_loss: 0.1242 - val_mae: 0.1242
Epoch 9/50
177/178 [============================>.] - ETA: 0s - loss: 0.1600 - mae: 0.1600
Epoch 00009: val_loss improved from 0.12425 to 0.11609, saving model to /home/ssac24/aiffel/anomaly_detection/kospi/mymodel.ckpt
178/178 [==============================] - 6s 35ms/step - loss: 0.1596 - mae: 0.1596 - val_loss: 0.1161 - val_mae: 0.1161
Epoch 10/50
177/178 [============================>.] - ETA: 0s - loss: 0.1697 - mae: 0.1697
Epoch 00010: val_loss did not improve from 0.11609
178/178 [==============================] - 6s 35ms/step - loss: 0.1698 - mae: 0.1698 - val_loss: 0.1892 - val_mae: 0.1892
Epoch 11/50
177/178 [============================>.] - ETA: 0s - loss: 0.1621 - mae: 0.1621
Epoch 00011: val_loss improved from 0.11609 to 0.11218, saving model to /home/ssac24/aiffel/anomaly_detection/kospi/mymodel.ckpt
178/178 [==============================] - 6s 35ms/step - loss: 0.1618 - mae: 0.1618 - val_loss: 0.1122 - val_mae: 0.1122
Epoch 12/50
177/178 [============================>.] - ETA: 0s - loss: 0.1524 - mae: 0.1524
Epoch 00012: val_loss did not improve from 0.11218
178/178 [==============================] - 6s 34ms/step - loss: 0.1530 - mae: 0.1530 - val_loss: 0.1173 - val_mae: 0.1173
Epoch 13/50
177/178 [============================>.] - ETA: 0s - loss: 0.1404 - mae: 0.1404
Epoch 00013: val_loss did not improve from 0.11218
178/178 [==============================] - 6s 34ms/step - loss: 0.1414 - mae: 0.1414 - val_loss: 0.1807 - val_mae: 0.1807
Epoch 14/50
177/178 [============================>.] - ETA: 0s - loss: 0.1304 - mae: 0.1304
Epoch 00014: val_loss did not improve from 0.11218
178/178 [==============================] - 6s 34ms/step - loss: 0.1302 - mae: 0.1302 - val_loss: 0.2328 - val_mae: 0.2328
Epoch 15/50
177/178 [============================>.] - ETA: 0s - loss: 0.1421 - mae: 0.1421
Epoch 00015: val_loss did not improve from 0.11218
178/178 [==============================] - 6s 33ms/step - loss: 0.1433 - mae: 0.1433 - val_loss: 0.1196 - val_mae: 0.1196
Epoch 16/50
177/178 [============================>.] - ETA: 0s - loss: 0.1389 - mae: 0.1389
Epoch 00016: val_loss improved from 0.11218 to 0.09901, saving model to /home/ssac24/aiffel/anomaly_detection/kospi/mymodel.ckpt
178/178 [==============================] - 6s 34ms/step - loss: 0.1384 - mae: 0.1384 - val_loss: 0.0990 - val_mae: 0.0990
Epoch 17/50
177/178 [============================>.] - ETA: 0s - loss: 0.1241 - mae: 0.1241
Epoch 00017: val_loss did not improve from 0.09901
178/178 [==============================] - 6s 33ms/step - loss: 0.1237 - mae: 0.1237 - val_loss: 0.1519 - val_mae: 0.1519
Epoch 18/50
177/178 [============================>.] - ETA: 0s - loss: 0.1267 - mae: 0.1267
Epoch 00018: val_loss did not improve from 0.09901
178/178 [==============================] - 6s 33ms/step - loss: 0.1276 - mae: 0.1276 - val_loss: 0.1147 - val_mae: 0.1147
Epoch 19/50
177/178 [============================>.] - ETA: 0s - loss: 0.1396 - mae: 0.1396
Epoch 00019: val_loss did not improve from 0.09901
178/178 [==============================] - 6s 33ms/step - loss: 0.1395 - mae: 0.1395 - val_loss: 0.1588 - val_mae: 0.1588
Epoch 20/50
177/178 [============================>.] - ETA: 0s - loss: 0.1359 - mae: 0.1359
Epoch 00020: val_loss did not improve from 0.09901
178/178 [==============================] - 6s 33ms/step - loss: 0.1356 - mae: 0.1356 - val_loss: 0.2043 - val_mae: 0.2043
Epoch 21/50
177/178 [============================>.] - ETA: 0s - loss: 0.1347 - mae: 0.1347
Epoch 00021: val_loss did not improve from 0.09901
178/178 [==============================] - 6s 34ms/step - loss: 0.1347 - mae: 0.1347 - val_loss: 0.1128 - val_mae: 0.1128

In [42]:

model.load_weights(checkpoint_path)

Out[42]:

<tensorflow.python.training.tracking.util.CheckpointLoadStatus at 0x7f5da86acb50>

학습 과정을 확인¶

• 안정적으로 Training loss 가 수렴하고
• Validation loss 가 발산하지 않음을 확인

In [43]:

fig = plt.figure(figsize=(12,8))
plt.plot(hist.history['loss'], label='Training')
plt.plot(hist.history['val_loss'], label='Validation')
plt.legend()

Out[43]:

<matplotlib.legend.Legend at 0x7f5e2c1dfb90>

분석 시 주의할 점!!¶

시계열 데이터를 window_size 만큼 밀어가면서 예측하는 모델을 만들었기때문에 train_gen의 길이가 원본 df의 길이보다 window_size 만큼 짧다. 예측 결과와 비교할 때는 scaled의 앞에서 window_size 만큼을 건너뛰어야한다.

In [44]:

# 예측 결과를 pred 로, 실적 데이터를 real로 받습니다
pred = model.predict(train_gen)
real = scaled[window_size:]

mae_loss = np.mean(np.abs(pred-real), axis=1)

# 샘플 개수가 많기 때문에 y축을 로그 스케일로 그립니다
fig, ax = plt.subplots(figsize=(9,6))
_ = plt.hist(mae_loss, 100, density=True, alpha=0.75, log=True)

In [47]:

import copy

test_df = copy.deepcopy(df.loc[window_size:]).reset_index(drop=True)
test_df['Loss'] = mae_loss

threshold = 0.3
test_df.loc[test_df.Loss>threshold]

Out[47]:

	Date	Open	High	Low	Close	Adj Close	Volume	Loss
734	2000-07-11	856.169983	863.530029	836.859985	836.859985	836.859985	850000.0	0.315499
738	2000-07-18	822.299988	829.190002	806.880005	812.330017	812.330017	317100.0	0.310300
946	2001-05-28	620.200012	622.669983	618.140015	618.469971	618.469971	434800.0	0.312367
1014	2001-09-03	545.020020	547.630005	539.169983	541.830017	541.830017	415200.0	0.315305
1016	2001-09-05	550.429993	553.609985	547.880005	551.909973	551.909973	788700.0	0.312333
...	...	...	...	...	...	...	...	...
5634	2020-05-28	2047.079956	2054.520020	2003.750000	2028.540039	2028.540039	1172100.0	0.345942
5638	2020-06-03	2108.550049	2156.550049	2107.689941	2147.000000	2147.000000	1152000.0	0.320892
5639	2020-06-04	2181.639893	2191.000000	2139.679932	2151.179932	2151.179932	1393300.0	0.498471
5646	2020-06-15	2114.409912	2129.669922	2030.819946	2030.819946	2030.819946	1071600.0	0.379944
5658	2020-07-01	2128.810059	2133.550049	2101.330078	2106.699951	2106.699951	1116200.0	0.322882

176 rows × 8 columns

그래프로 그려서 이상치를 찾아보기¶

Open, Close, Low , High 같은 데이터는 스케일이 비슷하니 한 번에 그려보고, Volume과 loss는 스케일이 다르니 각각 그래프를 그려보는것이 좋다. 기준치로 분류해낸 이상치들은 붉은 점으로 그려서 나타내는것이 좋다.

In [48]:

fig = plt.figure(figsize=(12,15))

# 가격들 그래프입니다
ax = fig.add_subplot(311)
ax.set_title('Open/Close')
plt.plot(test_df.Date, test_df.Close, linewidth=0.5, alpha=0.75, label='Close')
plt.plot(test_df.Date, test_df.Open, linewidth=0.5, alpha=0.75, label='Open')
plt.plot(test_df.Date, test_df.Close, 'or', markevery=[mae_loss>threshold])

# 거래량 그래프입니다
ax = fig.add_subplot(312)
ax.set_title('Volume')
plt.plot(test_df.Date, test_df.Volume, linewidth=0.5, alpha=0.75, label='Volume')
plt.plot(test_df.Date, test_df.Volume, 'or', markevery=[mae_loss>threshold])

# 오차율 그래프입니다
ax = fig.add_subplot(313)
ax.set_title('Loss')
plt.plot(test_df.Date, test_df.Loss, linewidth=0.5, alpha=0.75, label='Loss')
plt.plot(test_df.Date, test_df.Loss, 'or', markevery=[mae_loss>threshold])

Out[48]:

[<matplotlib.lines.Line2D at 0x7f5da0585350>]

/home/ssac24/anaconda3/envs/aiffel/lib/python3.7/site-packages/matplotlib/lines.py:191: FutureWarning: Using a non-tuple sequence for multidimensional indexing is deprecated; use `arr[tuple(seq)]` instead of `arr[seq]`. In the future this will be interpreted as an array index, `arr[np.array(seq)]`, which will result either in an error or a different result.
  return Path(verts[markevery], _slice_or_none(codes, markevery))