吃遍中國(guó)最地道的100道美食（權(quán)威榜單）的最佳路線（TSP問(wèn)題）

時(shí)間：2023-04-04 02:20:42來(lái)源：food欄目：餐飲美食新聞閱讀：

前言

大家好，我是吃貨，會(huì)寫代碼的那種。

中國(guó)最地道的美食，我相信榜單很多。

為什么中餐在國(guó)際上不夠有名，因?yàn)槊朗程嗔?，很難選出能讓大家信服的最最代表性的，以至于最后勝出并盛行的竟是“炒面”。

但是，刀削面，小面，蘭州拉面，擔(dān)擔(dān)面，油潑面。。。哪個(gè)會(huì)服？

中國(guó)農(nóng)業(yè)部推薦了中國(guó)最地道的100道鄉(xiāng)間美食，我是信服的。

美食鏈接：http://www.moa.gov.cn/ztzl/zgnmfsj/xwzx/201809/t20180913_6157258.htm

除港澳臺(tái)的31省市的地道美食地圖如下：

（圖片來(lái)自人民日?qǐng)?bào)）

作為一個(gè)吃貨，肯定要想吃完這些美食。

作為一個(gè)會(huì)寫代碼的吃貨，就要思考一個(gè)問(wèn)題: 我們要游歷全國(guó)，吃遍這些美食，最優(yōu)的路線是什么？

這就是一個(gè)（吃貨）旅行商問(wèn)題TSP。

TSP問(wèn)題

TSP問(wèn)題是英文縮寫是Travelling Tasting Salesman Problem（旅行商吃貨問(wèn)題）。

旅行商問(wèn)題就是“給定一系列城市和每對(duì)城市之間的距離，求解訪問(wèn)每一座城市一次并回到起始城市的最短回路”。

已知TSP算法最壞情況下的時(shí)間復(fù)雜度隨著城市數(shù)量的增多而成超多項(xiàng)式（可能是指數(shù)）級(jí)別增長(zhǎng)。

一般來(lái)說(shuō)，對(duì)于上萬(wàn)個(gè)城市數(shù)量的TSP問(wèn)題，常見(jiàn)的求解器都是可以搞定的。

TSP應(yīng)該在工業(yè)上除了明顯的物流場(chǎng)景，還有一些有趣的場(chǎng)景：比如PCB電路板的焊接點(diǎn)順序。

我們需要分解一下這個(gè)有趣問(wèn)題的任務(wù) （源碼和數(shù)據(jù)見(jiàn)文末）：

如何獲取城市以及美食清單？如何獲取城市距離？這里需要分情況，坐飛機(jī)去還是開(kāi)車？如何求解這個(gè)旅行商問(wèn)題，得到最后的路線圖？如何畫出路線圖？

獲取美食清單

獲取數(shù)據(jù)，大家想到的第一方法一定是爬蟲(chóng)，農(nóng)業(yè)部的頁(yè)面是htm，屬于入門爬蟲(chóng)級(jí)別。但是不要過(guò)于依賴爬蟲(chóng)技術(shù)。因?yàn)橛行r(shí)候?qū)懘a是效率最差的，尤其是當(dāng)你爬美食網(wǎng)址，容易走神。

htm網(wǎng)址，其實(shí)只要右鍵保存，所有的圖片自動(dòng)就會(huì)給你保存到文件夾了。還需要爬嗎？

無(wú)非是文字需要解析而已，我們要從網(wǎng)頁(yè)里面提取到美食的清單和對(duì)應(yīng)的城市以及圖片。

準(zhǔn)備城市數(shù)據(jù)

人民網(wǎng)推薦的榜單很全面的覆蓋了除港澳臺(tái)的所有省市，我們可以使用分享的數(shù)據(jù)集來(lái)獲取經(jīng)緯度信息。

file = ../data/china_cities.csv all_cities_df = pd.read_csv(file,encoding=gb18030) main_cities = all_cities_df[all_cities_df.sort_values(by=行政代碼)[行政代碼].astype(str).str.find(0100)>0] main_cities.columns = [id,name,lon,lat] main_cities[code] = main_cities[name].str[0:2] meal_city_df = pd.merge(left = meal_df,right=main_cities,on=code,how=left) meal_city_df.head()

我們需要對(duì)數(shù)據(jù)進(jìn)行聚合，因?yàn)槊總€(gè)城市可能有好幾道名菜，我們需要對(duì)城市進(jìn)行g(shù)roup，然后針對(duì)不同的列進(jìn)行不同的統(tǒng)計(jì)聚合。

比如：經(jīng)緯度信息，取第一個(gè)值即可id：我們可以計(jì)數(shù)，用于統(tǒng)計(jì)每個(gè)城市的榜上有名的菜個(gè)數(shù)。

另外，我們可以將每個(gè)城市的菜清單整理成List，這樣可以一起顯示。

meal_city_df_count = meal_city_df.groupby(code).agg({name:first,lon:first,lat:first,id:count}) meal_list = meal_city_df.groupby(name)[meal].apply(list).to_frame().reset_index() meal_list.columns = [name,meal_list] meal_list[meal_list] = meal_list[meal_list].str.join(",") meal_city_df_count = pd.merge(left = meal_city_df_count,right = meal_list,on=name) meal_city_df_count.tail()

我們可以初步看一下這100道名菜的分布。

求解TSP問(wèn)題

TSP作為經(jīng)典的優(yōu)化問(wèn)題，算法已經(jīng)很成熟了，這里我們直接調(diào)用OR-Tools進(jìn)行求解。

OR-Tools求解TSP問(wèn)題需要兩步。

第一步：準(zhǔn)備好城市距離數(shù)據(jù)

我們提供每?jī)蓚€(gè)城市的距離，這里有31個(gè)省市，最后就形成[31x31]的矩陣。

我們已經(jīng)獲取了每個(gè)城市經(jīng)緯度信息，計(jì)算飛行距離就比較容易了。經(jīng)緯度能確定一個(gè)點(diǎn)，兩個(gè)點(diǎn)之間的距離是一個(gè)很好的量化指標(biāo)。距離其實(shí)有很多計(jì)算方式，我們默認(rèn)的“距離”一般指的是歐式距離，對(duì)于小范圍的距離計(jì)算問(wèn)題，比如一個(gè)城市內(nèi)，可能歐式距離也可以求得比較理想的結(jié)果。但是對(duì)于大范圍距離，就會(huì)造成很大的誤差。因?yàn)椋旱厍蚴菆A的，球面上的兩點(diǎn)不適用歐式距離。飛行距離就近似于地球大圓距離，可以采用haversine_distances進(jìn)行計(jì)算。

還有一個(gè)數(shù)據(jù)就是定義出發(fā)城市，這里我隨手定義一個(gè)內(nèi)蒙古。

為什么可以隨手呢？

因?yàn)樽詈蠼鉀Q的路徑是一個(gè)閉環(huán)，從閉環(huán)的任意一點(diǎn)出發(fā)都可以。

from sklearn.metrics.pairwise import haversine_distances from math import radians def geo_distance(p1,p2): dis = haversine_distances([[radians(_) for _ in p1], [radians(_) for _ in p2]])[0][1]* 6371000/1000 # multiply by Earth radius to get kilometers return dis from ortools.constraint_solver import routing_enums_pb2 from ortools.constraint_solver import pywrapcp meal_city_df_count[loc] = meal_city_df_count.apply(lambda x: list([x[lon], x[lat]]),axis=1) xv,yv= np.meshgrid(meal_city_df_count[loc], meal_city_df_count[loc],indexing=ij) data = {} data[distance_matrix] = [[round(geo_distance(xv[i][j],yv[i][j])) for i in range(xv.shape[0])] for j in range(xv.shape[1])] data[num_vehicles] = 1 data[depot] = 2

第二步：原封不動(dòng)的復(fù)制官方定義的函數(shù)

OR-Tools的案例很好，以至于這些函數(shù)都不用改一行代碼，只需要傳入data，得到solution變量即可。

這里對(duì)函數(shù)進(jìn)行說(shuō)明：

RoutingIndexManager 函數(shù)就是創(chuàng)建路徑索引，無(wú)非就是給每個(gè)城市按照順序編個(gè)號(hào)碼而已。RoutingModel：創(chuàng)建路由模型，也就是實(shí)例化一個(gè)模型distance_callback：告訴模型如何計(jì)算使用城市的距離信息，采用RegisterTransitCallback函數(shù)傳入距離信息到模型。這樣模型就可以計(jì)算總距離了SetArcCostEvaluatorOfAllVehicles：計(jì)算其它成本，這里我們默認(rèn)距離就是成本SolveWithParameters：開(kāi)始計(jì)算，很快就出結(jié)果了，畢竟才31個(gè)城市。# Create the routing index manager. manager = pywrapcp.RoutingIndexManager(len(data[distance_matrix]), data[num_vehicles], data[depot]) # Create Routing Model. my_routing = pywrapcp.RoutingModel(manager) def distance_callback(from_index, to_index): """Returns the distance between the two nodes.""" # Convert from routing variable Index to distance matrix NodeIndex. from_node = manager.IndexToNode(from_index) to_node = manager.IndexToNode(to_index) return data[distance_matrix][from_node][to_node] transit_callback_index = my_routing.RegisterTransitCallback(distance_callback) # Define cost of each arc. my_routing.SetArcCostEvaluatorOfAllVehicles(transit_callback_index) # Setting first solution heuristic. search_parameters = pywrapcp.DefaultRoutingSearchParameters() search_parameters.first_solution_strategy = ( routing_enums_pb2.FirstSolutionStrategy.PATH_CHEAPEST_ARC) # Solve the problem. solution = my_routing.SolveWithParameters(search_parameters) # Print solution on console. if solution: print_solution(manager, my_routing, solution)

可以看到總的飛行距離10954公里，路線也已經(jīng)規(guī)劃好了。

Objective: 10954 km Route for vehicle 0: 2 -> 10 -> 17 -> 16 -> 3 -> 6 -> 9 -> 8 -> 14 -> 26 -> 4 -> 30 -> 0 -> 18 -> 23 -> 15 -> 20 -> 21 -> 11 -> 19 -> 12 -> 25 -> 27 -> 5 -> 1 -> 13 -> 24 -> 29 -> 22 -> 7 -> 28 -> 2

展示路線

這里我采用的plotly，畫出路線的基本思路就是將相鄰的城市依次連接起來(lái)。從圖中我們可以看出：

多數(shù)情況下，路徑就是從一個(gè)省到鄰省。從江蘇開(kāi)始后直飛東三省，之后再飛回上海。另外一個(gè)出人意料的地方就是從云南直飛新疆，而不是西藏，這點(diǎn)與我們的直覺(jué)不太吻合。真相就是：我們采用的是大圓距離，而且地圖上會(huì)有麥卡托投影的視覺(jué)誤差。

自駕版本

如果想要自駕過(guò)去，我們同樣可以采用這樣的方式進(jìn)行求解。

自駕版本的計(jì)算要稍微復(fù)雜一些，因?yàn)樾枰@取城市之間的實(shí)時(shí)路線規(guī)劃，然后進(jìn)行求解。

一般我會(huì)采用openstreetmap來(lái)進(jìn)行規(guī)劃，因?yàn)樗腁PI是免費(fèi)的。

我們這次再看一下自駕版本結(jié)果：

吃完老北京炸醬面（北京）之后，我們應(yīng)該直奔東北了，取嘗嘗殺豬菜了（黑龍江）吃飯鹽水鴨（江蘇）之后，之后的確是應(yīng)該去嘗蟹殼黃了(上海)這次擔(dān)擔(dān)面（四川）離風(fēng)干牛羊肉（西藏）更近了。bjective: 19258 km Route for vehicle 0: 2 -> 3 -> 30 -> 4 -> 26 -> 6 -> 9 -> 8 -> 14 -> 0 -> 18 -> 23 -> 15 -> 20 -> 21 -> 11 -> 19 -> 12 -> 1 -> 25 -> 27 -> 5 -> 24 -> 13 -> 29 -> 22 -> 7 -> 28 -> 17 -> 16 -> 10 -> 2

可以看到總的飛行距離19258 公里，接近飛行距離的兩倍，不過(guò)自駕的風(fēng)景應(yīng)該是飛行風(fēng)景的N倍吧。