Hello,. tulisan kali ini sedikit akan bermain dengan data mining, yaitu terkait bagaimana cara kita mengekstrak informasi dari suatu data set dan mentransformasikannya ke dalam bentuk yang mudah dipahami untuk keperluan analisa lebih lanjut. Bayangkan saja jika kita mempunyai ribuan atau mungkin jutaan data, maka menampilkan hanya data statistik beserta distribusinya akan menjadi tidak cukup untuk menjelaskan perilaku (behaviour) dari suatu data set. Sehingga menemukan pola (pattern) atau sejumlah data yang bersifat paling informatif dari sekumpulan data set akan menjadi hal yang saya kira sangat penting untuk dipelajari.
Seperti halnya pada tulisan-tulisan sebelumnya (seperti PSO, ICA, PCA dll), tulisan kali ini lebih menekankan pada bagaimana cara kita menggunakan suatu persamaan ketimbang cara mendapatkan atau menurunkannya (ini kerjaan orang matematik saya kira hehehe). Sedangkan, paper yang akan dibahas pada tulisan kali ini berjudul "Optimal Multi-scale Patterns in Time Series Streams" yang dipublikasi oleh Papadimitriou pada tahun 2006. Ini memang paper yang cukup lama, sekitar 10 tahun yang lalu, namun masih sangat menarik untuk dipelajari dan masih mungkin untuk dikembangkan. Data yang digunakan dalam tulisan ini adalah data pengukuran CO$_2$ di kota Bandung pada bulan Februari 2011. Untuk lebih jelasnya, silahkan baca tulisan ini. Selamat membaca. (program dapat didownload di akhir tulisan ini).
 |
| Konsentrasi CO$_2$ bulan Februari 2011 |
Menemukan pola dengan membagi data berdasarkan panjang window adalah masalah umum dan ada banyak sekali metode yang diusulkan untuk memecahkannya [1-4]. Jika panjang window dapat diestimasi secara tepat, pola dalam suatu data sekuensial akan dengan mudah diekstrak. Analisa sekuensial berdasarkan panjang window sering kali digunakan seperti menemukan time series motif [5]. analisa electroencephalograms (EEGs) [6], dan lain sebagainya. Sedangkan dalam tulisan kali ini, kita akan mencoba mengekstrak pola untuk data konsentrasi CO$_2$ yang diukur selama satu bulan menggunakan metode yang diusulkan oleh Papadimitriou [1].
Problem
Seperti yang terlihat pada gambar di bawah ini, tujuan dari tulisan ini adalah untuk mendapatkan panjang window yang tepat sehingga pola dari data konsentrasi CO$_2$ dapat diektrak dan dianalisa. Alasan penggunakan data CO$_2$ dalam contoh ini adalah karena periode dari data CO$_2$ sangat mudah ditebak yaitu 24 jam. Dengan mengetahui terlebih dahulu informasi dari periode CO$_2$, kita akan dengan mudah mengetahui dan memastikan apakah algoritma yang sedang kita pelajari ini sudah bekerja dengan tepat atau belum.
Locally Optimal Patterns
Jika anda membaca dengan teliti pada paper [1] di link ini http://www.cs.cmu.edu/afs/cs.cmu.edu/Web/People/spapadim/pdf/msbasis_sigmod06.pdf, anda mungkin akan sadar bahwa metode yang sedang kita pelajari ini hanya sebagian kecil dari yang diusulkan oleh paper tersebut. Karena tidak semua bagian dari paper itu saya pahami dengan baik, hanya sebagian kecil saja yaitu bab 3 dan 4 saja. Jadi tulisan kali ini akan fokus terutama pada bagaimana menemukan panjang window yang tepat menggunakan metode locally optimal patterns yang dijelaskan di bab 4 paper tersebut.
- Example Problem
Pada paper [1] tersebut, contoh masalah yang digunakan sebagai ilustrasi bagaimana metode yang diusulkan bekerja yaitu gelombang sinus dengan periode 50 dan berisi Gaussian noise dengan rataan 5 dan varian 0.5 dimana secara matematik dapat dituliskan sebagai berikut:
Dengan menggunakan Octave, persamaan di atas dapat dituliskan menjadi
close all
clear
clc
fs = 1;
t = 0:1/fs:2000;
x = sin(2*pi*t/50) + randn(size(t))/sqrt(2) + 5;
figure
plot(t,x)
xlim([0 max(t)])
grid on
xlabel('waktu (detik)')
ylim([0 9])
Hasilnya akan seperti di bawah ini.
- Algorithm
Lebih rinci terkait algoritma dapat dilihat di bawah ini,
----------------------------------------------------------------------------------------------------------
Algoritma LocalPattern $(\textbf{x} \in \mathbb{R}^m, \ \omega, \ k)$
----------------------------------------------------------------------------------------------------------
Gunakan time-delay coordinate matrix. $\textbf{X}^{(\omega)} \gets \textrm{delay}(\textbf{x}, \omega)$
Hitung SVD (Singular Value Decomposition). $\textbf{X}^{(\omega)}= \textbf{U}^{(\omega)} \ \Sigma^{(\omega)} \textbf{V}^{(\omega)}$
Dapatkan Singular Values-nya. $\Sigma^{(\omega)}=\textrm{diag} ([\sigma_1^{(\omega)} \ \sigma_2^{(\omega)} \ ... \ \sigma_r^{(\omega)}])$
Kemudian hitung power profile. $\pi^{(\omega)} = \sum_{i=k+1}^{\omega} (\sigma_i^{(\omega)})^2/\omega$
----------------------------------------------------------------------------------------------------------
Algoritma time-delay coordinate matrix dimana $\textbf{x} ≡ [x_1, \ x_2, \ ... \, x_m]^T$ yakni sebagai berikut
----------------------------------------------------------------------------------------------------------
Prosedur delay$(\textbf{x} \in \mathbb{R}^{m}, \ \omega)$
----------------------------------------------------------------------------------------------------------
$m' \gets \lceil m/\omega\rceil$ dan $n' \gets \omega$
Output adalah $\textbf{X}^{(\omega)} \in \mathbb{R}^{m' \times n'}$
$\textrm{for} \ t=1$ to $m' \ \textrm{do}$
$\textrm{Row} \ \textbf{X}_{(t)}^{(\omega)} \gets [x_{(t-1)n'+1}, \ x_{(t-1)n'+2}, \ ... \ , x_{tn'}]$
$\textrm{end for}$
----------------------------------------------------------------------------------------------------------
Dengan menggunakan Octave, algoritma di atas dapat dituliskan seperti berikut ini.
close all
clear
clc
fs = 1;
t = 0:1/fs:2000;
x = sin(2*pi*t/50) + randn(size(t))/sqrt(2) + 5;
figure
plot(t,x)
xlim([0 max(t)])
grid on
xlabel('waktu (detik)')
ylim([0 9])
% definisikan rentang panjang window yang ingin dicari
W = 1:400;
K = 1; % default K = 1, tapi bisa juga K = 1:2
power_profile = zeros(length(K),length(W));
for i = 1:length(W)
% time-delay coordinate matrix
m = length(x);
w = round(fs*W(i));;
m_abs = (m-mod(m,w))/w;
X = reshape(x(1:end-mod(m,w))/norm(x(1:end-mod(m,w))),w,m_abs)';
% local pattern
for k = K
[U S V] = svd(X,0);
power_profile(k,i) = sum(diag(S(k+1:end,k+1:end)).^2)/w;
end
fprintf(stderr,[num2str(i) '\n'])
end
figure
plot(repmat(W,size(power_profile,1),1)',power_profile','LineWidth',1.1)
xlim([26 400])
ylim([0 max(power_profile(:,max(find(W<=26))))*1.05])
xlabel('Panjang Window')
title('Power Profile')
Hasilnya akan seperti di bawah ini.
- Choosing the window
Jika kita membaca paper [1] dan memperhatikan gambar 5 di halaman 5. Kita akan menyadari bahwa gambar yang kita peroleh dari program di atas memiliki bentuk atau hasil yang mirip sekali dengan gambar 5 pada paper dengan sedikit perbedaan pada yaxis. Namun secara umum hasil yang kita dapat sudah tepat. Jika kita perlebar rentang dari xaxis menjadi antara 15 dan 400 seperti terlihat pada gambar di bawah ini. Ternyata, ada penurunan tajam yang lain tepatnya pada 25 yang mungkin sengaja tidak ditampilkan pada paper. Tapi ini tidak berarti ada kesalahan pada paper tersebut, namun sebaliknya paper tersebut menunjukkan bagaimana cara yang benar untuk memilih panjang window dengan tepat untuk kasus tersebut.
Rincian bagaimana memilih panjang window yang tepat dijelaskan seperti pada berikut ini:
- Hitung power profile sebagai fungsi panjang window.
- Temukan penurunan tajam yang memenuhi syarat $\omega \approx i\omega_{\circ}$. Cara termudah adalah dengan memulai memilih window terpendek.
- Namun, jika tidak ada penurunan tajam, itu berarti tidak ada komponen periodis dalam data tersebut.
- (tambahan aturan) gunakan nilai $k$ lebih dari satu angka misal $k = [1 \ 2]$ (lihat font berwarna merah pada kode Octave di atas) sebagai informasi tambahan saat kita memilih panjang window seperti pada gambar di bawah ini.
Dari gambar di atas, sebaiknya kita memilih window yang memiliki penurunan tajam ganda (penurunan tajam terjadi pada posisi yang sama untuk warna biru dan hijau). Maka nilai $\omega$ adalah 50, 100, 150, 200, 250, 300, 350 dan 400. Sehingga bisa disimpulkan bahwa panjang window terbaik adalah 50.
- Pattern Extraction
Untuk mengekstrak pola, cukup tambahkan kode berikut ini pada kode sebelumnya.
% Ekstrak data berdasarkan informasi panjang window terbaik
wbest = 50;
m = length(x);
m_abs = (m-mod(m,wbest))/wbest;
X = reshape(x(1:end-mod(m,wbest))/norm(x(1:end-mod(m,wbest))),wbest,m_abs)';
Xbar = mean(X);
mini = min(min(X));
maxi = max(max(X));
figure
subplot(121)
plot(t(1:length(Xbar)),Xbar,'LineWidth',1.1)
xlim([1 wbest])
axis square
ylim([mini maxi])
title('Rataan Pola dari X')
subplot(122)
plot(t(1:length(Xbar)),X(1,:),'LineWidth',1.1)
xlim([1 wbest])
axis square
ylim([mini maxi])
title('Salah satu Pola dalam X')
Hasilnya dapat dilihat pada gambar di bawah ini.
Data CO$_2$
Jika kita menggunakan data CO$_2$, kita akan mendapatkan hasil seperti di bawah ini. Terlihat dengan jelas bahwa periode untuk data CO$_2$ adalah 24. (silahkan dijalankan sendiri program ini untuk memastikan hasilnya)
Download
Program di atas dapat didownload di link ini. Selamat mencoba.
Problem
Seperti yang terlihat pada gambar di bawah ini, tujuan dari tulisan ini adalah untuk mendapatkan panjang window yang tepat sehingga pola dari data konsentrasi CO$_2$ dapat diektrak dan dianalisa. Alasan penggunakan data CO$_2$ dalam contoh ini adalah karena periode dari data CO$_2$ sangat mudah ditebak yaitu 24 jam. Dengan mengetahui terlebih dahulu informasi dari periode CO$_2$, kita akan dengan mudah mengetahui dan memastikan apakah algoritma yang sedang kita pelajari ini sudah bekerja dengan tepat atau belum.
 | |||
| Panjang Window yang tepat untuk data CO$_2$ adalah 24 jam |
Jika anda membaca dengan teliti pada paper [1] di link ini http://www.cs.cmu.edu/afs/cs.cmu.edu/Web/People/spapadim/pdf/msbasis_sigmod06.pdf, anda mungkin akan sadar bahwa metode yang sedang kita pelajari ini hanya sebagian kecil dari yang diusulkan oleh paper tersebut. Karena tidak semua bagian dari paper itu saya pahami dengan baik, hanya sebagian kecil saja yaitu bab 3 dan 4 saja. Jadi tulisan kali ini akan fokus terutama pada bagaimana menemukan panjang window yang tepat menggunakan metode locally optimal patterns yang dijelaskan di bab 4 paper tersebut.
- Example Problem
Pada paper [1] tersebut, contoh masalah yang digunakan sebagai ilustrasi bagaimana metode yang diusulkan bekerja yaitu gelombang sinus dengan periode 50 dan berisi Gaussian noise dengan rataan 5 dan varian 0.5 dimana secara matematik dapat dituliskan sebagai berikut:
$x_t = \displaystyle \sin \Big(\frac{2 \ \pi \ t}{50}\Big) + \mathcal{N}(5,0.5)$
Dengan menggunakan Octave, persamaan di atas dapat dituliskan menjadi
close all
clear
clc
fs = 1;
t = 0:1/fs:2000;
x = sin(2*pi*t/50) + randn(size(t))/sqrt(2) + 5;
figure
plot(t,x)
xlim([0 max(t)])
grid on
xlabel('waktu (detik)')
ylim([0 9])
Hasilnya akan seperti di bawah ini.
 |
| $x_t = \sin (2 \pi t/50) + \mathcal{N}(5,0.5)$ |
Lebih rinci terkait algoritma dapat dilihat di bawah ini,
----------------------------------------------------------------------------------------------------------
Algoritma LocalPattern $(\textbf{x} \in \mathbb{R}^m, \ \omega, \ k)$
----------------------------------------------------------------------------------------------------------
Gunakan time-delay coordinate matrix. $\textbf{X}^{(\omega)} \gets \textrm{delay}(\textbf{x}, \omega)$
Hitung SVD (Singular Value Decomposition). $\textbf{X}^{(\omega)}= \textbf{U}^{(\omega)} \ \Sigma^{(\omega)} \textbf{V}^{(\omega)}$
Dapatkan Singular Values-nya. $\Sigma^{(\omega)}=\textrm{diag} ([\sigma_1^{(\omega)} \ \sigma_2^{(\omega)} \ ... \ \sigma_r^{(\omega)}])$
Kemudian hitung power profile. $\pi^{(\omega)} = \sum_{i=k+1}^{\omega} (\sigma_i^{(\omega)})^2/\omega$
----------------------------------------------------------------------------------------------------------
Algoritma time-delay coordinate matrix dimana $\textbf{x} ≡ [x_1, \ x_2, \ ... \, x_m]^T$ yakni sebagai berikut
----------------------------------------------------------------------------------------------------------
Prosedur delay$(\textbf{x} \in \mathbb{R}^{m}, \ \omega)$
----------------------------------------------------------------------------------------------------------
$m' \gets \lceil m/\omega\rceil$ dan $n' \gets \omega$
Output adalah $\textbf{X}^{(\omega)} \in \mathbb{R}^{m' \times n'}$
$\textrm{for} \ t=1$ to $m' \ \textrm{do}$
$\textrm{Row} \ \textbf{X}_{(t)}^{(\omega)} \gets [x_{(t-1)n'+1}, \ x_{(t-1)n'+2}, \ ... \ , x_{tn'}]$
$\textrm{end for}$
----------------------------------------------------------------------------------------------------------
Dengan menggunakan Octave, algoritma di atas dapat dituliskan seperti berikut ini.
close all
clear
clc
fs = 1;
t = 0:1/fs:2000;
x = sin(2*pi*t/50) + randn(size(t))/sqrt(2) + 5;
figure
plot(t,x)
xlim([0 max(t)])
grid on
xlabel('waktu (detik)')
ylim([0 9])
% definisikan rentang panjang window yang ingin dicari
W = 1:400;
K = 1; % default K = 1, tapi bisa juga K = 1:2
power_profile = zeros(length(K),length(W));
for i = 1:length(W)
% time-delay coordinate matrix
m = length(x);
w = round(fs*W(i));;
m_abs = (m-mod(m,w))/w;
X = reshape(x(1:end-mod(m,w))/norm(x(1:end-mod(m,w))),w,m_abs)';
% local pattern
for k = K
[U S V] = svd(X,0);
power_profile(k,i) = sum(diag(S(k+1:end,k+1:end)).^2)/w;
end
fprintf(stderr,[num2str(i) '\n'])
end
figure
plot(repmat(W,size(power_profile,1),1)',power_profile','LineWidth',1.1)
xlim([26 400])
ylim([0 max(power_profile(:,max(find(W<=26))))*1.05])
xlabel('Panjang Window')
title('Power Profile')
Hasilnya akan seperti di bawah ini.
 |
| Power Profile dengan xaxis antara 26 - 400 |
- Choosing the window
Jika kita membaca paper [1] dan memperhatikan gambar 5 di halaman 5. Kita akan menyadari bahwa gambar yang kita peroleh dari program di atas memiliki bentuk atau hasil yang mirip sekali dengan gambar 5 pada paper dengan sedikit perbedaan pada yaxis. Namun secara umum hasil yang kita dapat sudah tepat. Jika kita perlebar rentang dari xaxis menjadi antara 15 dan 400 seperti terlihat pada gambar di bawah ini. Ternyata, ada penurunan tajam yang lain tepatnya pada 25 yang mungkin sengaja tidak ditampilkan pada paper. Tapi ini tidak berarti ada kesalahan pada paper tersebut, namun sebaliknya paper tersebut menunjukkan bagaimana cara yang benar untuk memilih panjang window dengan tepat untuk kasus tersebut.
 |
| Power Profile dengan xaxis antara 15 - 400 |
Rincian bagaimana memilih panjang window yang tepat dijelaskan seperti pada berikut ini:
- Hitung power profile sebagai fungsi panjang window.
- Temukan penurunan tajam yang memenuhi syarat $\omega \approx i\omega_{\circ}$. Cara termudah adalah dengan memulai memilih window terpendek.
- Namun, jika tidak ada penurunan tajam, itu berarti tidak ada komponen periodis dalam data tersebut.
- (tambahan aturan) gunakan nilai $k$ lebih dari satu angka misal $k = [1 \ 2]$ (lihat font berwarna merah pada kode Octave di atas) sebagai informasi tambahan saat kita memilih panjang window seperti pada gambar di bawah ini.
 |
| Power Profile menggunakan dua nilai $k$. |
Dari gambar di atas, sebaiknya kita memilih window yang memiliki penurunan tajam ganda (penurunan tajam terjadi pada posisi yang sama untuk warna biru dan hijau). Maka nilai $\omega$ adalah 50, 100, 150, 200, 250, 300, 350 dan 400. Sehingga bisa disimpulkan bahwa panjang window terbaik adalah 50.
- Pattern Extraction
Untuk mengekstrak pola, cukup tambahkan kode berikut ini pada kode sebelumnya.
% Ekstrak data berdasarkan informasi panjang window terbaik
wbest = 50;
m = length(x);
m_abs = (m-mod(m,wbest))/wbest;
X = reshape(x(1:end-mod(m,wbest))/norm(x(1:end-mod(m,wbest))),wbest,m_abs)';
Xbar = mean(X);
mini = min(min(X));
maxi = max(max(X));
figure
subplot(121)
plot(t(1:length(Xbar)),Xbar,'LineWidth',1.1)
xlim([1 wbest])
axis square
ylim([mini maxi])
title('Rataan Pola dari X')
subplot(122)
plot(t(1:length(Xbar)),X(1,:),'LineWidth',1.1)
xlim([1 wbest])
axis square
ylim([mini maxi])
title('Salah satu Pola dalam X')
Hasilnya dapat dilihat pada gambar di bawah ini.
Data CO$_2$
Jika kita menggunakan data CO$_2$, kita akan mendapatkan hasil seperti di bawah ini. Terlihat dengan jelas bahwa periode untuk data CO$_2$ adalah 24. (silahkan dijalankan sendiri program ini untuk memastikan hasilnya)
 |
| Power Profile untuk Data CO$_2$ |
 |
| Hasil ektraksi pola |
Program di atas dapat didownload di link ini. Selamat mencoba.
References
[1] Papadimitriou, Spiros, and Philip Yu. "Optimal multi-scale patterns in time series streams." Proceedings of the 2006 ACM SIGMOD international conference on Management of data. ACM, 2006.
[2] Lijffijt, Jefrey, Panagiotis Papapetrou, and Kai Puolamäki. "Size
matters: choosing the most informative set of window lengths for mining
patterns in event sequences." Data Mining and Knowledge Discovery 29.6
(2015): 1838-1864.
[3] Vlachos, Michail, et al. "Identifying similarities, periodicities and bursts for online search queries." Proceedings of the 2004 ACM SIGMOD international conference on Management of data. ACM, 2004.
[4] Patel, Pranav, et al. "Mining motifs in massive time series databases." Data Mining, 2002. ICDM 2003. Proceedings. 2002 IEEE International Conference on. IEEE, 2002.
[5] Chiu, Bill, Eamonn Keogh, and Stefano Lonardi. "Probabilistic discovery of time series motifs." Proceedings of the ninth ACM SIGKDD international conference on Knowledge discovery and data mining. ACM, 2003.
[6] Sörnmo, Leif, and Pablo Laguna. Bioelectrical signal processing in cardiac and neurological applications. Vol. 8. Academic Press, 2005.
No comments:
Post a Comment