極小値を探す(？)

python

z = f(x, y)としたときの、zの極小値を探そうという話。
といっても、fの関数の形はわからず、zの値は行列で表わされているとします。
とりあえず書けたのでブログに載せてみます。あってるかな。。。


http://old.nabble.com/Finding-local-minima-of-greater-than-a-given-depth-td18988309.html
これはScipyのメーリングリストのようですが、この中で極値について議論しています。
1変数の時、ある要素の前後が、その要素より大きい時にその値を極小値としています。


というわけで、2変数の場合は、周りの要素がその要素より大きいときに極小値とすればいいのかな？

（図 この時、x=1, y=1で極小値1を取る）


この考えに基づいて、関数local_minimaを書いてみました。
引数にzの情報が入った行列を渡し、

(行, 列, その時のzの値)

のタプルを返します。

行列で、四方の要素と比較する時は、最も外側の行と列の扱いが問題になります。
今回はその外側の行と列を無視しています。(外側の行と列に極小値があってもわからない)


（図 外側は四方に全てに値があるわけではない。というわけで注意が必要）


というわけでソースコードを示します。
試しにz = cos(x) + cos(y)について極小値を見ています。

#! /usr/bin/python
# -*- coding: UTF-8 -*-
u"""
行列の極小値を求める
"""
import numpy as np
import matplotlib.pyplot as plt


def local_minima(mat):
    u"""
    与えられた行列matの極小値を求める。
    返り値は、行列の行, 列, その位置での値である。
    """
    minima_ls = []
    # 行列の一番外側の行、列は調べない。
    # そのため、rangeが1, len-2となっている。
    for i in range(1, len(mat)-2):
        for j in range(1, len(mat[i])-2):
            if mat[i,j] < mat[i-1,j] and mat[i,j] < mat[i+1,j] and mat[i,j] < mat[i,j-1] and mat[i,j] < mat[i,j+1]:
                minima_ls.append((i, j, mat[i,j])) # 
    return minima_ls


def main():
    # 適当な行列の作成
    x = np.arange(0, 10, 0.1)
    y = np.arange(0, 10, 0.1)
    X, Y = np.meshgrid(x, y)
    # 適当な関数で、Z方向の値を作る
    Z = np.cos(X) + np.cos(Y)

    # Local Minimumを探す
    ls = local_minima(Z)
    for s in ls:
        y_index, x_index, z_value = s # 行列の行がz, 列がxに対応するので、y_index, x_indexの順番
        print u"x=%f, y=%f の時にz=%f" % (x[x_index], y[y_index], z_value)

    # 試しにプロットしてみる。
    plt.contourf(X, Y, Z)
    plt.show()


if __name__ == '__main__':
    main()

この結果は

c:\>python localminimum.py
x=3.100000, y=3.100000 の時にz=-1.998270
x=9.400000, y=3.100000 の時にz=-1.998828
x=3.100000, y=9.400000 の時にz=-1.998828
x=9.400000, y=9.400000 の時にz=-1.999386

となり、グラフは以下のようになりました。
たぶんあってる？