프로그래밍 코드를 짜다보면 수행 절차나 방법, 사용하는 메소드에 따라서 수행 시간이 차이가 나는 경우가 종종 있습니다. 그리고 성능이 중요해서 여러가지 방법을 테스트해보면서 가장 실행시간이 짧도록 튜닝하면서 최적화하기도 합니다. 

 

이번 포스팅에서는 Python에서 코드를 실행시켰을 때 소요된 시간을 측정하는 2가지 방법을 소개하겠습니다. 

 

(1) datetime.now() 메소드 이용해서 실행 시간 측정하기

(2) %timeit 로 실행 시간 측정하기

 

 

python : measuring the execution time of code snippets

 

먼저, 예제로 사용할 샘플 데이터셋으로서, 1억개의 값을 가지는 xarr, yarr 의 두개의 배열(array)를 만들어 보겠습니다. 그리고 배열 내 각 1억개의 값 별로 True/False 의 조건값을 가지는 cond 라는 배열도 난수를 생성시켜서 만들어보겠습니다. 

 

import numpy as np

## generating sample array with 100 million values
xarr = np.arange(100000000)
yarr = np.zeros(100000000)
cond = np.where(np.random.randn(100000000)>0, True, False)


cond[:10]
# [Out] array([False,  True,  True, False, False,  True,  True,  True,  
#              True, True])

 

 

위에서 만든 1억개의 원소를 가지는 배열을 가지고 조건값으로 True/False 블리언 값 여부에 따라서 True 조건값 이면 xarr 배열 내 값을 가지고, False 조건값이면 yarr 배열 내 값을 가지는 새로운 배열을 만들어보겠습니다. 이때 (1) List Comprehension 방법과, (2) NumPy의 Vectorized Operations 방법 간 수행 시간을 측정해서 어떤 방법이 더 빠른지 성능을 비교해보겠습니다.

(물론, Vectorized Operations이 for loop 순환문을 사용하는 List Comprehension보다 훨~씬 빠릅니다! 눈으로 직접 확인해 보겠습니다. )

 

## Let's compare the elapsed time between 2 methods 
## (list comprehension vs. vectorized operations)

## (1) List Comprehension
new_arr = [(x if c else y) for (x, y, c) in zip(xarr, yarr, cond)]

## (2) Vectorized Operations in NumPy 
new_arr = np.where(cond, xarr, yarr)

 

 

 

(1) datetime.now() 메소드 이용해서 실행 시간 측정하기

 

datetime 모듈은 날짜, 시간, 시간대(time zone) 등을 다루는데 사용하는 모듈입니다 datetime.now() 메소드는 현재의 로컬 날짜와 시간을 반환합니다. 실행 시간을 측정할 코드 바로 앞에 start_time = datetime.now() 로 시작 날짜/시간을 측정해놓고, 실행할 코드가 끝난 다음 줄에 time_elapsed = datetime.now() - start_time 으로 '끝난 날짜/시간'에서 '시작 날짜/시간'을 빼주면 '코드 실행 시간'을 계산할 수 있습니다. 

 

아래 결과를 비교해보면 알 수 있는 것처럼, for loop 순환문을 사용하는 List Comprehension 방법보다 NumPy의 Vectorized Operation이 약 38배 빠른 것으로 나오네요. 

 

## (1) -- measuring the elapsed time using datetime

## (a) List Comprehension
from datetime import datetime
start_time = datetime.now() 
list_comp_for_loop = [(x if c else y) for (x, y, c) in zip(xarr, yarr, cond)]
time_elapsed = datetime.now() - start_time 

print('Time elapsed (hh:mm:ss.ms) {}'.format(time_elapsed))
# Time elapsed (hh:mm:ss.ms) 0:00:17.753036

np.array(list_comp_for_loop)[:10]
# array([0., 1., 2., 0., 0., 5., 6., 7., 8., 9.])



## (b) Vectorized Operations in NumPy 
start_time = datetime.now() 
np_where_vectorization = np.where(cond, xarr, yarr)
time_elapsed = datetime.now() - start_time 

print('Time elapsed (hh:mm:ss.ms) {}'.format(time_elapsed))
# Time elapsed (hh:mm:ss.ms) 0:00:00.462215

np_where_vectorization[:10]
# array([0., 1., 2., 0., 0., 5., 6., 7., 8., 9.])

 

 

(2) %timeit 로 실행 시간 측정하기

 

다음으로 Python timeit 모듈을 사용해서 짧은 코드의 실행 시간을 측정해보겠습니다. timeit 모듈은 터미널의 command line 과 Python IDE 에서 호출 가능한 형태의 코드 둘 다 사용이 가능합니다. 

 

아래에는 Jupyter Notebook에서 %timeit [small code snippets] 로 코드 수행 시간을 측정해본 예인데요, 여러번 수행을 해서 평균 수행 시간과 표준편차를 보여주는 특징이 있습니다. 

 

## (2) measuring the elapsed time using timeit

## (a) List Comprehension
import timeit

%timeit list_comp_for_loop = [(x if c else y) for (x, y, c) in zip(xarr, yarr, cond)]
# 17.1 s ± 238 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)


## (b) Vectorized Operations in NumPy 
%timeit np_where_vectorization = np.where(cond, xarr, yarr)
# 468 ms ± 8.75 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

 

 

[Reference]

* Python datetime: https://docs.python.org/3/library/datetime.html

* Python timeit: "measuring the execution time of small code snippets"
   : https://docs.python.org/3/library/timeit.html

 

 

이번 포스팅이 많은 도움이 되었기를 바랍니다. 

행복한 데이터 과학자 되세요! 

 

728x90
반응형
Posted by Rfriend
,

Word counts 할 때 많이 사용하는 코드인데요, 이번 포스팅에서는 

  (1) 리스트에서 원소별 개수를 세서 {Key:Value} 쌍의 Dictionary를 만들고 

  (2) 원소별 개수를 세어놓은 Dictionary에서 개수 상위 n 개의 {Key:Value} 쌍을 가져오기

하는 방법을 소개하겠습니다. 

 

 

counts dictionary, getting top n

 

(1) 리스트에서 원소별 개수를 세서 {Key:Value} 쌍의 Dictionary를 만들기

 

먼저, 예제로 사용할 간단한 리스트를 만들어보겠습니다. 

 

## creating sample lists
my_list = ['a', 'f', 'a', 'b', 'a', 'a', 'c', 'b', 
           'c', 'e', 'a', 'c', 'b', 'f', 'c']
           
print(my_list)
# ['a', 'f', 'a', 'b', 'a', 'a', 'c', 'b', 'c', 'e', 'a', 'c', 'b', 'f', 'c']

 

 

다음으로, 원소별 개수를 세서 저장할 비어있는 Dictionary 인 counts={} 를 만들어놓고, for loop 순환문으로 리스트의 원소를 하나씩 순서대로 가져다가 Dictionary counts 의 Key 값에 해당 원소가 들어있으면 +1을 하고, Key 값에 해당 원소가 안들어 있으면 해당 원소를 Key 값으로 등록하고 1 을 값으로 입력해 줍니다.  

 

def get_counts(seq): 
    counts = {}
    for x in seq:
        if x in counts:
            counts[x] += 1
        else:
            counts[x] = 1
    return counts
    
 
counts = get_counts(my_list)


print(counts)
# {'a': 5, 'f': 2, 'b': 3, 'c': 4, 'e': 1}


## access value by key
counts['a']
# 5

 

 

 

(2) 원소별 개수를 세어놓은 Dictionary에서 개수 상위 n 개의 {Key:Value} 쌍을 가져오기

 

Dictionary를 정렬하는 방법에 따라서 두 가지 방법이 있습니다.

 

(a) sorted() 메소드를 이용해서 key=lambda x: x[1] 로 해서 정렬 기준을 Dictionary의 Value 로 하여 내림차순으로 정렬(reverse=True) 하고, 상위 n 개까지만 슬라이싱해서 가져오는 방법입니다. 

 

## way 1
## reference: https://rfriend.tistory.com/473
def top_n(count_dict, n=3):
    return sorted(count_dict.items(), reverse=True, key=lambda x: x[1])[:n]
    

## getting top 2
top_n(counts, n=2)
# [('a', 5), ('c', 4)]

 

 

 

(b) 아래는 dict.items() 로 (Key, Value) 쌍을 for loop 문을 돌리면서 (Value, Key) 로 순서를 바꿔서 리스트 [] 로 만들고 (list comprehension), 이 리스트에 대해서 sort(reverse=True) 로 Value 를 기준으로 내림차순 정렬한 후에, 상위 n 개까지만 슬라이싱해서 가져오는 방법입니다. 

 

## way2
## reference: https://rfriend.tistory.com/281
def top_n2(count_dict, n=3):
    val_key = [(v, k) for k, v in count_dict.items()]
    val_key.sort(reverse=True)
    return val_key[:n]
    
## getting top 2
top_n2(counts, n=2)
# [(5, 'a'), (4, 'c')]

 

 

[Reference]

* Dictionary 정렬: https://rfriend.tistory.com/473

* List 정렬: https://rfriend.tistory.com/281

 

 

이번 포스팅이 많은 도움이 되었기를 바랍니다. 

행복한 데이터 과학자 되세요. 

 

728x90
반응형
Posted by Rfriend
,

이번 포스팅에서는 Python pandas 의 DataFrame에서 

 

(1) 그룹별로 x 칼럼을 정렬 후 누적 비율을 구한 후 (calculating the cumulative proportion by groups)  

(2) 그룹별로 특정 분위수의 위치 구하기 (getting the indices for a specific quantile p by groups)

 

하는 방법을 소개하겠습니다. 

 

그룹별로 연산을 수행하므로 pandas.DataFrame.groupby().apply(UDF) 형식으로 구문을 작성할 거예요. 

 

 

[ pandas DataFrame에서 그룹별로 정렬 후 누적 비율을 구한 후에 --> 그룹별로 특정 분위수 위치 구하기 ]

pandas getting cumulative proportion and quantile

 

 

 

먼저, 예제로 사용하기 위해 그룹('grp') 칼럼별 값('x')을 가지는 간단한 pandas DataFrame을 만들어보겠습니다. 

 

import numpy as np
import pandas as pd

df = pd.DataFrame({
    'grp': ['a', 'a', 'a', 'a', 'a', 'b', 'b', 'b', 'b', 'b'], 
    'x': [3, 1, 2, 7, 4, 4, 2, 5, 9, 7]})
    
print(df)
#   grp  x
# 0   a  3
# 1   a  1
# 2   a  2
# 3   a  7
# 4   a  4
# 5   b  4
# 6   b  2
# 7   b  5
# 8   b  9
# 9   b  7

 

 

 

(1) 그룹별로 x 칼럼을 정렬 후 누적 비율을 구한 후 (calculating the cumulative proportion by groups)  

 

그룹별로 x 칼럼에 대한 누적 비율을 구하기 위해, 먼저 그룹별로 x 칼럼의 비율(proportion)을 계산해서 'prop' 라는 칼럼을 추가해보겠습니다. x_prop() 라는 사용자 정의 함수를 정의한 후, df.groupby('grp').apply(x_prop) 처럼 그룹에 apply() 메소드로 사용자 정의 함수를 적용해서 연산을 했습니다. 

 

## adding the proportion column by group
def x_prop(group):
    group['prop'] = group.x / group.x.sum()
    return group

df = df.groupby('grp').apply(x_prop)


print(df)
#   grp  x      prop
# 0   a  3  0.176471
# 1   a  1  0.058824
# 2   a  2  0.117647
# 3   a  7  0.411765
# 4   a  4  0.235294
# 5   b  4  0.148148
# 6   b  2  0.074074
# 7   b  5  0.185185
# 8   b  9  0.333333
# 9   b  7  0.259259


## checking the sanity
df.groupby('grp').prop.sum()
#      grp
# a    1.0
# b    1.0
# Name: prop, dtype: float64

 

 

앞에서 계산한 그룹별 x 칼럼의 비율 'prop'을 그룹별로 내림차순(descending order)으로 정렬해서 보면 아래와 같습니다. 

 

## sorting in descending order by prop 
df.sort_values(by=['grp', 'prop'], ascending=False)

#     grp	x	prop
# 8	b	9	0.333333
# 9	b	7	0.259259
# 7	b	5	0.185185
# 5	b	4	0.148148
# 6	b	2	0.074074
# 3	a	7	0.411765
# 4	a	4	0.235294
# 0	a	3	0.176471
# 2	a	2	0.117647
# 1	a	1	0.058824

 

 

pandas 의 cumsum() 메소드를 사용해서 그룹별 x칼럼의 비율 'prop'의 누적 합계 (cumulative sum) 인 'cum_prop' 를 그룹별로 계산해보겠습니다. 역시 비율 'prop'에 대해서 누적 합계(cum_sum)를 구하는 사용자 정의 함수 'cum_prop()'를 먼저 정의한 후에, 이를 df.groupby('grp').apply(cum_prop) 처럼 apply() 메소드에 사용자 정의함수를 적용해서 계산했습니다. 

 

## sorting in descending order by prop and calculating the cumulative sum of prop
def cum_prop(group):
    group['cum_prop'] = group.sort_values(
        by='prop', ascending=False).prop.cumsum()
    return group

df = df.groupby('grp').apply(cum_prop)


df.sort_values(by=['grp', 'cum_prop'])

#     grp	x	prop	        cum_prop
# 3	a	7	0.411765	0.411765
# 4	a	4	0.235294	0.647059
# 0	a	3	0.176471	0.823529
# 2	a	2	0.117647	0.941176
# 1	a	1	0.058824	1.000000
# 8	b	9	0.333333	0.333333
# 9	b	7	0.259259	0.592593
# 7	b	5	0.185185	0.777778
# 5	b	4	0.148148	0.925926
# 6	b	2	0.074074	1.000000

 

 

 

위의 예시는 간단한 편이므로 아래처럼 사용자 정의 함수를 정의하는 대신에 apply() 메소드 안에 바로 lambda 로 'prop'에 대해서 내림차순 정렬 후 누적 합계를 구하는 함수를 바로 써줘도 됩니다. 

 

## or, equivalentsly, using lambda function for cumulative proportion
df.groupby('grp').apply(lambda x: x.sort_values(by='prop', ascending=False).prop.cumsum())

# grp   
# a    3    0.411765
#      4    0.647059
#      0    0.823529
#      2    0.941176
#      1    1.000000
# b    8    0.333333
#      9    0.592593
#      7    0.777778
#      5    0.925926
#      6    1.000000
# Name: prop, dtype: float64

 

 

 

(2) 그룹별로 특정 분위수의 위치 구하기 (getting the indices for a specific quantile p by groups)

 

이제 위에서 구한 '그룹별 비율의 누적 합계('cum_prop')'에 대해서 pandas.Series.searchsorted(values, side='left') 메소드를 사용해서 특정 비율이 들어갈 위치를 구해보겠습니다.

 

비율에 대해 내림차순 정렬 후 누적 합계를 구한 값에 대해 특정 값이 들어갈 위치를 구하는 것이므로, 결과적으로 자료 크기 순서에 따른 위치값인 분위수(quantile) 를 구할 수 있게 됩니다. 인덱스가 '0'부터 시작하므로 위치를 구하기 위해서 반환받는 값에 '+1' 을 해주었습니다. 

 

그룹별로 특정 분위수의 위치를 구하고 싶으므로, 분위수를 구하는 사용자 정의 함수인 'quantile_idx()' 함수를 정의한 후에, 이를 df.groupby('grp').apply(quantile_idx, p) 처럼 apply() 메소드에 사용자 정의 함수와 매개변수 p를 입력해서 적용해주면 되겠습니다.

 

그룹별로 분위수 p=0.2, p=0.5, p=0.7 인 곳의 위치를 구해보니 잘 작동하는군요. 

 

## pandas.Series.searchsorted
## Series.searchsorted(value, side='left', sorter=None)[source]
## Find indices where elements should be inserted to maintain order.

def quantile_idx(group, p=0.5):
    group = group.sort_values(by='cum_prop', ascending=True)
    return group.cum_prop.searchsorted(p) + 1
    

## getting the index of quantile p=0.2 by group
df.groupby('grp').apply(quantile_idx, p=0.2)

# grp
# a    1
# b    1
# dtype: int64


## getting the index of quantile p=0.5 by group
df.groupby('grp').apply(quantile_idx, p=0.5)

# grp
# a    2
# b    2
# dtype: int64


## getting the index of quantile p=0.7 by group
df.groupby('grp').apply(quantile_idx, p=0.7)

# grp
# a    3
# b    3
# dtype: int64

 

 

[Reference]

* pandas.Series.searchsorted() method:  https://pandas.pydata.org/docs/reference/api/pandas.Series.searchsorted.html

 

 

이번 포스팅이 많은 도움이 되었기를 바랍니다. 

행복한 데이터 과학자 되세요! :-)

 

728x90
반응형
Posted by Rfriend
,

이번 포스팅에서는 

 (1) 그룹별로 x 칼럼을 기준으로 내림차순 정렬 후 (sorting by x in ascending order)

 (2) 그룹별로 y 칼럼의 첫번째 값, 마지막 값을 DataFrame에 칼럼 추가하기

하는 2가지 방법을 소개하겠습니다. 

 

(방법 1) pandas.DataFrame 의 transform('first', 'last') 메소드를 사용하는 방법

(방법 2) 그룹별 y의 첫번째 값, 마지막 값을 구해 DataFrame을 만들고, merge() 메소드로 합치는 방법

 

pandas DataFrame sort_values() groupby() transform('first') transform('last')

 

 

먼저, 예제로 사용할 간단한 DataFrame을 만들어보겠습니다. 

 

import numpy as np
import pandas as pd

df = pd.DataFrame({
    'grp': ['a', 'a', 'a', 'b', 'b', 'b'], 
    'x': [2, 3, 1, 4, 6, 5], 
    'y': [10, 20, 30, 40, 50, 60]
}) 

df

#   grp	x	y
# 0	a	2	10
# 1	a	3	20
# 2	a	1	30
# 3	b	4	40
# 4	b	6	50
# 5	b	5	60

 

 

 

(방법 1) pandas.DataFrame 의 transform('first', 'last') 메소드를 사용하는 방법

 

그룹별로 'x' 칼럼을 기준으로 내림차순으로 정렬하려면 df.sort_values(by=['grp', 'x']) 메소드를 사용합니다. 

 

## sorting by 'grp' and 'x' in ascnding order
df.sort_values(by=['grp', 'x'])

#   grp	x	y
# 2	a	1	30
# 0	a	2	10
# 1	a	3	20
# 3	b	4	40
# 5	b	5	60
# 4	b	6	50

 

 

그러면, 이제 x를 기준으로 내림차순 정렬한 후에 'grp' 그룹별로 y 칼럼의 첫번째 값('first')과 마지막 값('last')을 가져다가 기존의 df DataFrame에 새로운 칼럼을 추가해 보겠습니다. groupby('grp') 메소드로 'grp' 그룹별 연산을 하게 되고, transform('first')는 첫번째 값을 가져다가 DataFrame에 칼럼을 추가하며, transform('last')는 마지막 값을 가져다가 DataFrame에 칼럼을 추가합니다. 

 

## adding columns of the first and last value of y by group
df['y_first'] = df.sort_values(by=['grp', 'x'])\
    .groupby('grp').y.transform('first')
df['y_last'] = df.sort_values(by=['grp', 'x'])\
    .groupby('grp').y.transform('last')


 df
# 	grp	x	y	y_first	y_last
# 0	a	2	10	30	20
# 1	a	3	20	30	20
# 2	a	1	30	30	20
# 3	b	4	40	40	50
# 4	b	6	50	40	50
# 5	b	5	60	40	50

 

 

 

(방법 2) 그룹별 y의 첫번째 값, 마지막 값을 구해 DataFrame을 만들고, merge() 메소드로 합치는 방법

 

두번째 방법은 그룹별로 x를 기준으로 정렬 후 그룹별로 y 값의 첫번째 값과 마지막 값을 구해서 별도의 DataFrame을 만든 후에, 이를 원래의 DataFrame에 merge() 하는 것입니다. DB의 테이블을 join 하는 것과 유사한 방식이예요. 

 

## creating a sample DataFrame with 2 groups
df = pd.DataFrame({
    'grp': ['a', 'a', 'a', 'b', 'b', 'b'], 
    'x': [2, 3, 1, 4, 6, 5], 
    'y': [10, 20, 30, 40, 50, 60]
}) 


## making a DataFrame with the first and last values of y by groups
y_first = df.sort_values(by='x').groupby('grp').y.first()
y_last = df.sort_values(by='x').groupby('grp').y.last()

df_grp_fst_lst = pd.DataFrame({
    'y_first': y_first, 
    'y_last': y_last
})

df_grp_fst_lst
#     y_first	y_last
# grp
# a	      30	  20
# b	      40	  50

 

 

 

pd.merge(DataFrame1, DataFrame2, how='left', on='key') 방식으로 key를 기준으로 Left Join 하면 되겠네요. 

 

## merging df_grp_fst_lst to df DataFrame by left join on 'grp'
df2 = pd.merge(df, df_grp_fst_lst, how='left', on='grp') 
# or, equivalently: df2= df.merge(df_grp_fst_lst, how='left', on='grp')


df2
# 	grp	x	y	y_first	y_last
# 0	a	2	10	30	20
# 1	a	3	20	30	20
# 2	a	1	30	30	20
# 3	b	4	40	40	50
# 4	b	6	50	40	50
# 5	b	5	60	40	50

 

 

* pandas DataFrame merge(): https://rfriend.tistory.com/258

* pandas DataFrame transform(): https://rfriend.tistory.com/403

 

 

이번 포스팅이 많은 도움이 되었기를 바랍니다. 

행복한 데이터 과학자 되세요!  :-)

 

728x90
반응형
Posted by Rfriend
,

이번 포스팅에서는 시계열 데이터에서

  (1) TimeStamp 행별로 칼럼별 비율을 구하고 

  (2) 시도표 (time series plot) 를 그리기

하는 방법을 소개하겠습니다. 

 

먼저, 예제로 사용할 간단한 pandas DataFrame을 만들어보겠습니다. index 로 2000년 ~ 2021년까지의 년도를 사용하고, 성별로 'M', 'F'의 두 개의 칼럼에 포아송분포로 부터 난수를 발생시켜서 만든 도수(frequency)를 가지는 DataFrame 입니다. 

 

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt


## creating a sample DataFrame
ts = np.arange(2000, 2022) # from year 2000 to 2021

print(ts)
# [2000 2001 2002 2003 2004 2005 2006 2007 2008 2009 2010 2011 2012 2013
#  2014 2015 2016 2017 2018 2019 2020 2021]

np.random.seed(1) # for reproducibility
M = np.arange(len(ts)) + np.random.poisson(lam=10, size=len(ts))
F = np.arange(len(ts))[::-1] + np.random.poisson(lam=2, size=len(ts))

df = pd.DataFrame({'M': M, 'F': F}, index=ts)

df.head()
#       M	F
# 2000	9	21
# 2001	7	24
# 2002	9	20
# 2003	12	19
# 2004	13	18

 

 

(1) TimeStamp 행별로 칼럼별 비율을 구하기

 

먼저, pandas DataFrame 에서 합(sum)을 구할 때 각 TimeStamp 별로 칼럼 축(axis = 1) 으로 합을 구해보겠습니다. 

 

## calculating summation by rows
df.sum(axis=1).head()

# 2000    30
# 2001    31
# 2002    29
# 2003    31
# 2004    31
# dtype: int64

 

 

참고로, index 축으로 칼럼별 합을 구할 때는 df.sum(axis=0) 을 해주면 됩니다. sum(axis=0) 이 기본설정값이므로 df.sum() 하면 동일한 결과가 나옵니다. 

 

## summation by index axis
df.sum(axis=0) # default setting

# M    426
# F    274
# dtype: int64

 

 

pandas DataFrame에서 div() 메소드를 사용하면 각 원소를 특정 값으로 나눌 수 있습니다. 가령, 위의 예제 df DataFrame의 각 원소를 10으로 나눈다고 하면 아래처럼 df.div(10) 이라고 해주면 됩니다. (나누어주는 값 '10' 이 broadcasting 되어서 각 원소를 나누어주었음.)

 

## pd.DataFrame.div()
## : Get Floating division of dataframe and other, 
##   element-wise (binary operator truediv).
df.div(10).head()

#         M	F
# 2000	0.9	2.1
# 2001	0.7	2.4
# 2002	0.9	2.0
# 2003	1.2	1.9
# 2004	1.3	1.8

 

 

이제 df DataFrame의 각 원소를 각 원소가 속한 TimeStamp별로 칼럼 축(axis=1)으로 합한 값(df.sum(axis=1))으로 나누어주면 우리가 구하고자 하는 각 TimeStamp별 칼럼별 비율을 구할 수 있습니다. 

 

df.div(df.sum(axis=1), axis=0).head()

#         M	        F
# 2000	0.300000	0.700000
# 2001	0.225806	0.774194
# 2002	0.310345	0.689655
# 2003	0.387097	0.612903
# 2004	0.419355	0.580645

 

 

 

(2) 시도표 (time series plot) 를 그리기

 

pandas DataFrame 의 plot() 메소드를 사용하면 편리하게 시계열 도표를 그릴 수 있습니다. 이때 성별을 나타내는 칼럼 'M', 'F' 별로 선의 모양(line type)과 색깔(color) 을 style={'M': 'b--', 'F': 'r-'} 매개변수를 사용해서 다르게 해서 그려보겠습니다. ('M' 은 파란색 점선, 'F' 는 빨간색 실선)

 

df.div(df.sum(1), axis=0).plot(
    style={'M': 'b--', 'F': 'r-'}, 
    figsize=(12, 8), 
    title='Proportion Trend by Gender')

plt.show()

proportion trend plot by gender

 

[ Reference ]

* pandas.DataFrame.div(): https://pandas.pydata.org/docs/reference/api/pandas.DataFrame.div.html

 

 

이번 포스팅이 많은 도움이 되었기를 바랍니다. 

행복한 데이터 과학자 되세요!  :-)

 

728x90
반응형
Posted by Rfriend
,

GitHub의 Issue와 댓글은 GitHub Flavored Markdown (GFM) 의 형식으로 작성합니다. GFM을 이용하면 본문 내용이나 코드를 강조하거나, 테이블을 만들거나, 선을 그리거나, 링크를 삽입하는 것 등을 할 수 있습니다. 이번 포스팅에서는 

 

1. GitHub Flavored Markdown (GFM) 은 무엇인가? 

2. GFM 을 사용해 Syntax Highlight 하기

3. GFM 을 사용해 그림 첨부하기

 

하는 방법을 소개하겠습니다.

 

 

1. GitHub Flavored Markdown (GFM) 은 무엇인가?

 

GitHub Flavored Markdown 은 줄여서 GFM 이라고도 하며, 현재 GitHub.com과 GitHub Enterprise 의 사용자 컨텐츠를 지원하는 Markdown의 한 종류(dialect)입니다.

 

CommonMark Spec에 기반을 둔 GFM formal specification 은 이 구문과 문법을 정의합니다. GFM 은 CommonMark의 엄격한 전체집합(strict superset of CommonMark) 입니다. GitHub 사용자 컨텐츠에서 지원되지만 원본의 commonMark Spec 에는 구체화되어 있지않는 기능은 확장자(extensions) 로 알려져 있습니다. 

 

비록 GFM 이 넓은 범위의 인풋을 지원하지만, GitHub.com과 GitHub Enterprise는 보안과 웹사이트의 일관성을 보장하기 위해 GFM 이 HTML로 변환된 후에도 사후처리(post-processing)과 건전하게 만들기(sanitization)을 추가적으로 수행합니다. 

 

 

 

2. GFM 을 사용해 Syntax Highlight 하기

 

GFM을 사용하면 Python, R, Java, C 등의 프로그래밍 언어의 구문 (Syntax) 을 좀더 가독성을 높여줄 수 있도록 각 프로그래밍 언어의 구문에 맞게 Syntax Highlight 를 해줄 수 있습니다. 

 

아래에 "Open a pull request" 의 Write 코너에 Python과 R을 사용해서 간단한 구문을 작성해보았습니다. 

틸드 기호 세개를 써주고 그 다음에 이어서 프로그래밍 언어 이름을 써주면 됩니다. 

오른쪽 하단에 보면 "Styling with Markdown is supported" 라는 문구를 볼 수 있습니다. 

참고로 '#' 은 HTML 의 H1 (헤드라인 제일 큰 글자 크기) 를 의미하며, '-' 는 점 구분 포인트입니다. 

 

GitHub Flavored Markdown: Syntax Highlight

 

 

참고로, 틸드 기호(til'de) 기호는 키보드의 좌측 상단의 (~, ₩) 단추를 눌러주면 됩니다. 우측 중간의 작은 따옴표(') 아니예요. 조심하세요. 

 

 

GFM 으로 'Write'메뉴의 본문에서 Syntax Highlight 를 해서 문서를 작성했다면, 상단의 'Preview' 메뉴에서 미리보기를 하여 확인해볼 수 있습니다. 아래에 보니 Python과 R 언어의 Syntax Highlight 가 잘 되었음을 확인할 수 있네요. 

 

GitHub Flavored Markdown: Preview

 

 

3. GFM 을 사용해 그림 첨부하기

 

GitHub에 파일을 첨부하거나 그림을 삽입할 수도 있습니다. 아래의 화면캡쳐에서 보는 것처럼, 간단하게 첨부하고자 하는 파일을 아래의 네모 칸에 드래그 & 드랍(dragging & dropping) 하거나, 아니면 아래의 네모 칸을 클릭한 후에 탐색기에서 해당 파일을 찾아서 첨부를 해주면 됩니다. 

GitHub: Attach files by dragging & dropping, selecting or pasting them

 

 

아래의 예제 화면은 'R Friend Logo' 이미지 파일을 첨부해한 후에 'Preview'로 미리보기 확인을 해 본 것입니다. 이미지가 잘 삽입되었네요. 

 

 

 

 

[ Reference ]

- GitHub Flavored Markdown Spec: https://github.github.com/gfm/

 

 

이번 포스팅이 많은 도움이 되었기를 바랍니다. 

행복한 데이터 과학자 되세요! :-)

 

728x90
반응형
Posted by Rfriend
,

맥북은 Git 이 기본으로 설치되어 있습니다. 맥북에 설치된 Git 의 기본 설정하는 방법으로서, 

  (1) Git 사용자 이름 설정

  (2) Git 이메일 주소 설정

  (3) Git 출력되는 명령어의 색깔 설정

  (4) Git 기본설정 확인 

하는 방법을 소개하겠습니다. 

 

Git 기본설정 (git global setting)

 

 

(1) Git 사용자 이름 설정하기

$ git config --global user.name "YourFirstName LastName"

굳이 실명을 입력하지 않아도 되며, 별명을 입력해도 됩니다. 이름은 영어로 입력해 주세요. 

 

 

 

(2) Git 이메일 주소 설정하기

$ git config --global user.email "your_email@exmaple.com"

위에서처럼 Git 사용자 이름과 이메일 주소를 설정해주면 ~/.gitconfig 에 설정값이 저장됩니다. 

 

 

 

(3) Git 출력되는 명령어의 색깔 설정하기

$ git config --global color.ui auto

출력되는 명령어를 읽기 쉽도록 하기 위해서 색깔 설정을 자동으로 해줍니다. (color.ui auto)

 

 

 

(4) Git 기본 설정 확인하기

 

셸 명령어에서 cat 으로 ~/.gitconfig 에 저장되어 있는 설정값을 읽어와서 확인해볼 수 있습니다.  

$ cat ~/.gitconfig
[user]
	name = YourFirstName LastName
	email = your_email@example.com
[color]
	ui = auto
$
$

 

 

Git 기본 설정값을 지정하거나 혹은 수정하고자 할 때 vim ~/.gitconfig 편집모드에서 직접 입력해줄 수도 있습니다. vim 편집모드에서 사용자 정보와 색깔 설정을 입력/수정하고 나면  'esc + :wq!' 로 변경사항을 저장한 후에 vim 편집모드에서 나오면 됩니다. 

$ vim ~/.gitconfig


[user]
        name = YourFirstName LastName
        email = your_email@example.com
[color]
        ui = auto
~                                                                                                                                                      
~                                                                                                                                                      
:wq!

 

Github(https://github.com)에 회원가입을 한 후에, Github에서 SSH-Key를 생성하고 SSH public key를 등록하는 방법은 https://rfriend.tistory.com/603 를 참고하세요. 

 

많은 도움이 되었기를 바랍니다. 

행복한 데이터 과학자 되세요! :-)

 

728x90
반응형
Posted by Rfriend
,

이전 포스팅에서는 

 (i) 정상확률과정(stationary process)의 정의 (https://rfriend.tistory.com/691)

 (ii) 통계적 가설 검증을 통한 시계열 정상성(stationarity test) 여부 확인 (https://rfriend.tistory.com/694)

하는 방법을 소개하였습니다. 

 

ARIMA 모형과 같은 통계적 시계열 예측 모델의 경우 시계열데이터의 정상성 가정을 충족시켜야 합니다. 따라서 만약 시계열 데이터가 비정상 확률 과정 (non-stationary process) 이라면, 먼저 시계열 데이터 변환을 통해서 정상성(stationarity)을 충족시켜주어야 ARIMA 모형을 적합할 수 있습니다

 

이번 포스팅에서는 Python을 사용하여 

(1) 분산이 고정적이지 않은 경우 분산 안정화 변환 (variance stabilizing transformation, VST)

(2) 추세가 있는 경우 차분을 통한 추세 제거 (de-trend by differencing)

(3) 계절성이 있는 경우 계절 차분을 통한 계절성 제거 (de-seasonality by seaanl differencing)

하는 방법을 소개하겠습니다. 

 

 

[ 비정상확률과정을 정상확률과정으로 변환하기 (Transforming non-stationary to stationary process) ]

 

 

먼저 예제로 사용할 약 판매량 (drug sales) 시계열 데이터를 가져와서 pandas DataFrame으로 만들고, 시계열 그래프를 그려보겠습니다. 

 

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt 

## getting drug sales dataset 
file_path = 'https://raw.githubusercontent.com/selva86/datasets/master/a10.csv' 
df = pd.read_csv(file_path, 
                 parse_dates=['date'], 
                 index_col='date') 

df.head(12)

#               value
# date	
# 	            3.526591
# 1991-08-01	3.180891
# 1991-09-01	3.252221
# 1991-10-01	3.611003
# 1991-11-01	3.565869
# 1991-12-01	4.306371
# 1992-01-01	5.088335
# 1992-02-01	2.814520
# 1992-03-01	2.985811
# 1992-04-01	3.204780
# 1992-05-01	3.127578
# 1992-06-01	3.270523

## time series plot 
df.plot(figsize=[12, 8]) 
plt.title('Non-Stationary Process: Increasing Variance + Trend + Seasonality', 
          fontsize=16) 
plt.show()

non-stationary process

 

위의 시계열 그래프에서 볼 수 있는 것처럼, (a) 분산이 시간의 흐름에 따라 증가 하고 (분산이 고정이 아님), (b) 추세(trend)가 있으며, (c) 1년 주기의 계절성(seasonality)이 있으므로, 비정상확률과정(non-stationary process)입니다. 

 

KPSS 검정을 통해서 확인해봐도 p-value가 0.01 이므로 유의수준 5% 하에서 귀무가설 (H0: 정상 시계열이다)을 기각하고, 대립가설(H1: 정상 시계열이 아니다)을 채택합니다. 

 

## UDF for KPSS test 
from statsmodels.tsa.stattools import kpss 
import pandas as pd 

def kpss_test(timeseries): 
    print("Results of KPSS Test:") 
    kpsstest = kpss(timeseries, regression="c", nlags="auto") 
    kpss_output = pd.Series(
        kpsstest[0:3], index=["Test Statistic", "p-value", "Lags Used"] ) 
    
    for key, value in kpsstest[3].items(): 
        kpss_output["Critical Value (%s)" % key] = value 
    print(kpss_output)
    
    
## 귀무가설 (H0): 정상 시계열이다
## 대립가설 (H1): 정상 시계열이 아니다 <-- p-value 0.01

kpss_test(df)

# Results of KPSS Test:
# Test Statistic           2.013126
# p-value                  0.010000
# Lags Used                9.000000
# Critical Value (10%)     0.347000
# Critical Value (5%)      0.463000
# Critical Value (2.5%)    0.574000
# Critical Value (1%)      0.739000
# dtype: float64

 

 

 

(1) 분산이 고정적이지 않은 경우 분산 안정화 변환 (variance stabilizing transformation, VST)

 

분산이 고정적이지 않은 경우 멱 변환(Power Transformation)을 통해서 분산을 안정화(variance stabilization) 시켜줍니다. 분산이 고정적이지 않고 추세가 있는 경우 분산 안정화를 추세 제거보다 먼저 해줍니다. 왜냐하면 추세를 제거하기 위해 차분(differencing)을 해줄 때 음수(-)가 생길 수 있기 때문입니다. 

 

power transformation

 

 

원래의 시계열 데이터의 분산 형태에 따라서 적합한 멱 변환(power transformation)을 선택해서 정상확률과정으로 변환해줄 수 있습니다. 아래의 예제 시도표를 참고하세요. 

variance stabilizing transformation (power transfortion)

 

이번 포스팅에서 사용하는 예제는 시간이 흐릴수록 분산이 점점 커지는 형태를 띠고 있으므로 로그 변환(log transformation) 이나 제곱근 변환 (root transformation) 을 해주면 정상 시계열로 변환이 되겠네요. 아래 코드에서는 자연로그를 취해서 로그 변환을 해주었습니다. 

 

## Variance Stabilizing Transformation (VST) by Taking Logarithm
df_vst = np.log(df.value)

df_vst.head()

# date
# 1991-07-01    1.260332
# 1991-08-01    1.157161
# 1991-09-01    1.179338
# 1991-10-01    1.283986
# 1991-11-01    1.271408
# Name: value, dtype: float64


## plotting
df_vst.plot(figsize=(12, 8))
plt.title("Variance Stabilizing Transformation by taking Logarithm", 
          fontsize=16)
plt.show()

variance stabilizing transformation (VST)

 

 

위의 시도표를 보면 시간이 경과해도 분산이 안정화되었음을 알 수 있습니다.  KPSS 검정을 한번 더 해주면 아직 추세(trend)와 계절성(seasonality)가 남아있으므로 여전히 비정상확률과정을 따른다고 나옵니다. 

 

## 귀무가설 (H0): 정상 시계열이다  
## 대립가설 (H1): 정상 시계열이 아니다  <-- p-value 0.01 

kpss_test(df_vst)

# Results of KPSS Test:
# Test Statistic           2.118189
# p-value                  0.010000
# Lags Used                9.000000
# Critical Value (10%)     0.347000
# Critical Value (5%)      0.463000
# Critical Value (2.5%)    0.574000
# Critical Value (1%)      0.739000
# dtype: float64

 

 

 

(2) 추세가 있는 경우 차분을 통한 추세 제거 (de-trend by differencing)

 

차분(differencing)은 현재의 시계열 값에서 시차 t 만큼의 이전 값을 빼주는 것입니다. 

 

1차 차분 = Delta1_Z(t) = Z(t) - Z(t-1) 

2차 차분 = Delta2_Z(t) = Z(t) - Z(t-1) - (Z(t-1) - Z(t-2)) =  Z(t) - 2Z(t-1) + Z(t-2)

 

Python의 diff() 메소드를 사용해서 차분을 해줄 수 있습니다. 이때 차분의 차수 만큼 결측값이 생기는 데요, dropna() 메소드를 사용해서 결측값은 제거해주었습니다. 

 

## De-trend by Differencing
df_vst_diff1 = df_vst.diff(1).dropna()

df_vst_diff1.plot(figsize=(12, 8))
plt.title("De-trend by 1st order Differencing", fontsize=16)
plt.show()

de-trending by 1st order differencing

 

위의 시도표를 보면 1차 차분(1st order differencing)을 통해서 이제 추세(trend)도 제거되었음을 알 수 있습니다. 하지만 아직 계절성(seasonality)이 남아있어서 정상성 조건은 만족하지 않겠네요. 그런데 아래에 KPSS 검정을 해보니 p-value가 0.10 으로서 유의수준 5% 하에서 정상성을 만족한다고 나왔네요. ^^;

 

## 귀무가설 (H0): 정상 시계열이다  <-- p-value 0.10
## 대립가설 (H1): 정상 시계열이 아니다 

kpss_test(df_vst_diff1)

# Results of KPSS Test:
# Test Statistic            0.121364
# p-value                   0.100000
# Lags Used                37.000000
# Critical Value (10%)      0.347000
# Critical Value (5%)       0.463000
# Critical Value (2.5%)     0.574000
# Critical Value (1%)       0.739000
# dtype: float64

 

 

 

(3) 계절성이 있는 경우 계절 차분을 통한 계절성 제거 (de-seasonality by seaanl differencing)

 

아직 남아있는 계절성(seasonality)을 계절 차분(seasonal differencing)을 사용해서 제거해보겠습니다. 1년 12개월 주기의 계절성을 띠고 있으므로 diff(12) 함수로 계절 차분을 실시하고, 12개의 결측값이 생기는데요 dropna() 로 결측값은 제거해주었습니다. 

 

## Stationary Process: De-seasonality by Seasonal Differencing
df_vst_diff1_diff12 = df_vst_diff1.diff(12).dropna()

## plotting
df_vst_diff1_diff12.plot(figsize=(12, 8))
plt.title("De-seasonality by Seasonal Differencing", 
          fontsize=16)
plt.show()

de-seasonality by seasonal differencing

 

위의 시도표를 보면 이제 계절성도 제거가 되어서 정상 시계열처럼 보이네요. 아래에 KPSS 검정을 해보니 p-value 가 0.10 으로서, 유의수준 5% 하에서 귀무가설(H0: 정상 시계열이다)을 채택할 수 있겠네요.

 

## 귀무가설 (H0): 정상 시계열이다  <-- p-value 0.10
## 대립가설 (H1): 정상 시계열이 아니다 

kpss_test(df_vst_diff1_diff12)

# Results of KPSS Test:
# Test Statistic           0.08535
# p-value                  0.10000
# Lags Used                8.00000
# Critical Value (10%)     0.34700
# Critical Value (5%)      0.46300
# Critical Value (2.5%)    0.57400
# Critical Value (1%)      0.73900
# dtype: float64

 

이제 비정상 시계열(non-stationary process)이었던 원래 데이터를 (1) log transformation을 통한 분산 안정화, (2) 차분(differencing)을 통한 추세 제거, (3) 계절 차분(seasonal differencing)을 통한 계절성 제거를 모두 마쳐서 정상 시계열(stationary process) 로 변환을 마쳤으므로, ARIMA 통계 모형을 적합할 수 있게 되었습니다. 

 

이번 포스팅이 많은 도움이 되었기를 바랍니다. 

행복한 데이터 과학자 되세요!  :-)

 

728x90
반응형
Posted by Rfriend
,

이번 포스팅에서는 R의 data.table 패키지의 dcast() 함수를 사용해서 문자열(string)을 대상으로 데이터를 재구조화할 때 집계 함수 (aggregation function) 로서 

 

  (1) 문자열 원소의 개수 (length)

  (2) 문자열을 콤마로 구분해서 붙여쓰기

  (3) 첫번째 문자열만 가져오기

 

하는 방법을 소개하겠습니다. 

 

 

 

 

먼저 간단한 예제 data.table을 만들어보겠습니다. 

 

##---------------------------------
## R data.table dcast() for string
##---------------------------------

#install.packages("data.table")
library(data.table)

x1 <- c('g1', 'g1', 'g1', 'g1', 'g2', 'g2', 'g2', 'g2')
x2 <- c('v1', 'v2', 'v3', 'v3', 'v1', 'v2', 'v2', 'v3')
x3 <- c('a', 'b', 'c', 'd', 'e', 'f', 'g', 'h')

dt <- data.table(x1, x2, x3)
print(dt)
#    x1  x2  x3
# 1: g1  v1   a
# 2: g1  v2   b
# 3: g1  v3   c
# 4: g1  v3   d
# 5: g2  v1   e
# 6: g2  v2   f
# 7: g2  v2   g
# 8: g2  v3   h

 

 

R 의 data.table 패키지의 dcast() 함수를 사용해 데이터를 재구조화하면서 문자열을 대상으로 집계(value.var)를 할 때 집계 함수 (aggregation function) 을 명시적으로 적어주지 않으면 아래와 같은 경고 메시지가 발생합니다. 

 

Warning message: Aggregate function missing, defaulting to 'length'

 

이것은 문자열을 대상으로는 합계(sum), 최소값(min), 최대값(max), 평균(mean) 등의 숫자형을 대상으로 하는 요약통계량을 사용할 수 없기 때문입니다. 

 

##-- warning message
##: Aggregate function missing, defaulting to 'length'
dcast(dt, x1 ~ x2, 
      value.var = "x3")

# Aggregate function missing, defaulting to 'length'
#    x1  v1  v2  v3
# 1: g1   1   1   2
# 2: g2   1   2   1

 

 

따라서 dcast() 함수로 데이터를 재고조화시 문자열을 대상으로 집계를 한다면

  (1) 문자열 원소의 개수 (length)

  (2) 문자열을 콤마로 구분해서 붙여쓰기

  (3) 첫번째 문자열만 가져오기

와 같이 문자열에 맞는 집계함수를 지정해주어야 합니다. 

 

 

 

(1) dcast() 함수로 데이터셋 재구조화 시 문자열을 원소의 개수 (length) 로 집계

 

문자열 대상 집계일 때는 default 설정이 원소의 개수 (length) 이므로 위와 결과는 동일합니다만, 이번에는 경고 메시지가 안떴습니다. 

 

##-- (1) counting the number of values as an aggregation function for string values
dcast(dt, x1 ~ x2, 
      fun.aggregate = length, 
      value.var = "x3")

#    x1  1  2  3
# 1: g1  1  1  2
# 2: g2  1  2  1

 

 

 

(2) dcast() 함수로 데이터셋 재구조화 시 문자열을 콤마로 구분해서 붙여쓰기

 

dcast() 로 재구조화 시 하나의 셀 안에 여러개의 원소가 존재하게 될 경우, 이들 문자열 원소들을 콤마로 구분해서 옆으로 나란히 붙여서 집계하는 사용자 정의 함수를 fun.aggregate 매개변수란에 써주었습니다. 

 

##-- (2) concatenation as an aggregation function for string values
dcast(dt, x1 ~ x2, 
      fun.aggregate = function(x) if (length(x)==1L) x else paste(x, collapse=","), 
      value.var = "x3")

#    x1  1     2     3
# 1: g1  a     b   c,d
# 2: g2  e   f,g     h

 

 

 

(3) dcast() 함수로 데이터셋 재구조화 시 첫번째 문자열만 가져오기

 

dcast() 로 재구조화 시 하나의 셀 안에 여러개의 원소가 존재하게 될 경우, 이들 복수개의 원소들 중에서 첫번째 원소만 가져오는 사용자정의함수를 fun.aggregate 매개변수란에 작성해주었습니다. 

 

##-- (3) keeping the first value as an aggregation function for string values
dcast(dt, x1 ~ x2, 
      fun.aggregate = function(x) if (length(x)==1L) x else x[1], 
      value.var = "x3")
      
#    x1  1  2  3
# 1: g1  a  b  c
# 2: g2  e  f  h

 

 

이번 포스팅이 많은 도움이 되었기를 바랍니다. 

행복한 데이터 과학자 되세요!  :-)

 

728x90
반응형
Posted by Rfriend
,

주식을 하는 분들은 아마도 대표적인 시계열 데이터인 주가의 이동평균, 누적평균 그래프에 이미 익숙할 것입니다. 

 

이번 포스팅에서는 R의 zoo 패키지의 rollapply() 라는 window function

 

(1) Rolling Windows 를 사용해서 시계열 데이터의 이동 평균 구하기

     (average of time series using rolling windows)

(2) Expanding Windows 를 사용해서 시계열 데이터의 누적 평균 구하기

     (average of time series using expanding windows)

 

방법을 소개하겠습니다. 

 

 

[ 이동 평균 (average using Rolling Windows) vs. 누적 평균 (average using Expanding Windows) ]

moving average (rolling windows) vs. cumulative average (expanding windows) using R zoo rollapply() function

 

시계열 데이터를 전처리하고 분석할 때 Window Function 을 자주 사용하는데요, 

 - Rolling Windows : 특정 window width (예: 10분, 1시간, 1일 등) 를 유지한채 측정 단위시간별로 이동하면서 분석 

 - Expanding Windows : 처음 시작 시점은 고정한 채, 시간이 흐름에 따라 신규로 포함되는 데이터까지 누적해서 분석

하는 차이가 있습니다. 바로 위에 Rolling Windows 와 Expanding Windows 를 도식화 해놓은 자료를 보면 금방 이해가 될거예요. 

 

만약 시계열 데이터에 추세(trend) 나 계절성 (seasonality) 이 있다면 Rolling Windows 가 적당하며, 시계열 데이터에 추세나 계절성이 없이 안정적(stable) 이다면 Expanding Windows 를 사용해서 더 많은 데이터를 이용해서 요약 통계량을 계산하는게 유리할 수 있겠습니다. 

 

시계열 예측 모델링할 때는 Rolling Windows 를 사용해서 모델 성능을 검증합니다. 

 

 

R 의 zoo 패키지의 rollapply() 함수를 사용할 것이므로, zoo 패키지를 먼저 설치하고 임포팅합니다. 

그리고 예제로 사용할 간단한 시계열 데이터를 만들어보겠습니다. 추세와 노이즈가 있는 시계열 데이터 입니다. 

 

## ------------
## Wimdow functions in Time Series
## (1) Rolling window
## (2) Expanding window
## R zoo's rollapply(): https://www.rdocumentation.org/packages/zoo/versions/1.8-9/topics/rollapply
## ------------

install.packages("zoo")
library(zoo)

## generating a time series with trend and noise
set.seed(1) # for reproducibility
x <- rnorm(n=100, mean=0, sd=10) + 1:100

plot(x, type='l', 
     main="time series plot with trend and noise")

time series with trend and noise

 

 

 

(1) Rolling Windows 를 사용해서 시계열 데이터의 이동 평균 구하기

     (average of time series using rolling windows)

 

zoo 패키지의 rollapply() 함수에서

- width 매개변수는 'window width' 를 설정할 때 사용합니다.

- FUN 매개변수에는 원하는 함수를 지정해줄 수 있으므로 매우 강력하고 유연하게 사용할 수 있습니다. 아래 예에서는 평균(mean)과 최대값(max) 을 계산하는 함수를 사용해보았습니다.

- align 은 데이터의 기준을 정렬할 때 왼쪽("left"), 중앙("centered", default 설정), 오른쪽("right") 중에서 지정할 수 있습니다. 이때 align="left"로 설정해주면 자칫 잘못하면 미래의 데이터를 가져다가 요약 통계량을 만드는 실수 (lookahead) 를 할 수도 있으므로, 만약 예측 모델링이 목적이라면 lookahead 를 하는건 아닌지 유의해야 합니다. 

- partial=TRUE 로 설정하면 양쪽 끝부분에 window width 의 개수에 데이터 포인트 개수가 모자라더라도 있는 데이터만 가지고 부분적으로라도 함수의 통계량을 계산해줍니다. 

 

## (1) Rolling Windows

## (1-1) moving average
f_avg_rolling_win <- rollapply(
  data=zoo(x), 
  width=10, # window width
  FUN=function(w) mean(w), 
  # 'align' specifies whether the index of the result should be left-aligned 
  # or right-aligned or centered (default) 
  # compared to the rolling window of observations. 
  align="right", 
  # If 'partial=TRUE', then the subset of indexes 
  # that are in range are passed to FUN.
  partial=TRUE)

## (1-2) moving max
f_max_rolling_win <- rollapply(
  zoo(x), 
  10, 
  function(w) max(w), 
  align="right", 
  partial=TRUE)

plot(x, col="gray", lwd=1, type="l", main="Average and Max using Rolling Window")
lines(f_avg_rolling_win, col="blue", lwd=2, lty="dotted")
lines(f_max_rolling_win, col="red", lwd=2, lty="dashed")
legend("topleft", 
       c("Average with Rolling Windows", "Max with Rolling Windows"), 
       col = c("blue", "red"), 
       lty = c("dotted", "dashed"))

moving average and max using the rolling windows

 

 

 

 

(2) Expanding Windows 를 사용해서 시계열 데이터의 누적 평균 구하기

     (average of time series using expanding windows)

 

R 에서 zoo 패키지의 rollapply() 함수로 Expanding Windwos 를 사용하려면 width = seq_along(x) 를 지정해주면 누적으로 함수를 계산해줍니다. 

 

아래 예에서는 누적으로 평균과 최대값을 계산해서 시각화 한건데요, 우상향 하는 추세가 있는 시계열이다보니 누적으로 평균을 구하면 시계열 초반의 낮은 값들까지 모두 포함이 되어서 누적평균 값이 최근 값들을 제대로 따라가지 못하고 있습니다.

반면, 누적으로 최대값을 계산한 값은 중간에 소폭 값이 줄어들더라도 계산 시점까지 누적으로 최대값을 계산하므로, 항상 우상향하는 누적 최대값을 보여주고 있습니다.

(위의 (1)번의 이동평균, 이동최대값과 (2) 누적평균, 누적최대값을 비교해서 보세요.)

 

# (2) Expanding Windows

## (2-1) cumulative average
f_avg_expanding_win <- rollapply(
  data=zoo(x), 
  width=seq_along(x), # expanding windows
  FUN=function(w) mean(w), # average
  align="right", 
  partial=TRUE)

## (2-2) cumulative max
f_max_expanding_win <- rollapply(
  zoo(x), 
  seq_along(x), # expanding windows
  function(w) max(w), # max
  align="right", 
  partial=TRUE)

## plotting
plot(x, col="gray", lwd=1, type="l", main="Average and Max using Expanding Window")
lines(f_avg_expanding_win, col="blue", lwd=2, lty="dotted")
lines(f_max_expanding_win, col="red", lwd=2, lty="dashed")
legend("topleft", 
       c("Average with Expanding Windows", "Max with Expanding Windows"), 
       col = c("blue", "red"), 
       lty = c("dotted", "dashed"))

 

 

[ Reference ]

- R zoo's rollapply(): https://www.rdocumentation.org/packages/zoo/versions/1.8-9/topics/rollapply

 

 

이번 포스팅이 많은 도움이 되었기를 바랍니다. 

행복한 데이터 과학자 되세요!  :-)

 

728x90
반응형
Posted by Rfriend
,