読者です 読者をやめる 読者になる 読者になる

旅行好きなソフトエンジニアの備忘録

プログラミングや技術関連のメモを始めました

【Python】 粒子群最適化法による最適化

メタヒューリスティクス Python 最適化 数学

粒子群最適化法による大域的最適化を行うクラスを実装します。粒子群最適化法は

  1. ハイパーパラメータが少ないので扱いやすい
  2. 単峰性が強い関数の最適化に向いている

といった特徴を持っています。最適化問題は慣例として"目隠しした登山家が山の頂上をどうにかして目指す"ことに例えられます。この慣例に従えば、粒子群最適化法は"目隠しした沢山の登山家を山の様々な場所に派遣し、各々が自分が今までに見つけた最も高い標高と、他の登山家が今までに見つけた最も高い標高を頼りに移動先を決めて頂上を目指す"ようなイメージの手法です。自分が持っている情報だけでなく、他の人の情報も使って移動先を決めるため、局所解から抜け出せる余地を持っています。粒子群最適化法では登山家を粒子として表現するのですが、まずは粒子をParticleというクラスで表現します。

import sys
import numpy as np

class Particle(object):
    """粒子群最適化法の粒子を表すクラス"""

    def __init__(self, f, position, velocity, maxs=None, mins=None):
        """
        コンストラクタ
        @param f 最適化を行う関数
        @param position 粒子の初期位置
        @param velocity 粒子の初期速度
        @param maxs 粒子が移動可能な位置最大値
        @param mins 粒子が移動可能な位置最小値
        """
        self.__f = f
        self.__my_best_position = position
        self.__my_best_score = f(position)
        self.__my_position = position
        self.__my_velocity = velocity
        self.__maxs = maxs
        self.__mins = mins
        self.__WEIGHT_V = 0.8
        self.__WEIGHT_ME = 1.0
        self.__WEIGHT_US = 1.0

    @property
    def score(self):
        return self.__my_best_score

    @property
    def position(self):
        return np.array(self.__my_best_position)

    def move(self, best_position):
        """
        粒子を移動させるメソッド
        @param best_position 関数の最小値を与える粒子の位置
        """
        # 粒子の位置を更新する
        self.__my_position += self.__my_velocity
        # 範囲外に出た粒子は範囲内に収める
        if (self.__maxs is not None):
            max_out_of_range_index = self.__my_position > self.__maxs
            self.__my_position[max_out_of_range_index] = self.__maxs[max_out_of_range_index]
        if (self.__mins is not None):
            min_out_of_range_index = self.__my_position < self.__mins
            self.__my_position[min_out_of_range_index] = self.__mins[min_out_of_range_index]
        # 評価値を計算する
        score = self.__f(self.__my_position)
        # 最良解を更新する
        if (score < self.__my_best_score):
            self.__my_best_position = self.__my_position
            self.__my_best_score = score
        # 速度を更新する
        item1 = self.__WEIGHT_V*self.__my_velocity
        random_me = np.random.rand(self.__my_position.size)
        item2 = self.__WEIGHT_ME*random_me*(self.__my_best_position - self.__my_position)
        random_us = np.random.rand(self.__my_position.size)
        item3 = self.__WEIGHT_US*random_us*(best_position - self.__my_position)
        self.__my_velocity = item1 + item2 + item3
        # 範囲外に出た粒子の速度は0とする
        if (self.__maxs is not None): self.__my_velocity[max_out_of_range_index] = 0
        if (self.__mins is not None): self.__my_velocity[min_out_of_range_index] = 0

次に、粒子を制御するクラスをParticleSwarmOptimizationとして作成します。このクラスは単に粒子群から最適解に最も近い粒子の位置を探し、粒子の位置を更新しているだけです。

import sys
import Particle

class ParticleSwarmOptimization(object):
    """粒子群最適化法により最適化を行うクラス"""

    def __init__(self, particles):
        """
        コンストラクタ
        @param particles 粒子のリスト
        """
        self.__particles = particles
        self.__best_score = sys.float_info.max
        for particle in self.__particles:
            if (particle.score < self.__best_score):
                self.__best_score = particle.score
                self.__best_position = particle.position

    @property
    def best_score(self):
        return self.__best_score

    @property
    def best_position(self):
        return self.__best_position

    def position(self, index):
        return [particle.position[index] for particle in self.__particles]

    def update(self):
        """
        粒子の位置を更新するメソッド
        """
        for particle in self.__particles:
            particle.move(self.__best_position)

        for particle in self.__particles:
            if (particle.score < self.__best_score):
                self.__best_score = particle.score
                self.__best_position = particle.position

作成したParticleクラスとParticleSwarmOptimizationクラスを利用してn=2のRosenbrock関数の最適化を行います。
f(x)=100(x_1-x_0)^2+(1-x_0)^2 \ \ \ \ \ (-2.048 \lt x_i \lt 2.048)

import numpy as np
import matplotlib.pyplot as plt
from ParticleSwarmOptimization import ParticleSwarmOptimization
from Particle import Particle

if __name__ == '__main__':
    # 最適化する関数
    def f(x): return 100.0*(x[1]-x[0]**2)**2 + (1.0-x[0])**2
    # 変数が取り得る最大値
    maxs = np.array([2.048, 2.048])
    # 変数が取り得る最小値
    mins = np.array([-2.048, -2.048])
    # ばらまく粒子の個数
    PARTICLE_COUNT = 100
    # 粒子を作成する
    particles = []
    for i in range(PARTICLE_COUNT):
        position = (np.random.rand(2) - 0.5)*4.0
        velocity = (np.random.rand(2) - 0.5)/10.0
        particle = Particle(f, position, velocity, maxs, mins)
        particles.append(particle)
    # ParticleSwarmOptimizationクラスの生成
    pso = ParticleSwarmOptimization(particles)
    # 計算回数
    ITERATION = 100
    # 計算開始
    for i in range(ITERATION):
        pso.update()
        if i == 0 or i == 9 or i == 49 or i == 99:
            plt.figure()
            plt.xlim([-2.048, 2.048])
            plt.ylim([-2.048, 2.048])
            plt.grid()
            plt.title('i = ' + str(i))
            plt.scatter(pso.position(0), pso.position(1))
            plt.show()

i=0, i=9, i=49, i=99の時の粒子の様子を図示します。計算が進むにつれて、ばらばらに散らばっていた粒子が最適解である(1,1)に集まる様子を確認できます。 f:id:ni4muraano:20170125204342j:plain:w280f:id:ni4muraano:20170125204346j:plain:w280f:id:ni4muraano:20170125204351j:plain:w280f:id:ni4muraano:20170125204354j:plain:w280

メタヒューリスティクスとナチュラルコンピューティング

メタヒューリスティクスとナチュラルコンピューティング