机器学习 — 优化

优化问题

使用随机优化解决写作类问题:存在多种变量的影响,存在许多个可能的解,通过对题解打分,找到一个问题的最优解。

优化的主要思想:

  1. 找到影响结果的因素,比如这里旅行的航班价格、花费时间、租车费用等
  2. 将考虑到的主要因素根据权重组成,计算出总的成本
  3. 利用一定的算法找到成本最小时候的各个因素的取值

关键在于确定题解表示法和成本函数

计算最小成本的方法:

  1. 随机搜索
  2. 爬山法
  3. 模拟退火
  4. 遗传算法

在本节中要找到所有人乘坐航班的最佳班次表,也就是一个航班班次的列表

import time
import math
import random

people = [('Seymour','BOS'),
          ('Franny','DAL'),
          ('Zooey','CAK'),
          ('Walt','MIA'),
          ('Buddy','ORD'),
          ('Les','OMA')]
# Laguardia
destination='LGA'


# 解析schedule.txt为航班详情
flights = {}
for line in file('schedule.txt'):
    origin, dest, depart, arrive, price = line.split(',')
    flights.setdefault((origin, dest), [])
    
    # 将航班详情添加到航班列表中
    flights[(origin, dest)].append((depart, arrive, int(price)))
    
# 计算某一个时间在一天中的分钟数
def getminutes(t):
    x = time.strptime(t, '%H:%M')
    return x[3] * 60 + x[4]

# 一个表格的形式打印每个人来回的航班班次
def printschedule(r):
    """
    根据给定的航班班次列表,将每个人乘坐的航班详情打印成一个表格,直观
    
    Args:
        r:航班班次列表,也就是需要求出的解
    """
    for i in range(len(people)):
        name = people[i][0]
        origin = people[i][1]
        out = flights[(origin, destination)][r[2*i]]
        ret = flights[(destination, origin)][r[2*i+1]]
        print '%10s%10s %5s-%5s $%3s %5s-%5s $%3s' % (name, origin, out[0], out[1], out[2], ret[0], ret[1], ret[2])
# 同一个(origin, dest)航班班次
s = [1, 4, 3, 2, 7, 3, 6, 3, 2, 4, 5, 3]
printschedule(s)
   Seymour       BOS  8:04-10:11 $ 95 12:08-14:05 $142
    Franny       DAL 10:30-14:57 $290  9:49-13:51 $229
     Zooey       CAK 17:08-19:08 $262 10:32-13:16 $139
      Walt       MIA 15:34-18:11 $326 11:08-14:38 $262
     Buddy       ORD  9:42-11:32 $169 12:08-14:47 $231
       Les       OMA 13:37-15:08 $250 11:07-13:24 $171

成本函数

  1. 确定对成本影响的变量
  2. 根据变量计算出总成本

在这里影响成本的变量有:

  1. 每班次航班的价格
  2. 每个人在飞机上花费的总时间
  3. 在机场等待其他成员的时间
  4. 汽车租用时间,归还汽车时间必须早于租车时间,否则需要多付一天的租金

这里总成本使用金钱来衡量,每一分钟价值1元

# 计算总成本
def schedulecost(sol):
    # 总成本
    totalprice = 0
    # 最后到机场
    latestarrival = 0
    # 最早离开机场
    earliestdepart = 24 * 60
    
    # 计算出所有人来往航班价格总和,并得出所有人坐飞机最晚到的时间和最早坐飞机离开的时间
    for d in range(len(sol)/2):
        origin = people[d][1]
        # 得到往返航班
        outbound = flights[(origin, destination)][int(sol[2*d])]
        returnf = flights[(destination, origin)][int(sol[2*d+1])]
        
        # 计算所有往返航班价格之和
        totalprice += outbound[2]
        totalprice += returnf[2]
        
        # 记录最晚到达时间和最早离开时间
        if latestarrival < getminutes(outbound[1]):
            latestarrival = getminutes(outbound[1])
        if earliestdepart > getminutes(returnf[0]):
            earliestdepart = getminutes(returnf[0])
        
    # 每个人必须在机场等到最晚到达的航班的人到齐之后才可以离开机场
    # 每个人也必须和最早离开航班的那个人一起到达机场
    totalwait=0  
    for d in range(len(sol)/2):
        origin = people[d][1]
        outbound = flights[(origin,destination)][int(sol[d])]
        returnf = flights[(destination,origin)][int(sol[d+1])]
        totalwait += latestarrival-getminutes(outbound[1])
        totalwait += getminutes(returnf[0]) - earliestdepart
    
    # 如果离开时间大于到达时间,那么需要多付一天的费用
    if latestarrival < earliestdepart:
        totalprice += 50
        
    # 假设一分钟的价值是1元
    return totalprice + totalwait
print 'total cost: %d' % schedulecost(s)
total cost: 1623

随机搜索

不是一种非常好的优化算法,但是其他优化算法评价的基线

def randomoptimize(domain, costfunc):
    """
    随机搜索算法:
        1. 循环n次
        2. 每次产生一个随机解(这里就是一组航班班次的列表)
    随机算法不能办证找到的是最优解,只有循环次数足够大的时候,才能求出最接近的真实解
    
    Args:
        domain:解的每个值的取值范围
        costfunc:计算成本函数
        
    Returns:
        bestr:循环n次得到随机最优解
    """
    best = 9999999
    bestr = None
    
    for i in range(10000):
        # 创建一个随机解——航班班次列表
        r = [random.randint(domain[i][0], domain[i][1]) for i in range(len(domain))]
        
        # 成本
        cost = costfunc(r)
        
        if cost < best:
            best = cost
            bestr = r
    print bestr    
    return bestr
domain = [(0, 9)] * (len(people)*2)
print domain
s = randomoptimize(domain, schedulecost)
print schedulecost(s)
printschedule(s)
[(0, 9), (0, 9), (0, 9), (0, 9), (0, 9), (0, 9), (0, 9), (0, 9), (0, 9), (0, 9), (0, 9), (0, 9)]
[4, 7, 4, 6, 1, 7, 2, 9, 2, 8, 4, 8]
1322
   Seymour       BOS 12:34-15:02 $109 17:03-18:03 $103
    Franny       DAL 12:19-15:25 $342 15:49-20:10 $497
     Zooey       CAK  8:27-10:45 $139 16:33-18:15 $253
      Walt       MIA  9:15-12:29 $225 20:27-23:42 $169
     Buddy       ORD  9:42-11:32 $169 18:33-20:22 $143
       Les       OMA 12:18-14:56 $172 18:25-20:34 $205
# 爬山法,找到局部最优解
def hillclimb(domain, costfunc):
    """
    爬山法:
        1. 指定一组随机解
        2. 进行不断循环
        3. 在每次循环中,找到该解附近的所有解
        4. 计算所有相邻解中的最优解,作为新的解
        5. 回到3
        6. 直到相邻解中没有更优的解(成本都一样)跳出循环
    爬山法得到的是一组局部最优解
    
    Args:
        domain:解中每个值的取值范围
        costfunc:计算成本的函数
    """
    # 创建一个随机解
    sol = [random.randint(domain[i][0], domain[i][1]) for i in range(len(domain))]
    print domain
    print sol
    print '----------------'
    # 主循环
    while 1:
        # 创建相邻解的列表
        neighbors = []
        for j in range(len(domain)):
            # 在每一个方向上相对于原值偏离一点,作为一个相邻解

            if sol[j] > domain[j][0]:
                neighbors.append(sol[0:j] + [sol[j] - 1] + sol[j+1:])
            if sol[j] < domain[j][1]:
                neighbors.append(sol[0:j] + [sol[j] + 1] + sol[j+1:])
        
        # 在相邻解中寻找最优解
        current = costfunc(sol)
        best = current
        for j in range(len(neighbors)):
            cost = costfunc(neighbors[j])
            if cost < best:
                best = cost
                sol = neighbors[j]
        
        # 如果没有更好的解则跳出循环
        if best == current:
            break
        
        return sol
s = hillclimb(domain, schedulecost)
print schedulecost(s)
printschedule(s)
[(0, 9), (0, 9), (0, 9), (0, 9), (0, 9), (0, 9), (0, 9), (0, 9), (0, 9), (0, 9), (0, 9), (0, 9)]
[8, 1, 1, 0, 8, 9, 8, 0, 8, 9, 2, 2]
----------------
7078
   Seymour       BOS 18:34-19:36 $136  8:23-10:28 $149
    Franny       DAL  6:12-10:22 $230  6:09- 9:49 $414
     Zooey       CAK 18:35-20:28 $204 19:46-21:45 $214
      Walt       MIA 18:23-21:35 $134  6:33- 9:14 $172
     Buddy       ORD 18:48-21:45 $246 19:32-21:25 $160
       Les       OMA  9:15-12:03 $ 99  9:31-11:43 $210
# 模拟退火算法,随机选择一个位置的值进行改变,改变的多少也是随机的
def annelingoptimize(domain, costfunc, T=1000.0, cool=0.95, step=1):
    """
    模拟退火法:
        1. 得到一组随机解
        2. 不断循环,直到温度降到最低(这里是小于0.1)
        3. 每次循环中,随机选择一个方向(列表中任意一个位置的索引值),
            朝某一个方向(正、负)变化一个值(步长正负值范围内的一个随机值)
        4. 计算变化前后的成本,如果成本减少、或者一个随机值小于pow(math.e, -(eb-ea)/T)(这种情况下,
            即使成本没有减少,变化后的列表也有可能被选择作为新的列表,避免了局部最优),
            则将变化后的列表作为新的列表
        5. 降低温度,回到2
        6. 循环结束,返回最后得到的最优解
    Args:
        domain:解中每个值的取值范围
        costfunc:计算成本的函数
        T:初始温度
        cool:退火速度
        step:每次变化的步长
        
    Returns:
        vec:退火之后的最优解
    """
    # 随机初始化值
    vec = [random.randint(domain[j][0], domain[j][0]) for j in range(len(domain))]
    while T > 0.1:
        # 选择一个初始索引值
        i = random.randint(0, len(domain)-1)
        
        # 选择一个改变索引值方向
        direction = random.randint(-step, step)
        # 创建一个代表题解的新列表
        vecb = vec[:]
        vecb[i] += direction
        if vecb[i] < domain[i][0]:
            vecb[i] = domain[i][0]   
        elif vecb[i] > domain[i][1]:
            vecb[i] = domain[i][1]
            
            
        # 计算当前成本和新的成本
        ea = costfunc(vec)
        eb = costfunc(vecb)
        
        # 如果成本减少则选择较小的解
        # 如果成本没有减少,但是也会有机会选择较小的解,这样就避免了局部最优
        if (eb < ea or random.random() < pow(math.e, -(eb-ea)/T)):
            vec = vecb
        
        # 降低温度
        T = T * cool
    return vec
s = annelingoptimize(domain, schedulecost)
print schedulecost(s)
print s
printschedule(s)
2553
[1, 0, 0, 3, 1, 1, 0, 0, 1, 1, 1, 3]
   Seymour       BOS  8:04-10:11 $ 95  6:39- 8:09 $ 86
    Franny       DAL  6:12-10:22 $230 10:51-14:16 $256
     Zooey       CAK  8:27-10:45 $139  8:19-11:16 $122
      Walt       MIA  6:25- 9:30 $335  6:33- 9:14 $172
     Buddy       ORD  8:25-10:34 $157  7:50-10:08 $164
       Les       OMA  7:39-10:24 $219 11:07-13:24 $171
# 遗传算法
def geneticoptimize(domain, costfunc, popsize=50, step=1, mutprob=0.2, elite=0.2, maxiter=100):
    """
    遗传算法:
        1. 随机生成popsize个解组成的列表,作为刚开始的种群
        2. 循环maxiter次,遗传maxiter代
        3. 每次循环中,计算出种群中所有解的成本
        4. 将种群内的解按成本大小排序
        5. 取最优的(成本最小的)一些(elite*popsize决定有多少个可以胜出)胜出者作为新的种群
        6. 从新种群中随机选择一些进行变异或者交叉(修改题解的两种方法)之后添加到新种群,将新种群补齐为popsize个
        7. 回到3
        8. 循环结束后返回最后种群中的最优解
    
    Args:
        popsize:种群大小
        step:每次变异大小
        mutprob:新种群由变异得来的概率
        elite:种群中是优解,且被允许遗传给下一代的部分
        maxiter:需运行多少代
    
    Returns:
        遗传maxiter代之后的最优解
    """
    
    # 变异操作
    def mutate(vec):
        i = random.randint(0, len(domain)-1)
        if random.random() < 0.5 and vec[i] > domain[i][0]:
            return vec[0:i] + [vec[i]-step] + vec[i + 1:]
        elif vec[i] < domain[i][1]:
            return vec[0:i] + [vec[i]+step] + vec[i + 1:]
        else:
            return vec
        
    # 交叉操作
    def crossover(r1, r2):
        i = random.randint(1, len(domain)-2)
        return r1[0:i] + r2[i:]
    
    # 构造种群
    pop = []
    for i in range(popsize):
        vec = [random.randint(domain[i][0], domain[i][1]) for i in range(len(domain))]
        pop.append(vec)
        
    # 每一代中有多少个胜出者
    topelite = int(elite*popsize)
    
    # 主循环
    for i in range(maxiter):
        scores = [(costfunc(v), v) for v in pop]
        scores.sort()
        ranked = [v for (s, v) in scores]
        
        # 选出胜出者
        pop = ranked[0:topelite]
        
        # 添加变异和交叉之后,变成新的种群pop
        while len(pop) < popsize:
            if random.random() < mutprob:
                # 变异:选择胜出者进行变异,然后添加到新种群中
                c = random.randint(0, topelite)
                pop.append(mutate(ranked[c]))
                test = mutate(ranked[c])
            else:
                # 交叉:任选胜出者中的两个进行交叉
                c1 = random.randint(0, topelite)
                c2 = random.randint(0, topelite)
                pop.append(crossover(ranked[c1], ranked[c2]))
        # 打印当前最优解
        print scores[0][0]
        
    return scores[0][1]                
geneticoptimize(domain, schedulecost)
3093
3093
1946
1940
1931
1707
1563
1479
1479
1464
1246
1246
1246
680
680
554
554
554
447
399
399
399
366
366
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247
247





[4, 7, 3, 7, 4, 9, 2, 9, 4, 7, 4, 9]
posted @ 2017-03-27 20:29  lacker  阅读(619)  评论(0编辑  收藏  举报