Python基礎圖論演算法

摘要： 本來下學期在學姐的強力安利之下選了演算法這個課，但是取消了，於是在家打發時間看了edX上的一個法國人講的演算法網課。主要講一些基礎的圖論演算法，結合一個老鼠走迷宮的問題，用Python 寫寫程式。 image.png 個人感覺這個課對Python...

本來下學期在學姐的強力安利之下選了演算法這個課，但是取消了，於是在家打發時間看了edX上的一個法國人講的演算法網課。主要講一些基礎的圖論演算法，結合一個老鼠走迷宮的問題，用Python 寫寫程式。

個人感覺這個課對Python新手還是很友好的，因為他讓你寫的程式每一行都給你了註釋。不想花錢搞verified，又覺得過段時間記錄丟失有點虧，於是總結一下，留給自己看。

後來在回來的飛機上看了看5047 Intelligence System的PPT，發現這些演算法全在在5047的Week 2 lecture裡。我看了一個多星期的內容居然這個lecturer要用兩個小時加50頁ppt講完。所以基本上只是提到了概念，沒有做深入探討 (畢竟這篇文章加上程式4000字)。

一， 先說演算法

這裡涉及到的演算法有：DFS, BFS, Dijkstra, Brute force, Back tracking, Greedy。

正常的一個圖，由邊（Edge）和節點（Vertex）組成，可以用矩陣表示（Corresponding Matrix）。如果不考慮邊的長度，兩個點之間連線算1，不連線算0，這種叫adjacency matrix。這個矩陣橫座標代表從V1到Vn，縱座標也是V1到Vn，關於主對角線對稱（往返），並且對角線上全是0（自己到自己）。

Tree是在圖的基礎上得到的，要求cycle-free 。隨便定義一個點為root，這樣其他點就有了parent和children 的關係。因為cycle-free，一個children只能對應一個parent，這樣從root到每一個節點的path就被唯一確定。有點像9135裡網路的Routing Table。根據圖我們可以通過BFS或者DFS得到圖的Spinning Tree。這個過程叫做graph traversal ，這個結果不唯一。

DFS指的是deep first search，深度優先，就是一條路走到黑，沒路了再回頭。每連線一個節點，下一個節點選相鄰的沒連線的，如果相鄰的節點都已經連線過了，就返回，到上一步裡選一個相鄰的沒連線的，直到所有節點都連線到這個tree上。DFS產生的tree的特點是分支儘可能的少。

BFS指的是breadth first search，寬度優先，就是先把周圍的連線上，再去連線更遠的點。從一個節點開始，先連線和他相鄰的節點，再連線隔一個hop的到和他相鄰的，再連線隔兩個hop的到隔一個hop的節點上，直到所有點都連線到這個tree上。BFS生成的tree有一個特點，root到任意一個節點經過的hop都是最少的情況。

如果每個邊的長度都是1，那麼BFS產生的tree就是從root到任意一個點的最短距離的情況。但是如果考慮比較一般的情況，邊的長度是不一定相等的正數，這時候要找到一個點的最短距離就要用Dijkstra 演算法。

先考慮最簡單的情況，從v1到v3的距離。v2和 v3之間的距離如果大於2，顯然最短距離是v1到v3，如果小於2，最短距離是v1-v2-v3。這個最短路徑需要進一步的比較才能確定。但是，如果我們要找的是從v1到v2的最短距離，因為所有邊的距離都是正數，很顯然不管v2到v3的距離是多少，v1直接到v2一定是最短距離，而不是經過v3。所以我們可以得出一個結論，一個節點和他相鄰的節點裡距離最短的節點之間的路徑，就是這兩個節點之間的最短路徑。

對於比較複雜的情況，使用迭代的思想。從root開始，找相鄰節點中距離最短的，然後把已經確定的最短路徑上的節點和root節點當作一個整體，更新這個大的root節點到相鄰節點之間的距離為到實際root節點的距離，root整體節點不斷擴大，不斷迭代，直到包含所有的節點，就能找到每一個節點距離都最近的tree。

例子如上圖所示，能夠確定v1到v3一定是最短距離，於是把v1和v3當作一個整體，一個大的root節點，然後整體節點到相鄰節點v2的最短距離從5更新為4。這時候包含v1和v3的大的root節點相鄰的路徑到v2是最短（到v2更新為4，到v7更新為5），下一個被包進來的節點就是v2，以此類推，實現迭代，最終找到root到每個節點的最短路徑。

以上考慮的都是從一個節點到其他節點的距離問題。下面考慮一個更復雜的問題Traveling Salesman Problem，簡稱TSP。這個問題是在找從一個點開始，把圖上的所有節點走一遍的最短距離。

這個問題圖上如果有N個節點，那麼就要列舉計算(N-1)!種情況，這個數字非常大，所以這類問題沒有一個有效的演算法去解決。所謂Brute force就是計算出所有可能的路徑，找出最小的。Back tracking在此基礎上稍稍優化，記錄一個當前找到的最短路徑，在計算下一條路徑的時候如果算到一半發現總的距離已經比這個最短路徑大了，就不再進行這條路上後面的計算，返回，算下一條路。原理很簡單，只是稍作優化。

Greedy演算法就是說每一步都是區域性最優（locally optimal），對整體最有結果進行逼近，不一定能找到整體的最優結果，但是高效。放到這個TSP問題上，就是我在v1的時候走到v4是最短距離，在v4的時候走到v3是最短距離，從而認為v1-v4-v3-v2是最短距離，不再進行其他判斷。這個結果顯然不是最優的，但是高效。其他問題上也可以使用greedy演算法，比如要湊齊一定的錢用最少的硬幣。

二，再來寫程式

以上內容的原理其實很好理解，寫程式就比較難了。

在Python裡表示一個圖，使用的是字典（Dict）的巢狀。先把圖上點分好座標，作為外層字典的key，每個點的value值包含內層字典，表示其他點到這個點的距離。

maze_graph = {
(0,0): {(0,1):1,(1,0):1},
(0,1): {(0,2):1,(0,0):1},
(1,0): {(1,1):1,(0,0):1},
(1,1): {(1,2):1,(1,0):1},
(0,2): {(0,1):1,(1,2):1},
(1,2): {(0,2):1,(1,1):1}
}

寫DFS需要用到一個數據結構，LIFO Stack（last in first out）。最後放到這個stack裡（Push）的元素最先出來（Pop）。這個東西雖然Python裡有，可以直接用，但也可以用List實現。從root節點開始，每次先pop一個節點，當作走過的，再把相鄰的沒走過的節點push進這個stack，先push進去的會最後pop出來，這樣同一級的節點就會最後走。直到所有的節點都被pop，整個圖就按照DFS的方式走了一遍。

LIFO_list = list()
def LIFO_push(LIFO_list,element):
LIFO_list.append(element)
def LIFO_pop(LIFO_list):
return LIFO_list.pop(-1)

這裡定義一個list，用來存已經走過的節點，並且定義方法來加加上已經走過的節點，判斷一個點是否已經走過，即是否在這個list中。

def add_to_explored_vertices(explored_vertices,vertex):
explored_vertices.append(vertex)

def is_explored(explored_vertices,vertex):
return vertex in list(explored_vertices)

以下是實現這個DFS演算法的方法，引數是圖和root節點，返回值是已經走過的節點和每個節點的parent的字典。Python比Java在程式設計上的一個很明顯的改進就是支援一個方法可以有很多個返回值，讀取這些返回值時只要把變數按照順序排好就可以。

def DFS(maze_graph, initial_vertex) :
# explored vertices list
explored_vertices = list()
#LIFO stack
queuing_structure = list()
#Parent Dictionary
parent_dict = dict()
LIFO_push(queuing_structure,(initial_vertex,None)) # push the initial vertex
while len(queuing_structure) > 0: # while queuing_structure is not empty:
# current_vertex,parent = queuing_structure.pop()
# if the current vertex is not explored
# add current_vertex to explored vertices
# use parent_dict to map the parent of the current vertex
# for each neighbor of the current vertex in the maze graph:
# if neighbor is not explored:
# push the tuple (neighbor,current_vertex) to the queuing_structure
current_vertex,parent =queuing_structure.pop()
if not is_explored(explored_vertices,current_vertex):
add_to_explored_vertices(explored_vertices,current_vertex)
parent_dict.update({current_vertex:parent})
for neighbor in maze_graph[current_vertex]:
if not is_explored(explored_vertices,neighbor):
LIFO_push(queuing_structure,(neighbor,current_vertex))
return explored_vertices,parent_dict

寫BFS需要用到另外一個數據結構，FIFO Queue（first in first out），先放進去的元素後出來。這個和DFS的原理差不多，也是先pop出來一個節點，當作走過的，再把相鄰的沒走過的節點push進這個queue，因為先放進去的會先出來，這樣就會先走同一級的節點，直到所有的節點都走遍。

def BFS(maze_graph, initial_vertex) :
# explored vertices list
explored_vertices = list()
#LIFO stack
queuing_structure = list()
#Parent Dictionary
parent_dict = dict()
FIFO_push(queuing_structure,(initial_vertex,None)) # push the initial vertex to the while
len(queuing_structure) > 0: # while queuing_structure is not empty:
# current_vertex,parent = queuing_structure.pop()
# if the current vertex is not explored
# add current_vertex to explored vertices
# use parent_dict to map the parent of the current vertex
# for each neighbor of the current vertex in the maze graph:
# if neighbor is not explored:
# push the tuple (neighbor,current_vertex) to the queuing_structure
current_vertex,parent = FIFO_pop(queuing_structure)
if current_vertex not in explored_vertices:
explored_vertices.append(current_vertex)
parent_dict.update({current_vertex:parent})
FIFO_push(queuing_structure,(current_vertex,parent))
for neighbor in maze_graph[current_vertex]:
if neighbor not in explored_vertices:
FIFO_push(queuing_structure,(neighbor,current_vertex))
return explored_vertices,parent_dict
一下這個方法用來利用DFS或者BFS 的結果找到圖上從一個點到任意一個點的路徑。
def create_walk_from_parents(parent_dict,initial_vertex,target_vertex):
relist = list()
relist.append(target_vertex)
while(initial_vertex != target_vertex):
key = parent_dict[target_vertex]
relist.append(key)
target_vertex = key
relist.reverse()
if relist != []:
del relist[0]
return relist

initial_vertex = (0,0)
target_vertex = (0,2)
explored_vertices,parent_dict = DFS(maze_graph,initial_vertex)
route = create_walk_from_parents(parent_dict,initial_vertex,target_vertex)

在寫Dijkstra 演算法的時候，需要用到另一個數據結構，min-heap。Min-heap可以理解為一個字典，裡面都是（key, value）一對一對的。這個min-heap主要用兩個操作，add-or-replace和remove，add-or-replace就是在新增一個新的(key, value)的時候，如果這個key在min-heap裡已經有了，value只有小於min-heap中的，才能更新value。這樣可以保證value始終是最小值，並且key值不重複。因為Dijkstra需要找到最短路徑，我們就可以把節點作為key，距離作為value，不斷更新這個value，就能找到最短路徑。

但是在做這個演算法的時候，因為需要得到routing table，節點的parent也要加上。所以我們需要定義一個triplet，(vertex, weight, parent)，在比較value的時候，只比較weight。

先定義一個heap的pop方法，得到第一個元素。

# heap_pop function returns the first element of the list implementing the heap, providing the heap is not empty
def heap_pop(heap):
if heap != []:
vertex,weight,parent = heap.pop(0)
return (vertex, weight, parent)
else:
raise

因為要比較weight的值，即第二個元素，最後在sort中也要按照weight排序。定義一個方法，返回tuple中的第二個元素

def fun(s):
return s[1]

之後定義這個add_or_replace方法。

def heap_add_or_replace(heap,triplet):
i = 0
while i < len(heap):
if (heap[i])[0] == triplet[0]:
if (heap[i])[1] > triplet[1]:
heap[i] = triplet
heap.sort(key = fun)
return
i = i+1
heap.append(triplet)
heap.sort(key = fun)

def is_explored(explored_vertices,vertex):
return vertex in explored_vertices
def add_to_explored_vertices(explored_vertices,vertex):
explored_vertices.append(vertex)
def Dijkstra(maze_graph,initial_vertex):
# Variable storing the exploredled vertices vertexes not to go there again
explored_vertices = list()
# Stack of vertexes
heap = list()
#Parent dictionary
parent_dict = dict()
# Distances dictionary
distances = dict()
# First call
initial_vertex = (initial_vertex, 0, initial_vertex)#vertex to visit, distance from
heap_add_or_replace(heap,initial_vertex)
while len(heap) > 0:
# get the triplet (vertex, distance, parent) with the smallest distance from heap # if the vertex of the triplet is not explored:
# map the vertex to its corresponding parent
# add vertex to explored vertices.
# set distance from inital_vertex to vertex
# for each unexplored neighbor i of the vertex, connected through an edge of # add (i, distance + wi, vertex) to the heap
#
vertex,distance,parent = heap_pop(heap)
if not is_explored(explored_vertices,vertex):
parent_dict.update({vertex:parent})
add_to_explored_vertices(explored_vertices,vertex)
distances[vertex] = distance
#distances.update({vertex:(distance+d)})
for neighbor in maze_graph[vertex]:
if neighbor not in explored_vertices:
wi = (maze_graph[vertex])[neighbor]
heap_add_or_replace(heap,(neighbor,wi + distance,vertex))

return explored_vertices, parent_dict, distances

最後把圖的尺寸作為變數，畫一個計算時間隨圖的尺寸變數變化的關係圖。

from imports import maze
import time
import math
%matplotlib inline
import matplotlib.pyplot as plt
maze_size=[]
time_to_execute=[]
maze_shape=[(1,2),(2,3),(3,3),(5,3),(6,7),(8,9),(11,10),(12,15),(16,18),(19,22),(21,23)]
for i in range(len(maze_shape)):
length,width=maze_shape[i]
maze_graph=maze.generate_maze(length,width,0,True,False,0.5,5,"",0)
start=time.time()
Dijkstra(maze_graph[3],(0,0))
time_to_execute.append(time.time()-start)
maze_size.append(length*width)
plt.plot(maze_size,time_to_execute,color="blue")

可以大致看出來程式執行時間基本和節點的數量是一個二次曲線的關係。這裡涉及到演算法複雜度的討論。Dijkstra 演算法的複雜度主要依賴與最小堆的實現方法。我們這個最小堆算每一個節點的距離的時候線性搜尋了堆中的所有元素，所以複雜度是O(N^2). 經過對這個堆的演算法優化，讓他不是線性地去找最小值，可以進行優化，減少演算法複雜度，實現O(NlogN)。這個比較難懂，不做探討。

解決TSP問題的演算法需要計算每一條路徑，總路徑的條數是(N-1)!，因此Brute force演算法的複雜度就是O((N-1)!)。這個演算法的複雜度大於N的多項式，不能在多項式時間內求解，這種問題就叫做NP問題（Non-deterministic Polynomial）。

Python裡連線兩個list可以用“+”，也可以用append()，區別是“+”不會改變原來的list()。

def create_vertices_meta_graph(piece_of_cheese, player_location):
#
return piece_of_cheese + [player_location]
#
定義方法，把side的距離新增到adjacency矩陣裡。
import utils
def create_edge_weight_maze_graph(maze_graph,vertices):
adjacency_matrix={}
# for each initial_vertex in vertices:
# considere this vertex as source vertex
# use this source vertex and maze_graph to browse the graph with dijkstra algorithm
# for each vertex in vertices:
# use adjacency_matrix to store distances between source vertex and each vertex # remember to not store the distance from the source vertex to itself.
#
for iv in vertices:
source = iv
explored_vertices,parent_dict,dist_dict = utils.Dijkstra(maze_graph,source)
adjacency_matrix.update({source:{}})
for v in vertices:
if v != source:
adjacency_matrix[source].update({v:dist_dict[v]})
return adjacency_matrix

定義Brute force演算法的方法，找到最短路徑和最短距離。

def auxbf(current_walk,best_walk,adjacency_matrix,vertices,current_distance,best_distance):
# if length of current walk is larger than the order of the graph:
# if moreover the current distance is shorter than the best distance obtained so # update the value of the best distance
# update the corresponding best walk
# otherwise:
# for any vertex in the graph:
# if the vertex is not in the current walk:
# obtain potential values for best walk and best distance by recursively # if potential distance is shorter than best distance:
# update the value of the best distance
# update the corresponding best walk
#
if len(current_walk) > len(adjacency_matrix)-1:
if current_distance < best_distance:
best_distance = current_distance
best_walk = current_walk
else:
for v in adjacency_matrix:
if v not in current_walk:
walk_temp,dis_temp =
auxbf(current_walk+[v],best_walk,adjacency_matrix,vertices,current_distance+
adjacency_matrix[current_walk[-1]][v],best_distance)
if dis_temp < best_distance:
best_distance = dis_temp
best_walk = walk_temp

return best_walk,best_distance

呼叫上述三個方法進行brute force TSP問題的計算。

def bruteforceTSP(maze_graph,pieces_of_cheese,player_location):
# first we compute the vertices of the meta_graph:
vertices=create_vertices_meta_graph(pieces_of_cheese, player_location)
# then we create the adjacency matrix of the meta graph
adjacency_matrix = create_edge_weight_maze_graph(maze_graph,vertices)
# now we can start defining our variables
# current_distance is the length of the walk for the current exploration branch
current_distance = 0
# current_walk is a container for the current exploration branch
current_walk = [player_location]
# best_distance is an indicator of the shortest walk found so far
best_distance = float('inf')
# best_walk is a container for the corresponding walk
best_walk = []
# we start the exploration:
best_walk, best_distance =auxbf(current_walk,best_walk,adjacency_matrix,
pieces_of_cheese,current_distance,best_distance)
return best_walk, best_distance

在brute force的基礎上稍加改動，加入一個和當前最短距離的判斷，就能定義back tracking 方法

def auxbt(current_walk,best_walk,adjacency_matrix,vertices,current_distance,best_distance):
# First we test if the current walk have gone through all vertices
# if that is the case, we compare the current distance to the best distance
# and in the case it is better we update the best distance and the best walk
# if the current_walk is not finished we see if the current distance is lower than # if that is the case, for each possible vertex not explored,
# we add it and call ourself recursively
##
YOUR CODE HERE
#
if (len(current_walk) == len(vertices)):
if (current_distance < best_distance):
best_distance=current_distance
best_walk=current_walk
else:
if current_distance < best_distance:
for v in adjacency_matrix:
if v not in current_walk:
walk_temp,dis_temp =auxbt(current_walk+[v],best_walk,adjacency_matrix,
adjacency_matrix[current_walk[-1]][v],best_distance)
if dis_temp < best_distance:
best_distance = dis_temp
best_walk = walk_temp
return best_walk,best_distance

用time 計算兩個演算法對不同尺寸的圖的計算時間，並繪圖比較。

import time
%matplotlib inline
import matplotlib.pyplot as plt
maze_shape = [(3,3),(5,6),(6,7),(9,10),(11,13)]
poc = [2,3,4,5,6]
bruteforce_time = list()
backtracking_time = list()
maze_size=list()
for i in range(len(maze_shape)):
width,height = maze_shape[i]
_,_,_,maze_graph = maze.generate_maze(width,height,0,True,False,0.5,5,"",0)
pieces_of_cheese,player1_location,_ = maze.generate_pieces_of_cheese(
poc[i], width, height, False, None, None, False)
start = time.time()
bruteforceTSP(maze_graph,pieces_of_cheese,player1_location)
bruteforce_time.append(time.time()-start)
start = time.time()
backtrackingTSP(maze_graph,pieces_of_cheese,player1_location)
backtracking_time.append(time.time() - start)
maze_size.append(width*height)
plt.plot(maze_size,bruteforce_time,color="blue")
plt.plot(maze_size,backtracking_time,color="green")

分別是3x3的圖裡2個目標點，5x6的圖裡3個，6x7的圖裡4個，9x10的圖裡5個，11x13的圖裡6個的情況。可以看出back tracking在brute force的基礎上有所優化，尤其是目標點多了之後。

最後是greedy演算法，這個背景是老鼠要找迷宮中的所有乳酪，就每次都去找最近的，實現對總距離最短的近似實現。每一次的最近距離使用Dijkstra 演算法找。

def choose_target(distances,pieces_of_cheese):
distance = float('inf')
next_target = None
# for each vertex in pieces of cheese
#compare its distance from player location with distance variable value.
#update distance variable with the smallest value
#keep the nearest vertex from player location
# returns the nearest vertex to be next target
for vertex in pieces_of_cheese:
if distances[vertex] < distance or distance is float('inf') :
distance = distances[vertex]
next_target = vertex
return next_target

import utils
movements = list()
def movements_greedy_algorithm(maze_graph,pieces_of_cheese,player_location):
# get distances using Dijkstra algorithm
# get next_target using choose_target function
# use A_to_B function to get a list of movements that should be done to reach the nearest piece of cheese
# from player location
explored_vertices,parent_dict,distance_dict = utils.Dijkstra(maze_graph,player_location)
next_target = choose_target(distance_dict,pieces_of_cheese)
movements = utils.A_to_B(maze_graph,player_location,next_target,parent_dict)
return movements,next_target

最後，畫一個程式執行時間和圖的節點數的關係圖。