Tackling the TSP with a genetic algorithm

This is my first go at genetic algorithms, using PyGame for the visual representation. Is there a way to improve the speed of the code? Click to place a point, press to start looking for the shortest route, or press c to clear all.

import pygame, time, random

#Size of next generation, change number at the end if needed
resize_rate = 1.005
#Controls crossover appearence
cross_rate = 0.7
#Controls number of chromosomes created at the start, exponential
exponent = 4

black = (0, 0, 0)
white = (255, 255, 255)
red = (255, 0, 0)

genes = []

def draw_background(screen):
screen.fill(white)
def draw_connection(screen, pointlist):
pygame.draw.lines(screen, black, True, pointlist, 1)
def draw_best_connection(screen, pointlist):
pygame.draw.lines(screen, red, True, pointlist, 2)

def create_random_chromosomes(genes):
temp_chromosomes = []
chromosomes = []
for i in range(len(genes)**exponent):
genes_copy = genes[:]
random.shuffle(genes_copy)
temp_chromosomes.append(genes_copy)
for element in temp_chromosomes:
if element not in chromosomes:
chromosomes.append(element)
print(str(len(chromosomes)) + " random chromosomes created")
return chromosomes

def calculate_distance(chromosome):
distances = []
total_distance = 0
for i in range(len(chromosome)):
if i < (len(chromosome)-1):
distance_A_B = (((chromosome[i+1][0]-chromosome[i][0])**2)+(chromosome[i+1][1]-chromosome[i][1])**2)**0.5
distances.append(distance_A_B)
else:
distance_A_B = (((chromosome[i][0]-chromosome[0][0])**2)+(chromosome[i][1]-chromosome[0][1])**2)**0.5
distances.append(distance_A_B)
for i in range(len(distances)):
total_distance += distances[i]

def calculate_fitness(chromosomes):
fitness_scores = []
for i in range(len(chromosomes)):
fitness = 0
fitness = 1/calculate_distance(chromosomes[i])
fitness_scores.append(fitness)
return fitness_scores

def build_chromosome_fitness_dict(chromosomes, fitness_scores):
chrom_fit_dict = {}
for i in range(len(chromosomes)):
chrom_fit_dict[fitness_scores[i]] = chromosomes[i]
return chrom_fit_dict

def roulette_selection(chromosomes, fitnesses, new_size):
sum_fitness = sum(fitnesses)
rel_fitness = [fitness/sum_fitness for fitness in fitnesses]
#Generate probability intervals for each chromosome
probabilities = [sum(rel_fitness[:i+1]) for i in range(len(rel_fitness))]
#Select chromosomes
selected_chromosomes = []
for n in range(new_size):
r = random.random()
for (i, chromosome) in enumerate(chromosomes):
if r <= probabilities[i]:
selected_chromosomes.append(chromosome)
break
return selected_chromosomes

def find_closest(cities):
dist_a = (((cities[0][0]-cities[1][0])**2)+(cities[0][1]-cities[1][1])**2)**0.5
dist_b = (((cities[0][0]-cities[2][0])**2)+(cities[0][1]-cities[2][1])**2)**0.5
if dist_a < dist_b:
return cities[1]
else:
return cities[2]

#Greedy crossover selects the first city of one parent,
#compares the cities leaving that city in both parents,
#and chooses the closer one to extend the tour.
#If one city has already appeared in the tour, we choose the other city.
#If both cities have already appeared, we randomly select a non-selected city.
def crossover(chromosomes, rate):
new_chromosomes = []
#len -1 because index out of range error might occur otherwise
for i in range(len(chromosomes)-1):
#Crossover T/F
if random.random() < rate:
#Take 2 chromosomes
a = chromosomes[i]
b = chromosomes[i+1]
#Randomly select crossover location
cross_location = random.randint(1, len(chromosomes[i]))
#Create parts
c = a[:cross_location]
for i in range((len(a)-len(c))):
#City to start from
start_pos = c[-1]
#2 most plausible targets
target_0 = a[cross_location+i]
target_1 = b[cross_location+i]
if target_0 in c:
if target_1 in c:
not_yet_selected = []
for city in genes:
if city not in c:
not_yet_selected.append(city)
c.append(not_yet_selected[random.randint(0,(len(not_yet_selected)-1))])
else:
c.append(target_1)
elif target_1 in c:
c.append(target_0)
else:
cities = [start_pos, target_0, target_1]
closest = find_closest(cities)
c.append(closest)
new_chromosomes.append(c)
else:
new_chromosomes.append(chromosomes[i])
return new_chromosomes

pygame.init()
screen = pygame.display.set_mode((640,340))
clock = pygame.time.Clock()
done = False
first = True
while done == False:
for event in pygame.event.get():
if event.type == pygame.QUIT:
done = True
if event.type == pygame.KEYDOWN:
if event.key == pygame.K_UP:
if genes != []:
#Start Genetic Algorithm
print('Creating initial solutions')
chromosomes = create_random_chromosomes(genes)
#Draw most possible solutions
for i in range(len(chromosomes)):
draw_connection(screen, chromosomes[i])
#Main Logic Loop
print('Starting calculations')
while len(chromosomes) > 10:
print(str(len(chromosomes)) + ' chromosomes left')
new_size = int(len(chromosomes)/resize_rate)
#Calculate fitness scores
fitness_scores = calculate_fitness(chromosomes)
#Build dictionary to preserve order
chrom_fit_dict = build_chromosome_fitness_dict(chromosomes, fitness_scores)
#Roulette selection - select chromosomes for next generation
selected_chromosomes = roulette_selection(chromosomes, fitness_scores, new_size)
#Adopt newly selected chromosomes as population
chromosomes = selected_chromosomes
#Apply crossover mutation
mutated_chromosomes = crossover(chromosomes, cross_rate)
chromosomes = mutated_chromosomes
best_chrom = fitness_scores.index(max(chrom_fit_dict))
draw_best_connection(screen, chromosomes[best_chrom])
if event.key == pygame.K_c:
draw_background(screen)
genes = []
if first == True:
draw_background(screen)
first = False
#Get points
if pygame.mouse.get_pressed() == (1, 0, 0):
mouse_pos = pygame.mouse.get_pos()
genes.append(mouse_pos)
pygame.draw.rect(screen, black, [mouse_pos[0], mouse_pos[1], 2, 2], 1)
genes = list(set(genes))
pygame.display.flip()
clock.tick(30)
pygame.quit()

-
How many vertexes (cities) do you want to process? Do you need a genetic algorithm only, not an exact one? – cat_baxter Jul 3 '12 at 20:54
It isn't meant for real-world application so it doesn't really matter how many cities/points it can process, largest test I ran was with 30 different points, it executed in about 7-8 minutes and was most of the times I ran it 90% - 100% correct. I'm only interested in the genetic algorithm, this was a test project to learn how these work, so an exact result is not requiered. – Daquicker Jul 4 '12 at 12:57
Ok, I see. If you need the exact TSP solution (may be in order to check you genetic algorithm solution) let me know. It solves 30 points TSP in about 0.5 second. – cat_baxter Jul 4 '12 at 15:14
oh, to check solutions that would be interesting indeed. – Daquicker Jul 4 '12 at 17:13

Ok, so here is exact ATSP solver, just in order to check your genetic algorithm solution. It's based on pascal code implementing the branch-and- bound algorithm by Little, Murty, Sweeney, and Karel. It's very very fast (30 points in 0.5 sec!).

'''

Branch and bound Algorithm for Travelling Salesman Problem

based on

http://www.cs.sunysb.edu/~algorith/implement/syslo/distrib/processed/babtsp.p

'''

import time

INF = 100000000

best_cost = 0

def reduce(size, w, row, col, rowred, colred):
rvalue = 0
for i in xrange(size):
temp = INF
for j in xrange(size):
temp = min(temp, w[row[i]][col[j]])
if temp > 0:
for j in xrange(size):
if w[row[i]][col[j]] < INF:
w[row[i]][col[j]] -= temp
rvalue += temp
rowred[i] = temp
for j in xrange(size):
temp = INF
for i in xrange(size):
temp = min(temp, w[row[i]][col[j]])
if temp > 0:
for i in xrange(size):
if w[row[i]][col[j]] < INF:
w[row[i]][col[j]] -= temp
rvalue += temp
colred[j] = temp
return rvalue

def bestEdge(size, w, row, col):
mosti = -INF
ri = 0
ci = 0
for i in xrange(size):
for j in xrange(size):
if not w[row[i]][col[j]]:
minrowwelt = INF
zeroes = 0
for k in xrange(size):
if not w[row[i]][col[k]]:
zeroes += 1
else:
minrowwelt = min(minrowwelt, w[row[i]][col[k]])
if zeroes > 1: minrowwelt = 0
mincolwelt = INF
zeroes = 0
for k in xrange(size):
if not w[row[k]][col[j]]:
zeroes += 1
else:
mincolwelt = min(mincolwelt, w[row[k]][col[j]])
if zeroes > 1: mincolwelt = 0
if minrowwelt + mincolwelt > mosti:
mosti = minrowwelt + mincolwelt
ri = i
ci = j
return mosti, ri, ci

def explore(n, w, edges, cost, row, col, best, fwdptr, backptr):
global best_cost

colred = [0 for _ in xrange(n)]
rowred = [0 for _ in xrange(n)]
size = n - edges
cost += reduce(size, w, row, col, rowred, colred)
if cost < best_cost:
if edges == n - 2:
for i in xrange(n): best[i] = fwdptr[i]
if w[row[0]][col[0]] >= INF:
avoid = 0
else:
avoid = 1
best[row[0]] = col[1 - avoid]
best[row[1]] = col[avoid]
best_cost = cost
else:
mostv, rv, cv = bestEdge(size, w, row, col)
lowerbound = cost + mostv
fwdptr[row[rv]] = col[cv]
backptr[col[cv]] = row[rv]
last = col[cv]
while fwdptr[last] != INF: last = fwdptr[last]
first = row[rv]
while backptr[first] != INF: first = backptr[first]
colrowval = w[last][first]
w[last][first] = INF
newcol = [INF for _ in xrange(size)]
newrow = [INF for _ in xrange(size)]
for i in xrange(rv): newrow[i] = row[i]
for i in xrange(rv, size - 1): newrow[i] = row[i + 1]
for i in xrange(cv): newcol[i] = col[i]
for i in xrange(cv, size - 1): newcol[i] = col[i + 1]
explore(n, w, edges + 1, cost, newrow, newcol, best, fwdptr, backptr)
w[last][first] = colrowval
backptr[col[cv]] = INF
fwdptr[row[rv]] = INF
if lowerbound < best_cost:
w[row[rv]][col[cv]] = INF
explore(n, w, edges, cost, row, col, best, fwdptr, backptr)
w[row[rv]][col[cv]] = 0

for i in xrange(size):
for j in xrange(size):
w[row[i]][col[j]] = w[row[i]][col[j]] + rowred[i] + colred[j]

def atsp(w):
global best_cost
size = len(w)
col = [i for i in xrange(size)]
row = [i for i in xrange(size)]
best = [0 for _ in xrange(size)]
route = [0 for _ in xrange(size)]
fwdptr = [INF for _ in xrange(size)]
backptr = [INF for _ in xrange(size)]
best_cost = INF

explore(size, w, 0, 0, row, col, best, fwdptr, backptr)

index = 0
for i in xrange(size):
route[i] = index
index = best[index]
index = []
cost = 0

for i in xrange(size):
if i != size - 1:
src = route[i]
dst = route[i + 1]
else:
src = route[i]
dst = 0
cost += w[src][dst]
index.append([src, dst])
return cost, index

def main():

m = [
[INF, 514, 230, 92, 172, 201, 320, 205, 329, 285, 404, 128, 352, 357, 190, 269, 484, 237, 567, 209, 413, 404,
333, 480, 371, 175, 470, 281, 423, 328],
[452, INF, 359, 410, 450, 401, 315, 293, 244, 200, 275, 446, 469, 143, 278, 362, 399, 344, 544, 644, 335, 381,
248, 457, 626, 375, 551, 258, 397, 305],
[414, 360, INF, 90, 170, 42, 90, 46, 282, 238, 186, 126, 320, 198, 31, 267, 437, 190, 358, 623, 246, 345, 286,
271, 267, 16, 326, 222, 421, 119],
[324, 427, 269, INF, 80, 109, 359, 113, 263, 233, 312, 36, 260, 265, 98, 177, 418, 257, 475, 533, 321, 312, 281,
388, 279, 83, 378, 189, 331, 236],
[244, 382, 312, 226, INF, 335, 402, 246, 197, 153, 334, 262, 486, 466, 231, 315, 352, 297, 670, 453, 241, 538,
201, 583, 505, 309, 298, 415, 350, 431],
[448, 318, 355, 216, 268, INF, 445, 111, 240, 196, 288, 252, 385, 156, 96, 358, 395, 148, 691, 640, 348, 528,
244, 455, 373, 299, 566, 405, 393, 452],
[330, 437, 126, 216, 280, 119, INF, 123, 206, 315, 312, 252, 348, 275, 108, 324, 361, 267, 268, 539, 372, 255,
242, 181, 393, 93, 236, 132, 478, 29],
[495, 557, 521, 539, 313, 563, 469, INF, 406, 444, 543, 575, 274, 700, 440, 524, 442, 506, 729, 633, 554, 716,
442, 642, 788, 537, 611, 593, 641, 490],
[208, 256, 115, 168, 206, 157, 205, 49, INF, 148, 137, 204, 323, 313, 34, 118, 155, 100, 473, 417, 197, 460, 36,
386, 382, 131, 441, 337, 272, 234],
[252, 229, 159, 212, 250, 201, 249, 93, 44, INF, 181, 248, 367, 357, 78, 162, 199, 144, 517, 444, 241, 398, 48,
430, 426, 175, 485, 275, 197, 278],
[534, 568, 363, 453, 411, 250, 237, 254, 326, 344, INF, 225, 449, 406, 239, 366, 481, 398, 505, 547, 60, 375,
362, 418, 623, 224, 473, 252, 520, 266],
[445, 507, 334, 389, 469, 332, 393, 336, 413, 369, 276, INF, 224, 488, 321, 141, 392, 321, 659, 583, 336, 496,
417, 572, 601, 306, 629, 373, 295, 420],
[221, 283, 313, 313, 393, 314, 195, 263, 214, 170, 351, 349, INF, 426, 248, 325, 168, 314, 455, 359, 411, 442,
218, 368, 580, 288, 431, 319, 367, 216],
[309, 286, 216, 269, 307, 258, 306, 150, 101, 57, 132, 305, 326, INF, 135, 219, 256, 201, 574, 501, 192, 455,
105, 487, 483, 232, 542, 332, 254, 335],
[416, 554, 484, 398, 172, 507, 484, 15, 369, 325, 506, 434, 289, 638, INF, 487, 457, 469, 717, 625, 413, 554,
373, 359, 677, 481, 470, 431, 522, 478],
[514, 509, 193, 283, 352, 191, 252, 195, 306, 387, 379, 319, 469, 347, 180, INF, 461, 339, 518, 723, 249, 355,
342, 431, 460, 165, 488, 232, 154, 279],
[53, 115, 278, 145, 225, 254, 368, 212, 163, 119, 300, 181, 405, 258, 197, 281, INF, 263, 620, 262, 360, 457,
167, 533, 424, 228, 523, 334, 316, 381],
[300, 277, 207, 68, 148, 177, 297, 141, 92, 48, 229, 104, 328, 333, 126, 210, 247, INF, 543, 492, 289, 380, 96,
456, 347, 151, 446, 257, 245, 304],
[758, 471, 620, 710, 754, 593, 710, 597, 550, 671, 687, 746, 571, 614, 582, 668, 705, 650, INF, 930, 747, 599,
586, 404, 887, 567, 915, 476, 822, 523],
[39, 434, 269, 131, 211, 240, 359, 227, 178, 324, 315, 167, 391, 396, 212, 296, 333, 276, 452, INF, 375, 439,
214, 365, 410, 214, 509, 316, 450, 213],
[526, 513, 443, 496, 534, 485, 533, 377, 328, 284, 441, 165, 389, 641, 362, 306, 483, 428, 698, 487, INF, 535,
332, 611, 710, 459, 769, 412, 460, 459],
[784, 654, 691, 552, 604, 336, 781, 447, 576, 532, 624, 588, 721, 492, 432, 694, 731, 484, 1027, 976, 684, INF,
580, 791, 709, 635, 902, 741, 729, 788],
[435, 220, 576, 527, 607, 569, 450, 513, 464, 420, 495, 563, 689, 363, 498, 582, 619, 564, 718, 396, 555, 601,
INF, 631, 806, 543, 686, 478, 617, 479],
[354, 402, 216, 306, 350, 189, 306, 193, 146, 294, 283, 342, 167, 345, 178, 264, 301, 246, 358, 526, 343, 195,
182, INF, 483, 163, 511, 72, 418, 119],
[311, 288, 218, 79, 159, 188, 308, 152, 103, 59, 240, 115, 113, 344, 137, 221, 258, 11, 554, 472, 300, 391, 107,
462, INF, 162, 457, 268, 256, 315],
[431, 344, 186, 242, 187, 26, 276, 30, 266, 222, 314, 278, 304, 182, 15, 384, 421, 174, 544, 640, 374, 531, 270,
374, 399, INF, 348, 408, 419, 305],
[477, 454, 384, 437, 475, 426, 474, 318, 269, 225, 300, 473, 494, 168, 303, 387, 424, 369, 742, 669, 360, 623,
273, 655, 651, 400, INF, 500, 422, 503],
[282, 330, 144, 234, 280, 186, 234, 123, 74, 222, 211, 270, 366, 342, 108, 192, 229, 174, 286, 491, 271, 123,
110, 199, 411, 160, 470, INF, 346, 47],
[360, 399, 39, 129, 209, 81, 98, 85, 152, 277, 225, 165, 359, 237, 70, 270, 307, 229, 364, 569, 95, 201, 188,
277, 306, 55, 334, 78, INF, 125],
[301, 433, 97, 187, 267, 139, 187, 143, 177, 325, 283, 223, 319, 295, 128, 295, 332, 277, 239, 510, 343, 226,
213, 152, 364, 113, 423, 103, 449, INF]
]

start_time = time.time()
cost, path = atsp(m)
print "Cost = ", cost
print "Path = ", path
print "Time (s)", time.time() - start_time

if __name__ == "__main__":
main()

-