SciPyで数列の極大値や極小値を求める

時系列データを分析している中で、極大値や極小値を特定したいケースは稀にあります。

極大値/極小値というのは、要は局所的な最大値/最小値のことで、その値の周囲(前後)の値と比較して最大だったり最小だったりする要素のことです。(とても雑な説明。もっと正確な説明はWikipediaの極値のページを参照。)

例えば、[5, 4, 3, 2, 3, 4, 3, 2, 1] みたいな数列があった時、一番値が大きいのは先頭の5なので、これが最大値(極大値でもある)ですが、6番目の4もその近くだけ見ると、[2, 3, 4, 3, 2]となっていて前後の値より大きいので、この6番目の4が極大値です。

今の説明の通りのコードを書いて数列の前後の値と比較して判定したら極大値も極小値も見つかるのですが、それを一発でやってくれるメソッドがSciPyにあるよ、ってのが今回の記事です。

使うのはargrelmin/ argrelmax です。ドキュメントは以下。
scipy.signal.argrelmin — SciPy v1.10.0 Manual
scipy.signal.argrelmax — SciPy v1.10.0 Manual

minとmaxの違いは極小値か極大値かの違いなので、以下の説明はargrelmaxの方でやっていきますね。

ちょっと適当な数列を一個用意して実行してみます。(さっきの例に値を付け足して長くしたものです。)

import numpy as np
from scipy.signal import argrelmax

# サンプルデータを用意
data = np.array([5, 4, 3, 2, 3, 4, 3, 2, 1, 8, 1])

# 極大値のindex
print(argrelmax(data))
# (array([5, 9]),)

# 極大値の値
print(data[argrelmax(data)])
# [4 8]

indexが5(インデックスは0始まりなので6番目の要素)である4と、indexが9(同様に10番目の要素)である8が検出されました。

printした結果を見ていただくと分かる通り、argrelmaxはindexが入ったarrayを0番目の要素に持つタプル、という特殊な形で結果を返してくれます。慣れないとトリッキーに見えますが、それをそのまま使うと極値の値を取り出せるので便利です。

デフォルトでは、直前直後の値だけを見て極大値極小値が判定されますが、例えばノイズを含むデータなどでは実用上検出が多すぎることもあります。
その場合、order(デフォルト1、1以上の整数のみ指定可能)という引数を使うことで、前後何個の要素と比較するかを指定できます。

order=3 とすると、前の3個の値と、後ろの3個の値、合計6個より大きい場合に極大値として判定されます。

data = np.array([2, 1, 1, 4, 1, 1, 5, 1])

# index3の4 と index6の5が極大値として検出される。
print(argrelmax(data, order=1))
# (array([3, 6]),)

# order=3とすると、index6の5だけが極大値として検出される。
print(argrelmax(data, order=3))
# (array([6]),)

上記の例でもう一個着目して欲しい点があります。order=1の時、先頭の2は極大値として検出されませんでした。デフォルトでは、orderの値に関係なく前後に最低1個の要素がないと対象にならないようです。そして、order=3の場合も、後ろから2番目の5が検出されています。orderで指定した数に足りなくても前後に1個以上あれば良いようです。

この、端に関する挙動はmodeという引数で指定できます。デフォルトは”clip”で、これは両端の値は極値として扱われません。ここに”wrap”を指定すると、両端の値も対象になります。

# index=0の2も対象になった
print(argrelmax(data, order=1, mode="wrap"))
# (array([0, 3, 6]),)

もう一つ気をつけないといけないのは、ドキュメントに書かれている通り、前後の値より真に大きくないと極値として扱われません。以下のように前後の値と一致したらダメということです。

data = np.array([1, 1, 2, 3, 3, 2, 1])
print(argrelmax(data))
# (array([], dtype=int64),)

以上が1次元のデータに対する使い方になります。

さて、このargrelmin/ argrelmaxですが、2次元以上のデータに対しても使えます。ドキュメントには2次元の例が載っていますが、3次元でも4次元でもいけます。

2次元、要するに行列形式のデータに対して使ったら、上下左右、できれば斜めも考慮した8方向の隣接データと比較して極大値/極小値を出してくれるのかな?と期待したのですがそういう動きはしておらず、軸(axis)を1個固定してその軸に沿った1次元データとして取り出してそれぞれに対して極大値/極小値の検索をやるようです。方向はaxis引数(デフォルト0)で指定します。

ちょっとでたらめに作ったデータでやってみます。

data = np.array([
        [1, 4, 0, 9, 0, 3, 2],
        [4, 0, 0, 1, 2, 3, 7],
        [2, 9, 2, 0, 9, 0, 7],
        [0, 0, 7, 9, 6, 3, 1],
        [0, 4, 4, 7, 2, 8, 3]
    ])

# axis省略は0と同じ。
print(argrelmax(data))
# (array([1, 2, 2, 3, 3]), array([0, 1, 4, 2, 3]))

print(argrelmax(data, axis=1))
# (array([0, 0, 0, 2, 2, 3, 4, 4]), array([1, 3, 5, 1, 4, 3, 3, 5]))

結果の読み方が慣れないと分かりにくいですが、インデックスの1次元目の配列、2次元目の配列のタプルとして帰ってきてます。

要するに、 (array([1, 2, 2, 3, 3]), array([0, 1, 4, 2, 3])) というのは [1, 0], [2, 1], [2, 4], [3, 2], [3, 3]が極大値だった、という意味です。

そして、またこれも分かりにくいのですが、axis=0 の時、この行列の各列を取り出して極大値を探しています。[1, 0]の値は4 ですが、これはdata[:, 0] = [1, 4, 2, 0, 0] の極大値として検出されており、[3, 2]の7 は data[:, 2] = [0, 0, 2, 7, 4] の極大値として検出されています。

スライスした時に:(コロン)になる次元をaxis引数で指定していると考えたら良いでしょうか。

引数を省略してた時の挙動が想定と違うというか、axis=1を指定した時の方がデフォルトっぽい動きしていますね。こちらは、
(array([0, 0, 0, 2, 2, 3, 4, 4]), array([1, 3, 5, 1, 4, 3, 3, 5]))
が結果として帰りますが、こちらも同様に[0, 1], [0, 3], [0, 5], ….(略) が極大値として検出されています。そしてこれは data[0, :] = [1, 4, 0, 9, 0, 3, 2] の極大値です。

滅多に使わない関数ですしさらにこれを多次元データに使うというのも稀だと思うので、完璧に理解し記憶して使いこなすというよりも、必要になった時に挙動をテストしながら使うのが現実的ではないでしょうか。

コメントを残す

メールアドレスが公開されることはありません。 が付いている欄は必須項目です