先从这个网站找到一个简单并且有点不靠谱的一维PSO的python实现。特点是简单，容易看懂，并且用pygame将优化过程保存成图片，方便以后查看。改为多维的实现，能运行后发现调用libsvm进行交叉校验每次都会打印出信息，大大减慢了算法的速度！于是昨晚就调试怎么去掉这个了～
今天早上起来想到那个保存成图片的功能不错，多维算法可以指定其中一维来实现，大不了多搞几个这样的类给PSO，分别保存每一维的优化过程就是了！实现后发现有问题，蓝色的最优点居然会跑到fitness函数外面！幸亏这个保存图片的功能，不然还真发现不了这个bug呢！经过调试发现是因为复制列表时仅复制了引用，这样在种群继续进化时就直接导致最优点的值也变了！找到问题就好办了，用deepcopy复制列表就OK啦～
多维PSO的代码如下：
Python语言: 粒子群（PSO）算法
# Copyright (C) 2004, Maxime Biais <maxime@biais.org>
# Copyright (C) 2009, Ace Strong <acestrong@gmail.com>
#
# This program is free software; you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation; either version 2 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program; if not, write to the Free Software
# Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
#
# $Id: pso.py,v 2.0 2009/02/02 10:00:45 ace Exp $
from random import uniform
from copy import deepcopy
class PSO:
def __init__(self, pop_size, particle_size, particle_scope, C1, C2, w, max_w, min_w, max_iter, func, extra_args=None):
'''
pop_size: size of population
particle_size: dimension of each particle
particle_scope: value range of each particle's dimension
C1: cognitive factor following personal best
C2: social facotr following global best
w: const inertial factor
min_w: min range of inertial factor
max_w: max range of inertial factor
max_iter: max iteration of algrithm
func: fitness evaluation function
'''
# fitness function
self.func = func
# used for fitness function
self.extra_args = extra_args
self.pop = []
# converging factor
self.r = 0.729
self.particle_size = particle_size
# parse range of each dimension
self.min_range = particle_scope[::2]
self.max_range = particle_scope[1::2]
# 0: position, 1: velocity, 2: fitness
for i in xrange(pop_size):
self.pop.append(self.initParticle())
self.evaluate()
# global best
self.gbest = deepcopy(self.pop[0])
# personal best, each partical has one
self.pbest = deepcopy(self.pop)
# weight following personal best, cognitive factor
self.C1 = C1
# weight following global best, social facotr
self.C2 = C2
# weight following current speed, inertial factor
self.w = w
# used by un-const method
self.max_w = max_w
self.min_w = min_w
# max iteration number
self.max_iter = max_iter
# current iteration number
self.curr_iter = 0
def initParticle(self):
particle = []
# position, in particle_size dimension
position = []
for j in xrange(self.particle_size):
position.append(uniform(self.min_range[j], self.max_range[j]))
particle.append(position)
# velocity, in particle_size dimension
velocity = []
for j in xrange(self.particle_size):
velocity.append(uniform(-1, 1))
particle.append(velocity)
# fitness
particle.append(0)
return particle
def update_velocity(self):
# linear descending w
w = self.max_w - self.curr_iter*((self.max_w - self.min_w)/self.max_iter)
# # fixed w
# w = self.w
# # nonlinear descending w, concave function
# w = (self.max_w - self.min_w)*(self.curr_iter/self.max_iter)**2 \
# + (self.min_w - self.max_w)*(2*self.curr_iter/self.max_iter) \
# + self.max_w
# # nonlinear descending w, concave function
# w = self.min_w*(self.max_w/self.min_w)**(1/(1+10*self.curr_iter/self.max_iter))
i = 0
for p in self.pop:
for j in xrange(self.particle_size):
p[1][j] = w * p[1][j] + uniform(0, self.C1) * (self.pbest[i][0][j] \
- p[0][j]) + uniform(0, self.C2) * (self.gbest[0][j] - p[0][j])
i += 1
def evaluate(self):
for p in self.pop:
p[2] = self.func(p[0], self.extra_args)
def move(self):
i = 0
for p in self.pop:
for j in xrange(self.particle_size):
p[0][j] += self.r * p[1][j]
if p[0][j] > self.max_range[j] or p[0][j] < self.min_range[j]:
self.pop[i] = self.initParticle()
break
# if self.func(self.gbest[0],self.extra_args)!=self.gbest[2]:
# print "error in move!"
i += 1
def run(self, update_func=None):
for i in xrange(self.max_iter):
print "current iter = %d" % self.curr_iter
if update_func:
update_func()
self.update_velocity()
self.move()
self.evaluate()
# current best of personal bests
cpbest = deepcopy(self.pbest[0])
for k in xrange(len(self.pop)):
if self.pop[k][2] < self.pbest[k][2]:
self.pbest[k] = deepcopy(self.pop[k])
if self.pbest[k][2] < cpbest[2]:
cpbest = self.pbest[k]
if cpbest[2] < self.gbest[2]:
# print "changing gbest from %s" % self.gbest
self.gbest = cpbest
# print "to %s" % self.gbest
# if self.func(self.gbest[0],self.extra_args)!=self.gbest[2]:
# print "error in run!"
# print self.gbest
self.curr_iter += 1
def __str__(self):
ret = ""
for i in self.pop:
ret += str(i) + "\n"
return ret
import pygame
import time
class PygamePrinter:
'''draw the given dimension's evolution
'''
def __init__(self, pso, w=400, h=300, dimension=0, extra_args=None):
self.calls = 0
self.w = w
self.h = h
self._init_pygame()
self.pso = pso
self.dimension = dimension
self.extra_args = extra_args
def _init_pygame(self):
self.screen = pygame.display.set_mode((self.w, self.h), 0, 8)
self.backcolor = (0, 0, 0)
self.funccolor = (255, 255, 255)
self.partcolor = (255, 0, 0)
self.elitecolor = (0, 0, 255)
def draw_point(self, color, x, y, size=3):
pygame.draw.rect(self.screen, color, (x - size, y - size, \
size*2, size*2))
def p2p(self, x, y):
return (x + 1) * 200, y * 300
def draw_func(self):
for i in range(self.w):
x = i / (self.w / float((self.pso.max_range[self.dimension] \
- self.pso.min_range[self.dimension]))) \
+ self.pso.min_range[self.dimension]
y = self.pso.func((x, self.dimension), self.extra_args)
rh = y * (self.h / 2.) + (self.h / 2.)
self.draw_point(self.funccolor, i, rh, 1)
def _draw_xy(self, color, x, y):
self.draw_point(color, (x + self.pso.max_range[self.dimension]) * self.w \
/ float((self.pso.max_range[self.dimension] - self.pso.min_range[self.dimension])), \
(y + 1)*self.h/2)
def draw_pop(self):
for i in self.pso.pop:
self._draw_xy(self.partcolor, i[0][self.dimension], i[2])
i = self.pso.gbest
self._draw_xy(self.elitecolor, i[0][self.dimension], i[2])
def __call__(self):
#time.sleep(0.5)
self.screen.fill(self.backcolor)
self.screen.lock()
self.draw_func()
self.draw_pop()
self.screen.unlock()
pygame.display.flip()
if self.calls % 2 == 0:
# print "call=%d" % self.calls
pygame.image.save(self.screen, "pso-%d.bmp" % self.calls)
self.calls += 1
import math
def testFunc(arg, extra_args):
x = arg[0]
return math.cos(x) * math.exp(math.sin(x)) * math.sin(x) / 1.5
def test():
import math
# func = lambda x:math.cos(x*math.sin(x*0.3)-x) / 1.5
# func = lambda x:math.cos(x) * math.exp(math.sin(x)) * math.sin(x) / 1.5
p = PSO(15, 1, (-4.5, 4.5), C1=2, C2=2, w=0.5, max_w=0.95, min_w=0.4, max_iter=20, func=testFunc)
printer = PygamePrinter(p)
p.run(update_func=printer)
# p.run()
print p
print p.gbest
if __name__ == "__main__":
test()
# Copyright (C) 2009, Ace Strong <acestrong@gmail.com>
#
# This program is free software; you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation; either version 2 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program; if not, write to the Free Software
# Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
#
# $Id: pso.py,v 2.0 2009/02/02 10:00:45 ace Exp $
from random import uniform
from copy import deepcopy
class PSO:
def __init__(self, pop_size, particle_size, particle_scope, C1, C2, w, max_w, min_w, max_iter, func, extra_args=None):
'''
pop_size: size of population
particle_size: dimension of each particle
particle_scope: value range of each particle's dimension
C1: cognitive factor following personal best
C2: social facotr following global best
w: const inertial factor
min_w: min range of inertial factor
max_w: max range of inertial factor
max_iter: max iteration of algrithm
func: fitness evaluation function
'''
# fitness function
self.func = func
# used for fitness function
self.extra_args = extra_args
self.pop = []
# converging factor
self.r = 0.729
self.particle_size = particle_size
# parse range of each dimension
self.min_range = particle_scope[::2]
self.max_range = particle_scope[1::2]
# 0: position, 1: velocity, 2: fitness
for i in xrange(pop_size):
self.pop.append(self.initParticle())
self.evaluate()
# global best
self.gbest = deepcopy(self.pop[0])
# personal best, each partical has one
self.pbest = deepcopy(self.pop)
# weight following personal best, cognitive factor
self.C1 = C1
# weight following global best, social facotr
self.C2 = C2
# weight following current speed, inertial factor
self.w = w
# used by un-const method
self.max_w = max_w
self.min_w = min_w
# max iteration number
self.max_iter = max_iter
# current iteration number
self.curr_iter = 0
def initParticle(self):
particle = []
# position, in particle_size dimension
position = []
for j in xrange(self.particle_size):
position.append(uniform(self.min_range[j], self.max_range[j]))
particle.append(position)
# velocity, in particle_size dimension
velocity = []
for j in xrange(self.particle_size):
velocity.append(uniform(-1, 1))
particle.append(velocity)
# fitness
particle.append(0)
return particle
def update_velocity(self):
# linear descending w
w = self.max_w - self.curr_iter*((self.max_w - self.min_w)/self.max_iter)
# # fixed w
# w = self.w
# # nonlinear descending w, concave function
# w = (self.max_w - self.min_w)*(self.curr_iter/self.max_iter)**2 \
# + (self.min_w - self.max_w)*(2*self.curr_iter/self.max_iter) \
# + self.max_w
# # nonlinear descending w, concave function
# w = self.min_w*(self.max_w/self.min_w)**(1/(1+10*self.curr_iter/self.max_iter))
i = 0
for p in self.pop:
for j in xrange(self.particle_size):
p[1][j] = w * p[1][j] + uniform(0, self.C1) * (self.pbest[i][0][j] \
- p[0][j]) + uniform(0, self.C2) * (self.gbest[0][j] - p[0][j])
i += 1
def evaluate(self):
for p in self.pop:
p[2] = self.func(p[0], self.extra_args)
def move(self):
i = 0
for p in self.pop:
for j in xrange(self.particle_size):
p[0][j] += self.r * p[1][j]
if p[0][j] > self.max_range[j] or p[0][j] < self.min_range[j]:
self.pop[i] = self.initParticle()
break
# if self.func(self.gbest[0],self.extra_args)!=self.gbest[2]:
# print "error in move!"
i += 1
def run(self, update_func=None):
for i in xrange(self.max_iter):
print "current iter = %d" % self.curr_iter
if update_func:
update_func()
self.update_velocity()
self.move()
self.evaluate()
# current best of personal bests
cpbest = deepcopy(self.pbest[0])
for k in xrange(len(self.pop)):
if self.pop[k][2] < self.pbest[k][2]:
self.pbest[k] = deepcopy(self.pop[k])
if self.pbest[k][2] < cpbest[2]:
cpbest = self.pbest[k]
if cpbest[2] < self.gbest[2]:
# print "changing gbest from %s" % self.gbest
self.gbest = cpbest
# print "to %s" % self.gbest
# if self.func(self.gbest[0],self.extra_args)!=self.gbest[2]:
# print "error in run!"
# print self.gbest
self.curr_iter += 1
def __str__(self):
ret = ""
for i in self.pop:
ret += str(i) + "\n"
return ret
import pygame
import time
class PygamePrinter:
'''draw the given dimension's evolution
'''
def __init__(self, pso, w=400, h=300, dimension=0, extra_args=None):
self.calls = 0
self.w = w
self.h = h
self._init_pygame()
self.pso = pso
self.dimension = dimension
self.extra_args = extra_args
def _init_pygame(self):
self.screen = pygame.display.set_mode((self.w, self.h), 0, 8)
self.backcolor = (0, 0, 0)
self.funccolor = (255, 255, 255)
self.partcolor = (255, 0, 0)
self.elitecolor = (0, 0, 255)
def draw_point(self, color, x, y, size=3):
pygame.draw.rect(self.screen, color, (x - size, y - size, \
size*2, size*2))
def p2p(self, x, y):
return (x + 1) * 200, y * 300
def draw_func(self):
for i in range(self.w):
x = i / (self.w / float((self.pso.max_range[self.dimension] \
- self.pso.min_range[self.dimension]))) \
+ self.pso.min_range[self.dimension]
y = self.pso.func((x, self.dimension), self.extra_args)
rh = y * (self.h / 2.) + (self.h / 2.)
self.draw_point(self.funccolor, i, rh, 1)
def _draw_xy(self, color, x, y):
self.draw_point(color, (x + self.pso.max_range[self.dimension]) * self.w \
/ float((self.pso.max_range[self.dimension] - self.pso.min_range[self.dimension])), \
(y + 1)*self.h/2)
def draw_pop(self):
for i in self.pso.pop:
self._draw_xy(self.partcolor, i[0][self.dimension], i[2])
i = self.pso.gbest
self._draw_xy(self.elitecolor, i[0][self.dimension], i[2])
def __call__(self):
#time.sleep(0.5)
self.screen.fill(self.backcolor)
self.screen.lock()
self.draw_func()
self.draw_pop()
self.screen.unlock()
pygame.display.flip()
if self.calls % 2 == 0:
# print "call=%d" % self.calls
pygame.image.save(self.screen, "pso-%d.bmp" % self.calls)
self.calls += 1
import math
def testFunc(arg, extra_args):
x = arg[0]
return math.cos(x) * math.exp(math.sin(x)) * math.sin(x) / 1.5
def test():
import math
# func = lambda x:math.cos(x*math.sin(x*0.3)-x) / 1.5
# func = lambda x:math.cos(x) * math.exp(math.sin(x)) * math.sin(x) / 1.5
p = PSO(15, 1, (-4.5, 4.5), C1=2, C2=2, w=0.5, max_w=0.95, min_w=0.4, max_iter=20, func=testFunc)
printer = PygamePrinter(p)
p.run(update_func=printer)
# p.run()
print p
print p.gbest
if __name__ == "__main__":
test()
没有评论:
发表评论