概要

LSTM（Long Short-Term Memory）を使って、正弦波の振幅と位相を求めてみます。

他の手法との比較

正弦波の振幅と位相をScipyのカーブフィットで求める
https://qiita.com/nnn112358/items/35cfee57bd1d8f21147b

正弦波の振幅と位相をSciPyのFFTで求める
https://qiita.com/nnn112358/items/dfded784541045f988d4

正弦波の振幅と位相をLSTMで求める
https://qiita.com/nnn112358/items/1c9b077a7a413ee1dd9f

正弦波の振幅と位相をRNNで求める
https://qiita.com/nnn112358/items/63b53fb192b7f10d509a

LSTMの説明

LSTMは、Long Short-Term Memoryの略称で、日本語では「長・短期記憶」と訳されます。再帰型ニューラルネットワーク（RNN）の一種で、時系列データの処理に特に適した構造を持つニューラルネットワークです。

LSTMのアーキテクチャ

LSTMは、以下の4つの主要なコンポーネントで構成されています。
セル: 過去の情報を保持する記憶装置です。
忘却ゲート: 過去の情報をどの程度忘れるかを制御します。
入力ゲート: 新しい情報を取り入れるかどうかを制御します。
出力ゲート: セルの状態をどの程度出力に反映させるかを制御します。
これらのコンポーネントは、それぞれシグモイド関数と点乗積と呼ばれる数学的操作を用いて実装されています。

LSTMの利点
長期的な依存関係を学習できる: 勾配消失問題を解決することで、長期的な依存関係を学習することができます。
複雑な時系列データを処理できる: 長期的な依存関係を学習できるため、複雑な時系列データを処理することができます。
幅広いタスクに適用できる: 自然言語処理、時系列データ予測、ロボット制御など、幅広いタスクに適用することができます。
LSTMの欠点
計算コストが高い: 他のニューラルネットワークと比べて、計算コストが高くなります。
ハイパーパラメータが多い: 多くのハイパーパラメータを調整する必要があるため、モデルの訓練が難しい場合があります。

Python開発環境のインストール

Python環境として、Minicondaをインストールする手順を説明します。
WSLにUbuntu22.04 をインストールして使用しました。

MinicondaのLinuxへのインストール手順

Minicondaのインストーラをダウンロード:
まず、Minicondaの公式ウェブサイトからインストーラをダウンロードします。以下のコマンドで最新版のMinicondaインストーラをダウンロードします。
```
wget https://repo.anaconda.com/miniconda/Miniconda3-latest-Linux-x86_64.sh
```
インストーラの実行:
ダウンロードしたスクリプトを実行してMinicondaをインストールします。インストール先を指定する場合は、オプションでディレクトリを指定できます。
```
bash Miniconda3-latest-Linux-x86_64.sh
```
インストールプロセスが開始されると、いくつかのプロンプトが表示されます。これに従って進めてください。
- 使用許諾契約書が表示されたら、スペースキーを押して全文を表示し、yesと入力して同意します。
- インストールディレクトリの確認が表示されるので、デフォルトのままで良ければそのままエンターキーを押します。カスタムディレクトリにインストールしたい場合は、パスを入力してエンターキーを押します。
- インストールが完了したら、condaコマンドを使えるようにするために、環境変数の設定を行うかどうかのプロンプトが表示されます。これもyesと入力して進めます。
環境変数の設定:
インストールが完了したら、Minicondaのバイナリディレクトリを環境変数PATHに追加します。多くの場合、インストールスクリプトが自動的にこれを行いますが、手動で行う場合は以下のコマンドを実行します。
```
export PATH=~/miniconda3/bin:$PATH
```
この設定を永続化するために、~/.bashrc（または使用しているシェルの設定ファイル）に追加します。
```
echo 'export PATH=~/miniconda3/bin:$PATH' >> ~/.bashrc
source ~/.bashrc
```
インストールの確認:
condaコマンドを実行して、Minicondaのインストールが正しく行われたことを確認します。
```
conda --version
```
これでMinicondaのバージョンが表示されれば、インストールは成功です。
Condaの初期設定とアップデート:
インストール後、初期設定を行い、Conda自体とデフォルトパッケージをアップデートします。
```
conda init
conda update conda
```

これでMinicondaのインストールは完了です。

仮想環境の作成

Minicondaのインストールが完了したら、仮想環境を作成します。

仮想環境の作成:

 conda create -n new_env python=3.8 tensorflow keras matplotlib scikit-learn pykalman

仮想環境のアクティベート:

conda activate new_env

仮想環境のディアクティベート:

仮想環境を終了するには、以下のコマンドを実行します。

conda deactivate

サンプルソース

以下に、LSTMを用いて正弦波の振幅と位相を求める手順を示します。この手順では、LSTMモデルを用いて信号を予測し、その予測結果から振幅と位相を推定します。

$ python lstm_fit.py
Original amplitude: 2.0000
Train amplitude: 2.0952
Test amplitude: 2.0965
Original phase: 0.7854
Train phase: 0.7490
Test phase: 0.7157

データの準備:
正弦波にノイズを加えたデータを生成し、スケーリングしてLSTMモデルに入力できるように整形します。

lstm_fit.py

import numpy as np
import matplotlib.pyplot as plt
import tensorflow as tf
from sklearn.preprocessing import MinMaxScaler

# パラメータ設定
N = 1000  # サンプル数
T = 0.1  # サンプリング間隔
t = np.linspace(0.0, N*T, N)  # 時間軸の生成

# 既知の周期と周波数
known_period = 5.0  # 既知の周期
known_freq = 1.0 / known_period  # 周波数の計算

# 元の正弦波を生成
amplitude_original = 2.0  # 元の正弦波の振幅
phase_original = np.pi / 4  # 元の正弦波の位相
y_original = amplitude_original * np.sin(2.0 * np.pi * known_freq * t + phase_original)  # 元の正弦波データ

# ノイズを加える
np.random.seed(42)  # 再現性のためにランダムシードを設定
noise = np.random.normal(0, 0.5, N)  # 正規分布に従うノイズを生成
y_noisy = y_original + noise  # ノイズを加えたデータ

# データのスケーリング
scaler = MinMaxScaler(feature_range=(-1, 1))
y_noisy_scaled = scaler.fit_transform(y_noisy.reshape(-1, 1)).reshape(-1)

データの整形:
LSTMモデルに入力するための形状にデータを変換し、訓練データとテストデータに分割します。LSTMモデルに入力するために、データを適切な形状に整形します。

python lstm_fit.py

# LSTMの入力データを作成
def create_dataset(data, time_step=1):
    X, Y = [], []
    for i in range(len(data) - time_step - 1):
        a = data[i:(i + time_step)]
        X.append(a)
        Y.append(data[i + time_step])
    return np.array(X), np.array(Y)

time_step = 100
X, Y = create_dataset(y_noisy_scaled, time_step)

# LSTMに入力するための形状に変換
X = X.reshape(X.shape[0], X.shape[1], 1)

# 訓練とテストの分割
train_size = int(len(X) * 0.67)
test_size = len(X) - train_size
X_train, X_test = X[0:train_size], X[train_size:len(X)]
Y_train, Y_test = Y[0:train_size], Y[train_size:len(Y)]

LSTMモデルの構築と訓練
LSTM層と全結合層からなるモデルを構築し、訓練データでモデルを訓練します。

# LSTMモデルの構築
model = tf.keras.Sequential()
model.add(tf.keras.layers.LSTM(50, return_sequences=True, input_shape=(time_step, 1)))
model.add(tf.keras.layers.LSTM(50, return_sequences=False))
model.add(tf.keras.layers.Dense(1))
model.compile(optimizer='adam', loss='mean_squared_error')

# モデルの訓練
model.fit(X_train, Y_train, validation_data=(X_test, Y_test), epochs=20, batch_size=64, verbose=1)

予測と振幅・位相の推定:

LSTMモデルを使って信号を予測し、その結果から振幅と位相を推定します。訓練済みのLSTMモデルを使って信号を予測し、予測データを逆スケーリングして元のスケールに戻します。予測データから振幅と位相を推定し、その結果を表示します。
LSTMを用いて正弦波の振幅と位相を推定し、それをグラフで確認することができます。LSTMは特に時系列データの解析に強力であり、ノイズが多い場合や非線形なデータに対しても有効です。

# 予測
train_predict = model.predict(X_train)
test_predict = model.predict(X_test)

# 逆スケーリング
train_predict = scaler.inverse_transform(train_predict)
test_predict = scaler.inverse_transform(test_predict)

# フィッティング関数を定義
def fit_sine_wave(t, y):
    from scipy.optimize import curve_fit

    def sine_wave(t, A, f, phi):
        return A * np.sin(2 * np.pi * f * t + phi)

    guess_freq = 1.0 / known_period
    guess_amplitude = np.std(y) * 2**0.5
    guess_phase = 0
    guess_offset = np.mean(y)
    p0 = [guess_amplitude, guess_freq, guess_phase]

    popt, _ = curve_fit(sine_wave, t, y, p0=p0)
    return popt

# 予測データから振幅と位相を推定
train_amplitude, train_frequency, train_phase = fit_sine_wave(t[:len(train_predict)], train_predict.ravel())
test_amplitude, test_frequency, test_phase = fit_sine_wave(t[len(train_predict):len(train_predict) + len(test_predict)], test_predict.ravel())

# 結果の表示
print(f"Original amplitude: {amplitude_original:.4f}")
print(f"Train amplitude: {train_amplitude:.4f}")
print(f"Test amplitude: {test_amplitude:.4f}")

print(f"Original phase: {phase_original:.4f}")
print(f"Train phase: {train_phase:.4f}")
print(f"Test phase: {test_phase:.4f}")

# グラフ描画
plt.figure(figsize=(14, 8))

plt.subplot(2, 1, 1)
plt.plot(t[:len(train_predict)], y_noisy[:len(train_predict)], label='Noisy Train Data', alpha=0.5)
plt.plot(t[:len(train_predict)], train_predict, label='LSTM Train Predict', linestyle='--')
plt.xlabel('Time')
plt.ylabel('Amplitude')
plt.legend()

plt.subplot(2, 1, 2)
plt.plot(t[len(train_predict):len(train_predict) + len(test_predict)], y_noisy[len(train_predict):len(train_predict) + len(test_predict)], label='Noisy Test Data', alpha=0.5)
plt.plot(t[len(train_predict):len(train_predict) + len(test_predict)], test_predict, label='LSTM Test Predict', linestyle='--')
plt.xlabel('Time')
plt.ylabel('Amplitude')
plt.legend()

plt.tight_layout()
plt.show()