PythonでCFF-RLS(Constant Forgetting Factor - Recursive Least Square: 忘却機能つき逐次最小二乗法)アルゴリズムを実装する

今回はPythonで忘却機能つきのRLSアルゴリズムを実装してみます。

以前本ブログで書いたRLSアルゴリズムに忘却機能をつけたものですね。

RLSアルゴリズムは値が振動せず収束してくれるという特徴がありますが、例えば推定しようとするシステムが時変であり追従して欲しい場合にはそれがデメリットになります。

それを回避するために、今から前の時点の情報を忘れるかのように損失関数を改良したものが忘却機能つきRLSです。

こうすることで、古い情報を捨て新しい情報のみを使うようになるので、時変システムに対しても追従してくれるようになってくれます。

ちなみに、タイトルに書いているCFF(Constant Forgetting Factor)についてですが、もしかしたらあまり一般的な用語ではないかもしれません（笑）可変忘却機能をVFF(VariableForgetting Factor)というので、その対比でCFFとしてみました。

さて、CFF-RLSのコードは以下のようになります。途中で特性が変わるシステムに対してCFF-RLSにてパラメータ逐次推定を行っています。

[In]

import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline
np.random.seed(0)

[In]

class Plant1(object):
    def __init__(self):
        self.k0 = 0.0
        self.k1 = 0.0
        self.k2 = 0.0
        self.k3 = 0.0
        return
    def get(self, signal):
        self.k3 = self.k2
        self.k2 = self.k1
        self.k1 = self.k0
        self.k0 = signal
        return ((self.k0 + self.k1 + self.k2 + self.k3) * 0.25) + np.random.normal(loc=0,scale = 0.1)

[In]

class Plant2(object):
    def __init__(self):
        self.k0 = 0.0
        self.k1 = 0.0
        self.k2 = 0.0
        self.k3 = 0.0
        return
    def get(self, signal):
        self.k3 = self.k2
        self.k2 = self.k1
        self.k1 = self.k0
        self.k0 = signal
        return ((self.k0 + self.k1 + self.k2 + self.k3) * 0.5) + np.random.normal(loc=0,scale = 0.1)

[In]

# CFF: Constant Forgetting Factor
# RLS: Recursive Least Square
class CFF_RLS(object):
    def __init__(self, forget):
        self.len = 4
        self.theta = np.zeros((self.len,1))
        self.fai = np.zeros((self.len,1))
        self.cov = np.identity(self.len) * 100
        self.forget = forget
        return

    def update(self, x, y):
        self.fai[1:] = self.fai[0:-1]
        self.fai[0] = x
        gain = (self.cov @ self.fai) / (self.fai.T @ self.cov @ self.fai + self.forget)
        error = y - (self.fai.T @ self.theta)
        self.theta = self.theta + (gain @ error)
        self.cov = (self.cov - (gain @ self.fai.T @ self.cov)) / self.forget
        return

[In]

plant1 = Plant1()
plant2 = Plant2()
time = np.arange(600)
input_signal_array = []
output_signal_array = []
for i in time:
    input_signal = np.sin(2.0 * np.pi * i /20.0) + np.random.uniform(-1, 1)
    if i < 300:
        output_signal = plant1.get(input_signal)
    else:
        output_signal = plant2.get(input_signal)
    input_signal_array.append(input_signal)
    output_signal_array.append(output_signal)
plt.figure(figsize=(12, 4))
plt.subplot(1, 2, 1)
plt.plot(time, input_signal_array)
plt.title("Input signal")
plt.subplot(1, 2, 2)
plt.plot(time, output_signal_array)
plt.title("Output signal")
plt.show()

[Out]

[In]

log = []
rls = CFF_RLS(forget=0.97)
for t in time:
    rls.update(input_signal_array[t], output_signal_array[t])
    log.append(rls.theta)
        
for i in range(rls.theta.shape[0]):
    plt.plot(time, [log[j][i] for j in time], label="g_" + str(i))
plt.grid(True)
plt.legend()
plt.ylim(-0.1, 0.7)
plt.show()

[Out]

時刻（横軸）が300で特性を変えていますが、それに対してちゃんと速く追従してくれてるのが分かるかと思います。

今回は以上です！