Yup, tulisan kali ini sedikit akan membahas tentang algoritma LMS. Seperti yang Wikipedia katakan bahwa LMS adalah salah satu adaptive filter yang digunakan untuk meniru perilaku sebuah filter (yang tidak kita ketahui) dengan memprediksi koefisien filter tersebut berdasarkan kuadrat terkecil dari sinyal error. Untuk lebih lengkapnya silahkan merujuk ke Wikipedia di link berikut https://en.wikipedia.org/wiki/Least_mean_squares_filter.
Bayangkan bahwa kita mempunyai sinyal input $x(n)$ dan sebuah unknown system $\mathbf{\textrm{h}}(n)$ dimana output dari unknown system ini adalah $d(n)(=\mathbf{\textrm{h}}(n)*x(n))$. Dalam suatu pengujian misalnya, kita hanya mengetahui informasi $x(n)$ dan $d(n)$ saja, padahal kita membutuhkan informasi $h(n)$ juga. Jadi apakah yang harus kita lakukan? jawabanya sederhana, kita cukup mengestimasi nilai $\mathbf{\textrm{h}}(n)$ menggunakan $x(n)$ dan $d(n)$, nah saat itulah algoritma LMS diperlukan.
Algoritma LMS sangatlah mudah dan dapat diringkas seperti di bawah ini
Parameter:
$p$ = order dari filter estimasi $\hat{\mathbf{\textrm{h}}}(n)$
$\mu$ = step size atau ada yang bilang juga learning rate
Inisialisasi:
$\hat{\mathbf{\textrm{h}}}(0)$ = zeros($p$)
Komputasi:
for $n = 0, 1, 2, 3, ...$
$x(n)$ = $[x(n),x(n-1), ..., x(n-p+1)]^T$
$e(n) = d(n) - \hat{\mathbf{\textrm{h}}}(n)^Hx(n)$
$\hat{\mathbf{\textrm{h}}}(n+1) = \hat{\mathbf{\textrm{h}}}(n) + \mu e(n) x(n)$
Parameter:
$p$ = order dari filter estimasi $\hat{\mathbf{\textrm{h}}}(n)$
$\mu$ = step size atau ada yang bilang juga learning rate
Inisialisasi:
$\hat{\mathbf{\textrm{h}}}(0)$ = zeros($p$)
Komputasi:
for $n = 0, 1, 2, 3, ...$
$x(n)$ = $[x(n),x(n-1), ..., x(n-p+1)]^T$
$e(n) = d(n) - \hat{\mathbf{\textrm{h}}}(n)^Hx(n)$
$\hat{\mathbf{\textrm{h}}}(n+1) = \hat{\mathbf{\textrm{h}}}(n) + \mu e(n) x(n)$
Algoritma di atas dapat dituliskan ke dalam beberapa bahasa pemrograman yakni sebagai berikut
Octave
close all
clear
clc
% Assumsikan h(n) adalah unknown system
h = [0.25 1.27 -0.3 0.5 -1.0 0.75 0.75 2.5];
% x(n) adalah input signal
x = rand(1,1000);
% d(n) adalah output sinyal
d = filter(h,1,x);
%% Parameter untuk algoritma LMS
% Order dari filter estimasi h_est(n)
p = 15;
% Step size atau learning rate
mu = 0.1;
% Inisialisasi awal koefisien h_est(0)
h_est = zeros(1,p);
% Inisialisasi state h_est(n) / x(n) = [x(n),x(n-1), ..., x(n-p+1)]
h_stt = zeros(1,p);
% error
e = zeros(size(x));
tic
for n = 1:length(x)
h_stt = [x(n) h_stt(1:end-1)];
e(n) = d(n) - h_est*h_stt.';
h_est = h_est + mu*e(n)*h_stt;
end
toc
figure
subplot(211)
plot(e)
title('error')
subplot(212)
stem(h)
hold on
stem(h_est,'r.')
legend('original coef','estimated coef')
xlim([0 p])
clear
clc
% Assumsikan h(n) adalah unknown system
h = [0.25 1.27 -0.3 0.5 -1.0 0.75 0.75 2.5];
% x(n) adalah input signal
x = rand(1,1000);
% d(n) adalah output sinyal
d = filter(h,1,x);
%% Parameter untuk algoritma LMS
% Order dari filter estimasi h_est(n)
p = 15;
% Step size atau learning rate
mu = 0.1;
% Inisialisasi awal koefisien h_est(0)
h_est = zeros(1,p);
% Inisialisasi state h_est(n) / x(n) = [x(n),x(n-1), ..., x(n-p+1)]
h_stt = zeros(1,p);
% error
e = zeros(size(x));
tic
for n = 1:length(x)
h_stt = [x(n) h_stt(1:end-1)];
e(n) = d(n) - h_est*h_stt.';
h_est = h_est + mu*e(n)*h_stt;
end
toc
figure
subplot(211)
plot(e)
title('error')
subplot(212)
stem(h)
hold on
stem(h_est,'r.')
legend('original coef','estimated coef')
xlim([0 p])
Phyton
# sudo apt-get install python3-matplotlib
# sudo apt-get install python3-pip
# sudo pip3 install scipy
# sudo pip3 install numpy
import matplotlib.pyplot as plt
import scipy.signal as scpsig
import numpy as np
import time as tm
# Assumsikan h(n) adalah unknown system
h = np.array([0.25, 1.27, -0.3, 0.5, -1.0, 0.75, 0.75, 2.5])
# x(n) adalah input signal
x = np.random.rand(1000)
# d(n) adalah output sinyal
d = scpsig.lfilter(h,1,x)
## Parameter untuk algoritma LMS
# Order dari filter estimasi h_est(n)
p = 15
# Step size atau learning rate
mu = 0.1
# Inisialisasi awal koefisien h_est(0)
h_est = np.zeros(p)
# Inisialisasi state h_est(n) / x(n) = [x(n),x(n-1), ..., x(n-p+1)]
h_stt = np.zeros(p)
# error
e = np.zeros(x.size)
tic = tm.time()
for n in range(0,x.size):
h_stt = np.append(x[n], h_stt[:-1])
e[n] = d[n] - np.dot(h_est,h_stt)
h_est = h_est + mu*e[n]*h_stt
toc = tm.time()
print('%.4f seconds' % (toc-tic))
plt.figure()
plt.subplot(211)
plt.plot(range(1,e.size+1),e)
plt.title('error')
plt.subplot(212)
plt.stem(range(1,h_est.size+1),h_est,'-.',label='estimated coef')
plt.plot(range(1,h.size+1),h,'y.',label='original coef')
plt.legend()
plt.xlim(0,p)
plt.ylim(-2,3)
plt.show()
# sudo apt-get install python3-pip
# sudo pip3 install scipy
# sudo pip3 install numpy
import matplotlib.pyplot as plt
import scipy.signal as scpsig
import numpy as np
import time as tm
# Assumsikan h(n) adalah unknown system
h = np.array([0.25, 1.27, -0.3, 0.5, -1.0, 0.75, 0.75, 2.5])
# x(n) adalah input signal
x = np.random.rand(1000)
# d(n) adalah output sinyal
d = scpsig.lfilter(h,1,x)
## Parameter untuk algoritma LMS
# Order dari filter estimasi h_est(n)
p = 15
# Step size atau learning rate
mu = 0.1
# Inisialisasi awal koefisien h_est(0)
h_est = np.zeros(p)
# Inisialisasi state h_est(n) / x(n) = [x(n),x(n-1), ..., x(n-p+1)]
h_stt = np.zeros(p)
# error
e = np.zeros(x.size)
tic = tm.time()
for n in range(0,x.size):
h_stt = np.append(x[n], h_stt[:-1])
e[n] = d[n] - np.dot(h_est,h_stt)
h_est = h_est + mu*e[n]*h_stt
toc = tm.time()
print('%.4f seconds' % (toc-tic))
plt.figure()
plt.subplot(211)
plt.plot(range(1,e.size+1),e)
plt.title('error')
plt.subplot(212)
plt.stem(range(1,h_est.size+1),h_est,'-.',label='estimated coef')
plt.plot(range(1,h.size+1),h,'y.',label='original coef')
plt.legend()
plt.xlim(0,p)
plt.ylim(-2,3)
plt.show()
Julia
# Assumsikan h(n) adalah unknown system
h = [0.25; 1.27; -0.3; 0.5; -1.0; 0.75; 0.75; 2.5];
# x(n) adalah input signal
x = rand(1000);
# d(n) adalah output sinyal
d = filt(h,1,x);
## Parameter untuk algoritma LMS
# Order dari filter estimasi h_est(n)
p = 15;
# Step size atau learning rate
mu = 0.1;
# Inisialisasi awal koefisien h_est(0)
h_est = zeros(p);
# Inisialisasi state h_est(n) / x(n) = [x(n),x(n-1), ..., x(n-p+1)]
h_stt = zeros(p);
# error
e = zeros(size(x));
tic()
for n = 1:length(x)
h_stt = [x[n]; h_stt[1:end-1]];
e[n] = d[n] - (h_est'*h_stt)[];
h_est = h_est + mu*e[n]*h_stt;
end
toc()
using PyPlot
figure()
subplot(211)
plot(1:length(e),e)
title("error")
subplot(212)
stem(1:length(h_est),h_est,"-.",label="estimated coef")
plot(1:length(h),h,"y.",label="original coef")
xlim(0,p)
ylim(-2,3)
legend()
h = [0.25; 1.27; -0.3; 0.5; -1.0; 0.75; 0.75; 2.5];
# x(n) adalah input signal
x = rand(1000);
# d(n) adalah output sinyal
d = filt(h,1,x);
## Parameter untuk algoritma LMS
# Order dari filter estimasi h_est(n)
p = 15;
# Step size atau learning rate
mu = 0.1;
# Inisialisasi awal koefisien h_est(0)
h_est = zeros(p);
# Inisialisasi state h_est(n) / x(n) = [x(n),x(n-1), ..., x(n-p+1)]
h_stt = zeros(p);
# error
e = zeros(size(x));
tic()
for n = 1:length(x)
h_stt = [x[n]; h_stt[1:end-1]];
e[n] = d[n] - (h_est'*h_stt)[];
h_est = h_est + mu*e[n]*h_stt;
end
toc()
using PyPlot
figure()
subplot(211)
plot(1:length(e),e)
title("error")
subplot(212)
stem(1:length(h_est),h_est,"-.",label="estimated coef")
plot(1:length(h),h,"y.",label="original coef")
xlim(0,p)
ylim(-2,3)
legend()
Cara menggunakan matplotlib di Julia dapat dilihat di link berikut http://gear-second19.blogspot.jp/2016/07/jupyter-notebook-aplikasi-web-untuk.html.
Selamat mencoba!
