ga

ga.py

# -*- coding: utf-8 -*-

import numpy as n
import random as r
import pylab as p
import copy

# Parameters ##################################################################

# how many cities?
CITIES = 50
# max number of generations
GENERATIONS = 1000
# total population size
POP_SIZE = 200
# crossover and mutation rate
CROSSOVER_RATE = 0.6
MUTATION_RATE = 0.1
# how many neighbors to change in each cluster?
TO_CHANGE_RATE = 0.1
# how many representatives to vote on? (=K in Kmeans)
NUM_REPRESENTATIVES = 3
# dna: [ ELITE_RATE | OP_RATE | RANDOM ]
ELITE_RATE = 0.2
OP_RATE = 0.6

# Operators ###################################################################

def crossover(i_mom, i_dad):
    """ Crossover by edge recombination. Based on
    http://en.wikipedia.org/wiki/Edge_recombination_operator and Pyevolve's
    Crossovers.G1DListCrossoverEdge method
    """
    g_mom = i_mom.copy()
    g_dad = i_dad.copy()
    sisterl = []
    brotherl = []

    def get_edges(ind):
        edg = {}
        ind_list = ind['dna']
        
        for i in xrange(len(ind_list)):
            a, b = ind_list[i], ind_list[i-1]
            if a not in edg:
                edg[a] = []
            else:
                edg[a].append(b)
            if b not in edg:
                edg[b] = []
            else:
                edg[b].append(a)
        return edg
    
    def merge_edges(edge_a, edge_b):
        edges = {}
        for value, near in edge_a.items():
            for adj in near:
                if (value in edge_b) and (adj in edge_b[value]):
                    edges.setdefault(value, []).append(adj)
        return edges
    
    def get_edges_composite(mom, dad):
        mom_edges = get_edges(mom)
        dad_edges = get_edges(dad)
        return (mom_edges, dad_edges, merge_edges(mom_edges, dad_edges))
    
    mom_edges, dad_edges, merged_edges = get_edges_composite(g_mom, g_dad)

    for c, u in (sisterl, set(g_mom['dna'])), (brotherl, set(g_dad['dna'])):
        curr = None
        for i in xrange(len(g_mom['dna'])):
            curr = r.choice(tuple(u)) if not curr else curr
            c.append(curr)
            u.remove(curr)
            d = [v for v in merged_edges.get(curr, []) if v in u]
            if d:
                curr = r.choice(d)
            else:
                s = [v for v in mom_edges.get(curr, []) if v in u]
                s += [v for v in dad_edges.get(curr, []) if v in u]
                curr = r.choice(s) if s else None
    # garante que sempre haverá um 0 no início do dna
    pos0sister = sisterl.index(0)
    sisterl[pos0sister] = sisterl[0]
    sisterl[0] = 0
    pos0brother = brotherl.index(0)
    brotherl[pos0brother] = brotherl[0]
    brotherl[0] = 0
    
    sister = g_mom.copy()
    brother = g_dad.copy()
    sister['dna'] = sisterl
    brother['dna'] = brotherl
    
    return (sister, brother)

def mutate(guy):
    """ Mutation by sublist reversing """
    inicio = r.choice(range(1,len(guy['dna'])-1))
    fim = r.choice(range(inicio, len(guy['dna'])-1))
    # invertemos a ordem dessa sublista
    aux = guy.copy()
    foo = aux['dna'][inicio:fim+1]
    foo.reverse()
    # trocamos a sublista antiga pela invertida
    aux['dna'][inicio:fim+1] = foo[:]
    return aux

def fitness(guy):
    return 1. / guy['score']

def score(guy):
    """ Scoring based on the sum of distances of all valid paths """
    def dist(c1, c2):
        p1 = cities[c1]
        p2 = cities[c2]
        return n.sqrt((p2[0] - p1[0])**2 + (p2[1] - p1[1])**2)
        #return n.abs(p2 - p1)
    
    s = 0
    for i in xrange(len(guy['dna'])-1):
        c1 = guy['dna'][i]
        c2 = guy['dna'][i+1]
        s = s + dist(c1, c2)
    
    return s

def new_guy():
    dna = range(1, CITIES)
    r.shuffle(dna)
    
    d = {'dna': [0] + dna,
         'fitness': .0,
         'score': .0,
         'parents': []}
    return d.copy()

# PCA #########################################################################

def pca(data):
    data = n.array(data, dtype=float)
    
    # normalizamos a matriz de dados (X = X - mean) e dividimos pelo d.p.
    #     X = (X - mean) / dp
    for i in xxrange(data.shape[1]):
        # adiciono um valor irrisorio 0.001 no denominador para nao
        # dar divisao por zero
        data[:,i] = (data[:,i] - data[:,i].mean())/(data[:,i].std()+0.001)
    
    # calculamos a matriz de covariância de X
    matriz_cov = n.cov(data, bias=1, rowvar=0)
    
    # calculamos os autovetores e autovalores e ordenamos em ordem decresc.
    autovalores, autovetores = n.linalg.eig(matriz_cov)
    args = n.argsort(autovalores)[::-1]
    autovalores = autovalores[args]
    autovetores = autovetores[args]
    
    # calculamos os componentes principais para todos os dados
    dados_finais = n.dot(autovetores.T, data.T)
    principais = dados_finais.T
    
    return principais

# Cluster #####################################################################

# Util ########################################################################

def ellipse(x, y, a, b, angle, steps):
    """ Returns the points for an ellipse centered in x, y with size a, b """
    beta = -angle * (n.pi / 180)
    sinbeta = n.sin(beta)
    cosbeta = n.cos(beta)

    alpha = n.linspace(0, 360, steps).T * (n.pi / 180)
    sinalpha = n.sin(alpha)
    cosalpha = n.cos(alpha)

    X = x + (a * cosalpha * cosbeta - b * sinalpha * sinbeta)
    Y = y + (a * cosalpha * sinbeta + b * sinalpha * cosbeta)

    ell = []
    for i in xrange(steps):
        ell.append([X[i], Y[i]])

    return n.array(ell)

# Main loop ###################################################################

num_elite = int(POP_SIZE * ELITE_RATE)
num_op = int(POP_SIZE * OP_RATE)
num_rand = int(POP_SIZE - num_elite - num_op)
print num_elite, num_op, num_rand
cities = ellipse(0, 0, 1, 1, 0, CITIES+1)

clusters = []

# inicializa a população com novos indivíduos aleatórios
pops = []
pop = []
for i in xrange(POP_SIZE):
    pop.append(new_guy())

for generation in xrange(GENERATIONS):
    # copia a elite ('num_elite' primeiros) para a nova população
    elite = []
    new_pop = []

    for i in xrange(num_elite):
        elite.append(copy.deepcopy(pop[i]))
        new_pop.append(copy.deepcopy(pop[i]))

    # aplica operadores de crossover e mutação apenas na elite, criando novos
    for i in xrange(num_op/2):
        # crossover: aplicado a dois elementos da elite, escolhidos ao acaso
        mom = r.choice(elite)
        dad = r.choice(elite)
        sis = None
        bro = None
        if r.random() < CROSSOVER_RATE:
            (sis, bro) = crossover(mom, dad)
        else:
            sis = copy.deepcopy(mom)
            bro = copy.deepcopy(dad)

        # mutation
        if r.random() < MUTATION_RATE:
            sis = mutate(sis)
            bro = mutate(bro)

        # store parents
        #sis['parents'] = [dad, mom]
        #bro['parents'] = [mom, dad]
        
        # store new guys in the new pop
        new_pop.append(sis)
        new_pop.append(bro)

    # o restante de new pop é obtido criando-se novos aleatórios
    for i in xrange(num_rand):
        ne = new_guy()
        new_pop.append(ne)

    # calcula o custo de cada indivíduo
    for i in xrange(POP_SIZE):
        sc = score(new_pop[i])
        new_pop[i]['score'] = sc

    # atualiza o fitness de cada indivíduo
    for i in xrange(POP_SIZE):
        fi = fitness(new_pop[i])
        new_pop[i]['fitness'] = fi

    # sobrescreve a população antiga pela nova
    pop = new_pop[:]

    # clusteriza
    # components = pca([guy['dna'] for guy in pop])
    # points = []
    # for i in xrange(len(components)):
    #     coords = components[i][:2]
    #     ref = pop[i]
    #     points.append(Point(coords, ref))
    # clusters = kmeans(points, 2, .05)

    # escolhe representante

    # *** TODO ***

    # ordenamos a população em função de seu fitness (maior -> menor)
    pop.sort(key=lambda x: x['fitness'], reverse=True)

    pops.append({'generation': generation,
                 'pop': pop,
                 'best': pop[0],
                 'best fitness': pop[0]['fitness'],
                 'fitness avg': sum([x['fitness'] for x in pop]) / len(pop),
                 'fitness min': min([x['fitness'] for x in pop]),
                 'fitness max': max([x['fitness'] for x in pop]),
                 'score avg': sum([x['score'] for x in pop]) / len(pop),
                 'score min': min([x['score'] for x in pop]),
                 'score max': max([x['score'] for x in pop])})


    print '*' * 10
    print 'generation:             ', generation
    print 'best                    ', pop[0]['dna']
    print 'best fitness:           ', pop[0]['fitness']
    print 'max fitness:            ', max([x['fitness'] for x in pop])

#
# PLOT #1: FITNESS
#
p.figure()
x = []
y = []
yerr_max = []
yerr_min = []
ymin = []
ymax = []
for po in pops:
    x.append(po['generation'])
    y.append(po['fitness avg'])
    ymin.append(po['fitness min'])
    ymax.append(po['fitness max'])
    yerr_max.append(ymax)
    yerr_min.append(ymin)

p.plot(x, ymin, 'g-')
p.plot(x, ymax, 'r-')
p.plot(x, y, '-')
p.xlabel('Generation')
p.ylabel('Fitness Min/Avg/Max')
p.grid(True)

p.figure()
x = []
y = []
yerr_max = []
yerr_min = []
ymin = []
ymax = []
for po in pops:
    x.append(po['generation'])
    y.append(po['score avg'])
    ymin.append(po['score min'])
    ymax.append(po['score max'])
    yerr_max.append(ymax)
    yerr_min.append(ymin)

p.plot(x, ymin, 'g-')
p.plot(x, ymax, 'r-')
p.plot(x, y, '-')
p.xlabel('Generation')
p.ylabel('Score Min/Avg/Max')
p.grid(True)


p.show()