优化算法--以Python实现（2）

优化算法最大的困难之处在于，其一，将问题以合理的方式形式化的表达出来；即《集体智慧编程》中所说：当处理类似问题时，我们有必要明确潜在的题解将如何表达。其二，构造合理的成本函数。



一种通用的方法就是以数字序列来表示题解。如前面的旅游问题中，我们以题解来表达没人的航班次。



但在某些情况下，以数字化的方式来表达题解不是很直接，需要进行一些变通。



例如，宿舍安排问题。

有5间宿舍，每间分割为2个隔间，可安排两人居住。由10名学生来竞争居住，每个学生尤其首选和次选。但如下的数据结构所示，并不是每个学生都能安排进他们的首选宿舍。Bacchus的安排就很成问题。

# chr05 -- dorm.py
import random
import math
#-------------------------------------------------------------------------------
# 宿舍，每个宿舍有两个隔间
dorms = ['Zeus', 'Athena', 'Hercules', 'Bacchus', 'Pluto']
# 学生，首选，次选
prefs=[
    ('Toby',    ('Bacchus', 'Hercules')),
    ('Steve',   ('Zeus', 'Pluto')),
    ('Karen',   ('Athena', 'Zeus')),
    ('Sarah',   ('Zeus', 'Pluto')),
    ('Dave',    ('Athena', 'Bacchus')), 
    ('Jeff',    ('Hercules', 'Pluto')), 
    ('Fred',    ('Pluto', 'Athena')), 
    ('Suzie',   ('Bacchus', 'Hercules')), 
    ('Laura',   ('Bacchus', 'Hercules')), 
    ('James',   ('Hercules', 'Athena')) ]

最困难的是以合理的形式化方式将问题表达为可以优化算法求解的形式。这关系到优化算法的可行性。以下domain的定义，即体现了这一点。

domain = [ (0, (len(dorms)*2)-i-1) for i in range( len(dorms)*2 ) ]

在该问题中，是有5间宿舍，而每个宿舍可以有2个隔间，因此将形式化的序列定义为元组的列表较为合理。实际上是把一个可安排的宿舍位看做是一个槽。这样的槽又需满足一些限制条件，我们的目标是要以最为优化的方式，在尽量满足所有人的偏好的同时，填满这10个槽。因此，关于5个宿舍的安排，最终变成了10个槽位的竞争。

可以由以下这个打印宿舍安排的函数来更好的了解domain和槽的概念。

#-------------------------------------------------------------------------------
# 打印宿舍安排序列
def printsolution(vec):
    slots = []
    # 为每个宿舍建两个槽
    for i in range( len(dorms) ): slots += [i, i]
    # 遍历每一个学生的安置情况
    for i in range(len(vec)):
        x = int(vec[i])
        # 从剩余槽中选择
        dorm = dorms[slots[x]]
        # 输出学生及其被分配的宿舍
            print( prefs[i][0], dorm )
        # 删除该槽
        del slots[x]

slots实际是0-4的一列表。vec[i]为第i个学生在slot中的位置，实际上槽是固定的，变化的是对学生的分配。

成本计算函数与之类似：

#-------------------------------------------------------------------------------
# 成本函数
def dormcost(vec):
    cost = 0
    # 建立一个槽序列
    slots = []
    for i in range( len(dorms) ): slots += [i, i]
    # 遍历每个学生
    for i in range(len(vec)):
        x = int(vec[i])
        dorm = dorms[slots[x]]
        pref = prefs[i][1]
        # 首选成本为0，次选成本为1
        if pref[0] == dorm: cost +=0
        elif pref[1] == dorm: cost += 1
        else: cost += 3
        del slots[x]
    return cost

定义成本为：如果是首选宿舍，则成本为0；次选为1，否则为3.容易看出，最优成本是0，即每个人都安排了其首选宿舍。虽然我们可以从以上定义的数据结构中看出这是不可能的，但知道最优成本为0可以让我们更容易的评价当前安排序列的优劣，以及距最优序列的远近如何。



成本函数定义是关键：

我们可以使用在《优化算法1》中定义的各优化函数。例如，我们调用randomoptimize函数。

#----------------------------------------------------------------
# 随机搜索
def randomoptimize(domain, costfunc):
    best = 999999999
    bestr = None
    # 随机进行1000次猜测
    for i in range(1000):
        # 创建一个随机群
        r = [ random.randint(domain[idx][0], domain[idx][1])
              for idx in range(len(domain)) ]
        # 计算成本
        cost = costfunc(r)
        # 与目前为止的最优解比较
        if cost < best:
            best = cost
            bestr = r
    return bestr

我们的测试代码：

#================================================================
# test: dorm.py
#================================================================
if __name__ == '__main__':
    print('=========================')
    print('Test: dorm.py')
    print('=========================')
    domain = [ (0, (len(dorms)*2)-i-1)
               for i in range( len(dorms)*2 ) ]
    li =optimization.randomoptimize(domain, dormcost)
    print( dormcost(li) )
    print( li )

可以看出，实际上，randomoptimize不过是根据用户给定的domain，产生一系列的随机序列，使用用户给的成本函数计算成本，以此来判断一个好的序列。

因此，成本函数是关键。

关键在于，如何定义合适的domain，和如何在成本函数中，使用randomoptimize随机产生的序列进行有效的成本计算，如何将randomoptimize的随机序列变得有意义（当然，randomoptimize产生的随机序列是可通过domain来控制的）。



具体到本例中，如果把randomoptimize产生的随机序列理解做：“学生在10个槽位中所站的那个槽位的编号”就大错特错了。因为randomoptimize产生的随机序列中的元素是会有重复值的。一个槽位不可能安排两个学生来住……



要理解randomoptimize随机序列的含义，我们可以通过成本函数dormcost对入内的序列的使用方式来看。在上面的代码中，可以看出，dormcost并不是将入内的序列当做槽位的标号来使用，而是如前所说，入内的vec中，vec[i]为第i个学生在slot中的位置，而slot才是真真的槽位，宿舍号从0-4，每个宿舍两个槽位，实际上槽是固定的，变化的是对学生的分配。slot中存放的是真正的宿舍号，而我们所要求的是一个序列，该序列表示某人在宿舍槽位序列（固定，为：[0, 0, 1, 1, 2, 2, 3, 3, 4, 4]）中所占的位置，如若入内vec为： [9, 0, 2, 0, 3, 4, 0, 2, 0, 0]，则表示第一个人Toby的宿舍将被安排在slot【9】所决定的宿舍，即：4号宿舍，Bacchus。



现在才终于真正理解生成domain的方法了：

domain = [ (0, (len(dorms)*2)-i-1)
               for i in range( len(dorms)*2 ) ]

domain是用于在randomoptimize中产生随机序列的，即dormcost的入内参数，即，宿舍槽位slot的编号。因为共有10个槽位，那么槽位的编号应该是0-9,所以我们所生成的domain的最大值应该是9.而在dormcost中，使用掉一个槽位，我们就将该槽位从slot中删除，因此整个slot就少了一个元素，slot的最大下标随之变化，由之前的9变到了8.所以，在生成的domain【2】中，最大值是8，其范围是0-8，以此类推，故domain是：

[(0, 9), (0, 8), (0, 7), (0, 6), (0, 5), (0, 4), (0, 3), (0, 2), (0, 1), (0, 0)]



所以，此前的优化算法都是可以通用的，但在使用之前应该定义好所需要的domain以及成本函数。并不是随便调用就可以了。



测试了下，随机搜索算法最大能达到成本为5，而“遗传算法--以爬山法产生初始种群”可以到达2（太强了。。。）：

2
[4, 0, 1, 0, 0, 3, 3, 1, 1, 0]
Toby Hercules
Steve Zeus
Karen Athena
Sarah Zeus
Dave Athena
Jeff Pluto
Fred Pluto
Suzie Bacchus
Laura Bacchus
James Hercules

这是在看《集体智慧编程》时的收获，但《集体智慧编程》并没有对该算法做更多的解释。当写完这个算法后，我一度百思不得其解，陷入函数randomoptimize产生的是随机序列，是有重复的，可能一个序列中，出现1的次数远远超过2，而1号宿舍最多只有2个隔间……书上并没有告诉我，作者是怎么得到该算法的，也并没有做如我上述的讲解，书中往往是：“让我们在XX.py中加入函数XX”……并没有告诉我怎么想出这个函数，该函数与其他函数背后有联系吗？是什么？彼此之间怎么协作？（实际上，我看到的很多发表在专业期刊上的论文也是这个风格：让我们这样来解决：……然后咔咔咔就出来了，完全不知为何如此。。。）。就如刘未鹏博客所说一样：我们要的不是相对论，而是诞生相对论的那个大脑。我们要的不是金蛋，而是下金蛋的那只鸡。





合理构造题解序列

再看网络可视化问题。

我们要把社交网站上的人际关系绘制成一个脉络清晰的网络图，比如，该网络图节点分布合理，交叉线少。对前一问题，我们可以利用质点弹簧算法解决。但质点弹簧算法对问题二无能为力。

于是我们采用优化算法来尽量的减少交叉线。



我们所要做的全部，就是定义一个良好的成本函数，并尝试让它返回的值尽可能的小。很显然的成本函数是通过计算彼此相交的线。

那么重点就在于怎么表示题解了。



由于我们要计算的彼此相交的线，而这些线是通过将“有社交关系的两个人”表示成两点并将之相连而形成的，而我们所要做的是，将这些点合理安排使之直接线段相交最少。即是，我们要求的是一组关于这些点的坐标，这些坐标使得相互直接的有效连接（有社交关系的两点）彼此相交最少。

因此，我们将题解表示做关于坐标点的集合。

sol = [120, 200, 250, 125 ...]

表示：Charlie位于（120,200）。。。

这样成本函数只要对彼此相交的线段进行计数即可。



定义的数据结构如：

#chr05 -- socialnetwork.py
import math
people=['Charlie','Augustus','Veruca','Violet','Mike','Joe','Willy','Miranda']
links=[('Augustus', 'Willy'), 
       ('Mike', 'Joe'), 
       ('Miranda', 'Mike'), 
       ('Violet', 'Augustus'), 
       ('Miranda', 'Willy'), 
       ('Charlie', 'Mike'), 
       ('Veruca', 'Joe'), 
       ('Miranda', 'Augustus'), 
       ('Willy', 'Augustus'), 
       ('Joe', 'Charlie'), 
       ('Veruca', 'Augustus'), 
       ('Miranda', 'Joe')]

定义成本函数为：

（原来关于两线相交也是有公式的。。。基本思路是，计算线条的分数值，若两线分数值位于0-1之间，则其彼此相交）

#------------------------------------------------------------------------
# 成本函数
def crosscount(v):
    # 将数字序列v转化为：people:(x,y)形式的字典
    loc = dict( [(people[i], (v[i*2], v[i*2+1]))
                 for i in range(0, len(people))] )
    total = 0
    
    # 遍历每一对连线
    for i in range(len(links)):
        for j in range(i+1, len(links)):
            # 获取坐标位置
            (x1, y1), (x2, y2) = loc[links[i][0]], loc[links[i][1]]
            (x3, y3), (x4, y4) = loc[links[j][0]], loc[links[j][1]]
            den = (y4 - y3)*(x2 - x1)-(x4 - x3)*(y2 - y1)
            # 如果两线平行，则den == 0
            if den == 0: continue
            # ua,ub为两线的分数值
            ua = ((x4-x3)*(y1-y3)-(y4-y3)*(x1-x3))/den
            ub = ((x2-x1)*(y1-y3)-(y2-y1)*(x1-x3))/den
            # 若两线的分数值介于0-1之间，则相交
            if ua > 0 and ua < 1 and ub > 0 and ub < 1:
                total += 1
    return total

#=========================================================================
# test: socialnetwork.py
#=========================================================================
if __name__ == '__main__':
    domain = [(10, 370)]*(len(people) *2)
    sol = optimization.randomoptimize (domain, crosscount)
    print (sol)
    print (crosscount(sol))

因此，在利用优化算法时，对题解的数字序列表示时，应灵活变通。



我想我有必要看看数据挖掘方面的Book。。。



By Kewing.