优化算法-随机搜索、爬山法、模拟退火法、遗传算法1

随机搜索、爬山法、模拟退火法、遗传算法1

案例1：行程优化问题

毕业多年，曾经上学时的好友们都已分散到不同城市工作生活。近日你联系组织了这群生活在不同城市的好友们，组团去日本旅行。约定到达东京羽田国际机场(Tokyo)后，在机场租车，一起出发去酒店。旅行结束后，一起成乘车去机场，还车，然后各自坐飞机离开日本返程。

你该如何协调大家乘坐的去程、返程飞机时间，使得大家在机场等待的时间最短？

假设你统计到出行人员及出发城市为:
（老高：杭州），（老张：成都），（老刘：武汉），（老宋：北京），（老何：上海），（老牛：上海），（老李：杭州），（老赵：厦门），（老钱：广州），（老曾：深圳）。

假设你也已经从航程网站上爬取了这几个城市的、对应出发日/返程日的飞机航班和价格表，并整理成了数据表：
格式：出发城市，到达城市，起飞时间，落地时间，机票价格，飞机信息1，飞机信息2
HangZhou,Tokyo,09:55,13:50,3188,CA929,中国国航
Tokyo,HangZhou,10:00,12:50,3351,NH929,全日空航空
HangZhou,Tokyo,08:10,12:30,3227,HO3139,吉祥航空
Tokyo,HangZhou,13:55,21:25,2683,HU440,海南航空
ShangHai,Tokyo,17:35,21:30,2381,CA157,中国国航

思路逻辑图：

Python代码

1、定义成本函数

#!/usr/bin/env python
# -*- coding: utf-8 -*-
#@Time: 2019-11-11 14:43
#@Author: gaoll

import time
import random
import math

#组团旅行,行程优化问题
'''
为people中各人员,寻找航班。每个人的航班可选阈来自航班表,其中去程的可选域为出发城市为该人员所在城市,到达城市为Tokyo；返程的解可选域出发城市和目的城市正好相反。
'''
#出行人员及其出发城市
people = [('Gao','HangZhou'),('Zhang','ChengDu'),('Liu','WuHan'),('Song','BeiJing'),('He','ShangHai'),('Niu','ShangHai'),('Li','HangZhou'),('Zhao','xiaMen'),('Qian','GuangZhou'),('Zeng','ShenZhen')]

#东京羽田国际机场
destination = 'Tokyo'

#加载飞机行程、价格数据
schedule_file = 'schedule.txt'  #HangZhou,Tokyo,09:55,13:50,3188,CA929,中国国航
flights = {}

#航班数据在这里就是可行解域
for line in open(schedule_file,'r'):
	if len(line)<1:
		continue
	origin,dest,depart_time,arrive_time,price,flight_info1,flight_info2 = line.strip().split(',')
	flights.setdefault((origin,dest),[])  

	flights[(origin,dest)].append((depart_time,arrive_time,int(price))) #构建每个成员去程、返程的航班可行解域

#比如为Gao选一个去程飞机,即在flights[('HangZhou','Tokyo')]中选一个出发时间、票价飞机行程。所以flights[('HangZhou','Tokyo')]构成Gao的去程行程的解行域。
#为Gao选一个返程飞机,即在flights[('Tokyo','HangZhou')]中选一个出发时间、票价飞机行程。flights[('Tokyo','HangZhou')]构成Gao的返程行程的解行域。
#这里用数字序列来指代航程,如Gao的去程选择第K个航班,即航班为flights[('HangZhou','Tokyo')]中第K个位置元素,flights[('HangZhou','Tokyo')][k]。
#这样Gao的去程航班的解域,就从flights[('HangZhou','Tokyo')]列表,变成了一个数字区间[0,n],其中n=len(flights[('HangZhou','Tokyo')])。
#由上,则sol=[0,0,0,0,0,0,0,0,0,0]即为行程问题的一组解。它表示每个人的去程、返程都乘坐各自出发城市到东京的第一班飞机。这样问题就转成了寻找sol这样的一组最优解。

#定义成本函数
'''
1、考查哪些变量可以衡量行程好坏
（1）机票价格。所有航班的总票价。如果有土豪不关心票价,可以用加权平均,把他对机票价格的权重调低或设为0。
（2）等待时间。10人要在机场聚齐才租车出发,所以,去程时先到的人在机场等待后到达的时间。还有返程时,一同出发去机场还车,飞机有先起飞、后起飞,等待飞机起飞的时间也要算上。机场等待时间按每小时60块。
（3）飞行时长。每个人在飞机上花费的总时间。飞行时长超过8小时后,每小时折算成30块人民币的时间成本。
（4）出发时间。早于早上6点出发或者晚于晚上24点出发的航班不是好的航班,影响睡眠,要折算成成本。每小时折算成80块人民币的健康成本。
（5）汽车租用费用。集体租车,必须在约定的时刻还车,（假设还车时间就是到达机场时间）如果还车晚了1-2小时会产生一天的租车费用,也要算入成本中。租车一天300块。
2、想办法把影响变量组合成一个值
 与机票价格统一,将集中变量都折算成价钱来定义成本函数
'''
def gethours(t):
	x = t.split(':')
	return float(x[0]) + float(float(x[1])/60.0)

def schedulecost(sol):
	totalprice = 0 
	latestarrival = 0 
	earliestdepart = 24*60

	#根据这个行程解,计算行程成本
	for d in range(len(people)):
		name = people[d][0]
		origin = people[d][1]
				
		k1 = sol[2*d] #sol中第 2*d个位置表示people[d]的去程选择的是第几个航班。 destination ='Tokyo'为全局变量,最开始做了设置。
		outbound = flights[(origin,destination)][k1] 
		k2 = sol[2*d+1] #sol中第 2*d+1个位置表示people[d]的返程选择的是第几个航班
		returnf = flights[(destination,origin)][k2]

		outbound_depart_time = outbound[0] #该成员去程出发时间
		outbound_arrive_time = outbound[1] #此成员去程到达机场的时间
		outbound_price = outbound[2]   #去程机票
		return_depart_time = returnf[0] #此成员返程飞机起飞的时间
		return_arrive_time = returnf[1] #该成员返程飞机到达时间
		returnf_price = returnf[2] #返程机票

		#成本累加上该成员的机票价格
		totalprice += outbound_price
		totalprice += returnf_price  
		
		#该成员去程、返程在飞机上的总时间
		flight_time = (gethours(outbound_arrive_time) - gethours(outbound_depart_time)) +(gethours(return_arrive_time) - gethours(return_depart_time))
		
		totalprice +=max(flight_time-8,0)*30 #飞行时长超过8小时后,每小时按30块成本算

		if gethours(outbound_depart_time) < gethours('06:00'):
			totalprice += (gethours('06:00') - gethours(outbound_depart_time))*80 #健康成本

		#记录最晚或最早到达时间,这决定了大家在机场等待最后到达的人,或者等待飞机起飞的时间
		if latestarrival < gethours(outbound_arrive_time):
			latestarrival = gethours(outbound_arrive_time)

		if earliestdepart > gethours(return_depart_time):
			earliestdepart = gethours(return_depart_time)

	#计算成员的总共等待时间,计入成本
	for d in range(len(people)):
		name = people[d][0]
		origin = people[d][1]

		k1 = sol[2*d] #sol中第 2*d个位置表示people[d]的去程选择的是第几个航班。 destination ='Tokyo'为全局变量,最开始做了设置。
		outbound = flights[(origin,destination)][k1] 
		k2 = sol[2*d+1] #sol中第 2*d+1个位置表示people[d]的返程选择的是第几个航班
		returnf = flights[(destination,origin)][k2]

		outbound_arrive_time = outbound[1] #此成员去程到达机场的时间
		return_depart_time = returnf[0] #此成员返程飞机起飞的时间

		totalprice += (latestarrival - gethours(outbound_arrive_time))*60 #去程时,这个成员在机场等待最后一个人到达的时间成本
		totalprice += (gethours(return_depart_time) - earliestdepart)*60 #返程时,这个成员在机场等待飞机起飞的时间成本

	#是否要多加一天的租车成本
	if earliestdepart > latestarrival:
		totalprice +=300

	return totalprice

#因为把行程题解转成了[0,1,2,..]这样的列表来表示行程。但它并不方便直接解读。我们需要把优化后的解再用方便直观的行程表展示出来。
def printschedule(sol):
	for d in range(len(people)):
		name = people[d][0]
		origin = people[d][1]

		k1 = sol[2*d] #sol中第 2*d个位置表示people[d]的去程选择的是第几个航班。 destination ='Tokyo'为全局变量,最开始做了设置。
		outbound = flights[(origin,destination)][k1] 
		k2 = sol[2*d+1] #sol中第 2*d+1个位置表示people[d]的返程选择的是第几个航班
		returnf = flights[(destination,origin)][k2]

		outbound_depart_time = outbound[0] #该成员去程出发时间
		outbound_arrive_time = outbound[1] #此成员去程到达机场的时间
		outbound_price = outbound[2]   #去程机票
		return_depart_time = returnf[0] #此成员返程飞机起飞的时间
		return_arrive_time = returnf[1] #该成员返程飞机到达时间
		returnf_price = returnf[2] #返程机票

		print('%s,%s,%5s-%5s,¥%4s,%5s-%5s,¥%4s'%(name,origin,outbound_depart_time,outbound_arrive_time,outbound_price,return_depart_time,return_arrive_time,returnf_price))

2、优化算法

#!/usr/bin/env python
# -*- coding: utf-8 -*-
#@Time: 2019-11-11 16:43
#@Author: gaoll

import time
import random
import math
import numpy as np

#优化算法

#随机搜索 --创建1000个随机解
#domain是解的取值区域，它是由二元组构成的列表,[(0,9),(0,8),(0,10),..]。第k个位置上的二元组(m,n)指定了第k个变量的最小值、最大值为m、n，即取值区间范围。
def randomoptimize(domain,costf):
	best = 99999999999
	best_sol = None
	for i in range(1000):
		#创建一个随机解
		r = [random.randint(domain[k][0],domain[k][1]) for k in range(len(domain))] 

		#计算这个随机解的成本
		cost = costf(r)

		#此随机解与目前为止的最优解进行比较
		if cost < best:
			best = cost
			best_sol = r	

	return best_sol

#爬山法 --从一个随机解开始，在其临近的解集中寻找最好的题解。
def hillclimb(domain,costf):
	#创建一个随机解
	sol = [random.randint(domain[k][0],domain[k][1]) for k in range(len(domain))]

	best = costf(sol)
	#主循环，为了防止找到局部最优解时，爬山法的迭代次数太大，我这里设了迭代次数不超过1000次
	t = 0 
	while t<1000:
		#创建相邻解的列表,相邻解含义：设第k个位置的解为i, 该解的取值区间为domain[k],即[domain[k][0], domain[k][1]]这个区间的数据。
		#那个这个区间中与i相邻的数值有，i-1,i+1。如果i-1或者i+1溢出了domin[k],则去掉不算入内。
		#于是由第k个位置得到的相邻值与解中其他位置的值组成新的解。对其他位置的解同样操作，就得到一组相邻解，从中选出最优解。
		t+=1
		neighbors = []
		for j in range(len(domain)):
			i = sol[j]  #sol当前解中第j个位置的值
			min = domain[j][0] #第j个位置可行域的最小值
			max = domain[j][1] #第j个位置可行域的最大值

			if i-1 >=min:
				new_sol = sol[0:j] + [i-1] +sol[j+1:] #用i-1替换第j个位置的值i,构成一个新的解，即相邻解
				neighbors.append(new_sol) #将此相邻解加入相邻解列表中

			if i+1 <=max:
				new_sol = sol[0:j] + [i+1] + sol[j+1:]
				neighbors.append(new_sol)

		#寻找相邻解中的最优解
		current = best
		for i in range(len(neighbors)):
			sol_tmp = neighbors[i]
			cost = costf(sol_tmp)
			if cost < best:
				best = cost
				sol = sol_tmp
		if best == current:
			break  #如果这组相邻解中并没有使成本函数变少，那么退出循环，已经找到了一个局部最优解

	return sol

#模拟退火算法
#以一个随机解开始，随机选择当前解中的一个位置，朝着某个方向变化。如果新解的成本值更低，新的题解会成为当前解，如果新解的成本值更高，则由概率决定会不会成为当前解。
#概率的计算公式：P=math.power(e,-(highcost-lowcost)/t), t为当前的温度。设初始温度T，冷却系数cool，迭代n次后，当前的温度t = T*power(cool,n)
#若P>random.random()则接受更低的解作为新解。随着退火过程的进行，算法越来越不可能接受较差的解。
def annealingoptimize(domain,costf,T=10000.0,cool=0.95,step=1):
	#随机初始值
	sol = [random.randint(domain[k][0],domain[k][1]) for k in range(len(domain))]
	ea = costf(sol)

	while T >0.1:
		#选择一个索引值，即解中的某个位置
		i = random.randint(0,len(domain)-1)

		#选择一个改变索引值的方向
		dir = random.randint(-step,step)

		new_sol = sol[:] 
		new_sol[i] += dir  #新解是对原解改变一个位置上的值形成

		min = domain[i][0] #第i个索引位置可行解最小值
		max = domain[i][1] #第i个索引	位置可行解最大值
		if sol[i] + dir < min:
			new_sol[i] = min
		if sol[i] + dir > max:
			new_sol[i] = max

		#计算并比较成本
		eb = costf(new_sol)
		P = np.power(math.e,-(eb-ea)/T)
		if (eb<ea or P >random.random()):
			sol = new_sol

		#降低温度
		T = T*cool

	return sol



#遗传算法
#随机生成一组解称为种群。根据种群的成本函数，选取成本最小的TopN解称为精英解。对精英解进行修改形成一个全新的解组，称为下一代。
#对精英解修改的方法有：变异(对某个位置的值增加或减小)、交叉(截取两个精英解的各一部分组成新解)
#运行几代后，选取最后一代cost最小的作为最优解
def geneticoptimize(domain,costf,popsize=50,step=1,mutprob=0.2,elite=0.2,maxiter=100):
	#变异操作
	def mutate(sol):
		i = random.randint(0,len(domain)-1)
		min = domain[i][0] #第i个索引位置可行解最小值
		max = domain[i][1] #第i个索引	位置可行解最大值
		new_sol = sol[:]
		if random.random() >0.5 and sol[i] - step >min: 
			new_sol[i] = new_sol[i] - step
		elif sol[i] + step < max:
			new_sol[i] = new_sol[i] + step
		return new_sol

	#交叉操作
	def crossover(sol1,sol2):
		i = random.randint(1,len(domain)-2)
		new_sol = sol1[0:i] + sol2[i:]
		return new_sol

	#构造初始种群,大小popsize=50
	pop = []
	for i in range(popsize):
		sol = [random.randint(domain[k][0],domain[k][1]) for k in range(len(domain))]
		pop.append(sol)

	#选择精英，比例elite=0.2
	topelite = int(np.floor(popsize*elite))

	#主循环,maxiter最多进行多少代
	for i in range(0,maxiter):
		scores = [(costf(sol),sol) for sol in pop]
		scores.sort()

		pop = [scores[k][1] for k in range(0,topelite)] #选出精英解作为新的种群
		#由精英种群经过变异、交叉填充足种群大小popsize。变异和交叉操作的比例由参数 mutprob决定
		while len(pop) < popsize:
			if random.random() < mutprob:
				#变异
				c = random.randint(0,topelite-1) #这里pop是随着新解的加入而增加的，不能用len(pop),会选出非精英解.
				mutate_c = mutate(pop[c]) #变异操作
				pop.append(mutate_c)
			else:
				#交叉
				k1 = random.randint(0,topelite-1)
				k2 = random.randint(0,topelite-1)
				cross_sol = crossover(pop[k1],pop[k2]) #交叉操作
				pop.append(cross_sol)
		#print(scores[0][0])
	best_sol = pop[0]
	return best_sol

3、运行

if __name__ == "__main__":
	import optimization
	from optimization import randomoptimize
	from optimization import hillclimb
	from optimization import annealingoptimize
	from optimization import geneticoptimize
	domain = [(0,9)]*(len(people)*2)

	#调用随机搜索算法
	sol1 = randomoptimize(domain,schedulecost)
	print(schedulecost(sol1))
	printschedule(sol1)

	#调用爬山法
	sol2 = hillclimb(domain,schedulecost)
	print(schedulecost(sol2))
	printschedule(sol2)

	#调用模拟退火算法
	sol3 = annealingoptimize(domain,schedulecost)
	print(schedulecost(sol3))
	printschedule(sol3)

	#调用遗传算法
	sol4 = geneticoptimize(domain,schedulecost)
	print(schedulecost(sol4))
	printschedule(sol4)

4、展示这四种算法的一个结果：

#展示这四种算法的一个结果：
cost: 60764.0
schedule:
Gao,HangZhou,13:51-11:71,¥2581,19:19-17:15,¥3112
Zhang,ChengDu,11:45-14:16,¥2542,14:34-17:22,¥3163
Liu,WuHan,16:30-15:24,¥2892,11:07-10:02,¥2532
Song,BeiJing,18:54-16:78,¥2718,12:36-10:53,¥2749
He,ShangHai,17:35-21:30,¥2381,08:06-09:18,¥2743
Niu,ShangHai,15:50-19:40,¥1865,08:06-09:18,¥2743
Li,HangZhou,13:51-11:71,¥2581,08:20-10:05,¥2833
Zhao,xiaMen,23:56-21:33,¥2632,17:22-14:03,¥2921
Qian,GuangZhou,07:38-10:15,¥2994,13:19-13:15,¥3464
Zeng,ShenZhen,23:11-21:07,¥3218,08:27-09:55,¥2679

 

Done