Por: Jesús Mager

Se presenta la familia de los algoritmos evolutivos y se trabaja el caso de un algoritmo genético binario para la solución de un problema no lineal. Se trata de encontrar el volumen mínimo de una lata data su fórmula y una serie de restricciones. A través del experimento se localiza que el resultado es variable, y no cuenta un óptimo local, pero esto no implica un óptimo global. Se analizan diferentes métodos para crear nuevas generaciones y seleccionar individuos. 

Introducción

Los algoritmos genéricos han sido ampliamente usados como herramientas de búsqueda para problemas de optimización no lineales.  El concepto de de algoritmo genético fue usado por primera vez por John Holland. En esta práctica se va a usar un algoritmo genético con codificación binaria. 

Estos algoritmos son bioinspirados, basados en los principios de la genética en la naturaleza. Se puede codificar un atributo en cadenas binarias para usarla como un cromosoma.  Cada uno de los bits puede ser interpretado a su vez como un gen, donde a cambios en este, también habrá cambios en el fenotipo, o solución. Por consiguiente se debe de elegir una representación del las soluciones y generar una población de cadenas a azar. Al elegir una cadena de bits, esta cadena tiene una longitud y tendrá la capacidad de utilizar máximo y mínimo número de soluciones. Esta cantidad marcará el la cota superior e inferior de soluciones posibles [1]

La evaluación de una solución significa calcular las funciones objetivos y evaluar las violaciones a las restricciones. Por lo anterior se debe definir una métrica para combinar ambas cuestiones, llamada aptitud (fitness). Una condición de terminación también debe ser definida y revisada. A continuación la población de soluciones es modificada con operaradores genéticos que generará una nueva generación. La reproducción de la población tiene los siguientes objetivos: 

  • Encontrar buenas soluciones(que mejoran el promedio)
  • Hacer varias copias de una buena solución
  • Eliminar malas soluciones, de tal manera que puedan insertarse las copias de las soluciones buenas.

Para resolver estas tareas se tienen un conjunto de métodos. Entre los más comunes son la selección por torneo, la selección proporcional y la selección por ranking. En la selección por torneo se selecciona a dos soluciones y la mejor es tomada y agregada a los especímenes seleccionados. Se repite el ciclo hasta que los de espacios de selección se encuentran llenos[1].    

La asignación proporcional las soluciones son copias asignadas, según el número que corresponde a su aptitud usando 

\[\frac{f_i}{f_{media}}.\]

Se puede usar el concepto de una rueda de la fortuna, donde la rueda se divide en en N divisiones, donde el tamaño de cada una es marcada por la proporción antes descrita. Dado que para sacar la media de fitness hay que recorrer toda la lista, el costo computacional de este método es superior al del torneo. Existen varios métodos para aproximarse a la media. Este método también tiene un problema de escalabilidad.

Ahora un ejemplo

El problema que se quiere solucionar es minimizar el volumen de un  cilindro para lo cual se necesitarán únicamente dos parámetros: diámetro d y altura h. Se va a requerir que los cilindros tengan al menos 300ml. El planteamiento del problema de programación no lineal es el siguiente:

 

\[     \min d(d,h) = c (\frac{\pi d^2}{2}+\pi dh),\\     \text{s.a.} \\               g_1(d,h) = \frac{\pi d^2h}{4} \geq 300, \\           d_{\min} \leq d \leq d_{\max}, \\                h_{\min} \leq h \leq {h_\max}  \]

El parámetro c es el costo del material por centímetro cúbico y las variable de decisión son d y h están acotados por

 

\[ [d_{\min}, d_{\max}] y [d_{\min}, d_{\max}] \]

A continuación se expondrá el código para implementar el algoritmo con las siguientes características paso a paso:

Se pueda ajustar una nueva longitud del genoma, número de individuos por población, número de poblaciones

Defino una clase que va a adquirir una serie de atributos con los cuales se generaran los tipos de datos correspondientes. Se puede ver también un método que permite la transformación de valores enteros a binarios y así pode trabajar con el algoritmo.

class Lata:
    def __init__(self, n=10, maxsol=10, gensize=5, term=100, crossf=0.5, prob=0.3):
        self.n = n
        self.maxsol = maxsol
        self.pob = []
        self.gen = 0
        self.gensize = gensize
        self.term = term
        self.crossf = crossf
        self.prob = prob
        self.best = [100000000000000,0]

    # Genera la población inicial, según el tamaño de maz sol. 
    def gen_population(self):
        # Genera pares aleatorios
        dv = np.random.uniform(low=0, high=self.maxsol, size=self.n)
        hv = np.random.uniform(low=0, high=self.maxsol, size=self.n)
        # Posteriormente los converte en binarios y llena con ceros
        print("Población generada: ")
        for i in range(len(dv)):
            d = "{0:b}".format(int(dv[i]))
            h = "{0:b}".format(int(hv[i]))
            d = d.zfill(self.gensize)
            h = h.zfill(self.gensize)
            print(d,h)
            self.pob.append([d,h])

Ajustar el criterio de terminación, tal que si el fitness se alcanza termine el proceso evolutivo, o bien al termino del número de generaciones reporte el individuo hallado

Nuestra criterio de termino es son dos, puede ser por número de generaciones, que es posible ajustar desde el inicio de clase o elegir el primer elemento factible. Para determinar la factibilidad tenemos la función g y para caluclar el fitness la función f. Esas se enumeran a continuación. 

    def f(self, pair):
        d = 0.65*(float(int(pair[0],2)))
        h = float(int(pair[1],2))
        res = ((pi*pow(d,2))/2) + pi*d*h
        return res

    def g(self,pair):
        print(pair)
        d = float(int(pair[0],2))
        h = float(int(pair[1],2))
        res = ((pi*pow(d,2)*h)/4)-300
        return res

Cambiar el porcentaje de cruza y mutación, iniciar con 50% y 50% respectivamente además de Iniciar con cruza en un punto 

En el código se puede ver el valor selr.crossf, que representa el punto en el cual se inicia la cruza de dos pares de valores. Estos se modifican según este valor de punto flotante y se transforma en posiciones de los arreglos. Posteriormente se intercambian.

   def cross(self, pair1, pair2):
        a1 = pair1[0]
        b1 = pair1[1]

        a2 = pair2[0]
        b2 = pair2[1]

        por = self.crossf


        size = len(a1)
        fin = int(size * por)
        ini = size - fin

        na1 = a1[:ini] + b2[fin+1:]
        nb1 = b1[:ini] + a2[fin+1:]

        na2 = a2[:ini] + b1[fin+1:]
        nb2 = b2[:ini] + a1[fin+1:]
        return [[na1, nb1], [na2,nb2]]

Iniciar con 10% para porcentaje de mutación por individuo

En este método se va a realizar la mutación por probabilidad. Se recorrerá cada elemento de la cadena y se generará números aleatorios con un porcentaje de éxito definido por el usuarios. En este caso utilizamos el 10%

    # Función que genera la multación de los individuos
    def mutt(self, pair):
        new = []
        for word in pair:
            newgen = []
            for l in word:
                r = np.random.choice(2 , p=[1-self.prob, self.prob])
                if r:
                    if l == '1':
                        newgen.append('0')
                    else:
                        newgen.append('1')
                else:
                    newgen.append(l)
            new.append("".join(newgen))
        print(" -- ", str(new))
        return new

Dibujar el cilindro-lata, con los valores hallados de "h" y "d" respectivamente

Para poder visualizar mejor el cilindro generado por el algoritmo se va a graficar este mismo. La figura generada se muestra a continuación, al igual que el código necesario para ejecutarse. 

#Método que grafica los valores encontrados en forma de cilindros
    def drawCyl(self): 

        d  = int(self.best[1][0],2)
        h  = int(self.best[1][1],2)
        center = 0

        fig = plt.figure()
        ax = fig.add_subplot(111, projection='3d')
        x=np.linspace(center-d, center+d, 100)
        z=np.linspace(0, h, 100)
        Xc, Zc=np.meshgrid(x, z)
        Yc = np.sqrt(d**2 -(Xc - center) ** 2) + center 
        # Draw parameters
        rstride = 20
        cstride = 10
        ax.plot_surface(Xc, Yc, Zc, alpha=0.2, rstride=rstride, cstride=cstride)
        ax.plot_surface(Xc, (2*center-Yc), Zc, alpha=0.2, rstride=rstride, cstride=cstride)

        ax.set_xlabel("X")
        ax.set_ylabel("Y")
        ax.set_zlabel("Z")

        ax.set_title("Mejor cilindro encontrado por el algoritmo evolutivo")

        plt.show()

¿Y cómo mejora la solución de la lata a través de las iteraciones?

Conclusion

Los algoritmos evolutivos permiten resolver problemas de programación no lineal de manera eficiente, pero no garantizan encontrar la solución óptima. Si bien la optimización clásica si logra encontrar esta solución óptima global, el costo de cálculo es muy alto. Si bien, en ocasiones es necesario encontrar el óptimo global, en una gran parte de las aplicaciones bastaría con encontrar una solución aceptable, que puede ser traducida como un óptimo local. 

El control de la composición de la población de la siguiente generación, la forma de codificación de los resultados y en general la aleatoriedad introducida por la mutación son parámetros controlados por el investigador, que le permitirá acercarse a los resultados estimados de una manera más eficiente. 

El código completo

# License: GPL v3+
# Copyright: Jesús Mager 2017
import numpy as np from math import pi, pow import matplotlib.pyplot as plt from mpl_toolkits.mplot3d import Axes3D import mpl_toolkits.mplot3d.art3d as art3d from matplotlib.patches import Circle class Lata: def __init__(self, n=10, maxsol=32, gensize=5, term=10, prob=0.3): self.n = n self.maxsol = maxsol self.pob = [] self.gen = 0 self.gensize = gensize self.term = term self.prob = prob self.best = [100000000000000,0] np.random.seed() # Genera la población inicial, según el tamaño de maz sol. def gen_population(self): # Genera pares aleatorios dv = np.random.uniform(low=0, high=self.maxsol, size=self.n) hv = np.random.uniform(low=0, high=self.maxsol, size=self.n) # Posteriormente los converte en binarios y llena con ceros print("Población generada: ") for i in range(len(dv)): d = "{0:b}".format(int(dv[i])) h = "{0:b}".format(int(hv[i])) d = d.zfill(self.gensize) h = h.zfill(self.gensize) print(str(i),":",str(d),str(h)) self.pob.append([d,h,i]) #Ciclo principal def calc(self): self.data = [] while True: print(" * Probar fitess y condiciones") for pair in self.pob: res = self.g(pair) cond = self.f(pair) # ¿Si es el mejor elemento, guardar? if res > 0: if self.best[0] > cond: self.best = [cond, pair] # Si no se ha alcanzado la condición de término if not self.term: if res > 0: print("Solución encontrada 1: ", str(self.best),end=" ") print(res, cond) return pair else: if self.gen > self.term: print("Solución encontrada 2: ", str(self.best), end=" ") return self.best # Una revisada la condición de término se pasa a # crear la nueva generación print(" * Crear nueva generación") new_gen = [] i = 0 count = 0 #Contador de la nueva generación # Genera nueva población por cruza while i < ((len(self.pob))/4): # Selecciona aleatoria oponentes r = np.random.uniform(0, len(self.pob), size=2) j = int(r[0]) k = int(r[1]) #Primer torneo! if self.f(self.pob[j]) < self.f(self.pob[k]): winner1 = j else: winner1 = k # Selecciona aleatoriamente los contrincantes r = np.random.uniform(0, len(self.pob), size=2) j = int(r[0]) k = int(r[1]) #Segundo torneo if self.f(self.pob[j]) < self.f(self.pob[k]): winner2 = j else: winner2 = k # Ahora hacemos la cruza de los mejores individuos elem = self.cross(self.pob[winner1], self.pob[winner2]) elem[0].append(count) count += 1 elem[1].append(count) count += 1 print(elem[0]) print(elem[1]) # Agrego a la nueva generación new_gen.append(elem[0]) new_gen.append(elem[1]) i += 1 i = 0 # Ahora vamos a mutar! while count < len(self.pob)-1: # Seleccionamos la primer víctima r = np.random.randint(0,len(self.pob)) e = self.mutt(self.pob[r][:2]) e.append(count) print(e) new_gen.append(e) count += 1 # Seleccionamos la primer víctima r = np.random.randint(0,len(self.pob)) e = self.mutt(self.pob[r][:2]) e.append(count) new_gen.append(e) print(e) count += 1 i += 1 print(new_gen) self.pob = new_gen # Transformamos a binarios :) if self.best[1] == 0: x1 = np.random.uniform(low=0, high=self.maxsol) x2 = np.random.uniform(low=0, high=self.maxsol) d = "{0:b}".format(int(x1)) d = d.zfill(self.gensize) h = "{0:b}".format(int(x2)) h = h.zfill(self.gensize) self.best[0] = self.f([d,h]) self.best[1] = [d,h] # Agregamos el mejor de la generación a la próxima best = self.best[1] best.append(count) self.pob.append(best) self.gen += 1 self.data.append([self.gen, self.best]) print("Mejor solución de la generación: ", str(self.best)) print("Volumen:", str(300+self.g(self.best[1]))) print("d:", str(int(self.best[1][0],2))) print("h:", str(int(self.best[1][1],2))) input(); print("-------------------- Nueva Generación -------------------") # Cálculo del la función de fitness def f(self, pair): d = 0.65*(float(int(pair[0],2))) h = float(int(pair[1],2)) res = ((pi*pow(d,2))/2) + pi*d*h return res # Cálculo de la restrucción def g(self,pair): print(pair) d = float(int(pair[0],2)) h = float(int(pair[1],2)) res = ((pi*pow(d,2)*h)/4)-300 return res # Se realiza la cruza tomando el punto de inicio self.crossf def cross(self, pair1, pair2): a1 = pair1[0] b1 = pair1[1] a2 = pair2[0] b2 = pair2[1] # Generamos nuevos números fin = np.random.randint(2,7) # Punto de cruza na1 = a1[:fin] + b2[fin:] nb1 = b1[:fin] + a2[fin:] fin = np.random.randint(2,7) # Punto de cruza 2 na2 = a2[:fin] + b1[fin:] nb2 = b2[:fin] + a1[fin:] return [[na1, nb1], [na2,nb2]] # Función que genera la multación de los individuos def mutt(self, pair): new = [] print(pair[0]) p1 = np.array(list(pair[0])) p2 = np.array(list(pair[1])) np.random.shuffle(p1) np.random.shuffle(p2) print(p1) pn1 = "".join(p1) pn2 = "".join(p2) new.append(pn1) new.append(pn2) return new ) #Método que grafica los valores encontrados en forma de cilindros def drawCyl(self): d = int(self.best[1][0],2) h = int(self.best[1][1],2) center = 0 fig = plt.figure() ax = Axes3D(fig, azim=30, elev=30) x=np.linspace(center-d, center+d, 100) z=np.linspace(0, h, 100) Xc, Zc=np.meshgrid(x, z) Yc = np.sqrt(d**2 -(Xc - center) ** 2) + center # Draw parameters rstride = 20 cstride = 10 ax.plot_surface(Xc, Yc, Zc, alpha=0.2, rstride=rstride, cstride=cstride) ax.plot_surface(Xc, (2*center-Yc), Zc, alpha=0.2, rstride=rstride, cstride=cstride) ax.set_xlabel("X") ax.set_ylabel("Y") ax.set_zlabel("Z") ax.set_title("Mejor cilindro encontrado por el algoritmo evolutivo") plt.show() def draw(self): x = [i[0] for i in self.data] y = [int(i[1][0]) for i in self.data] print(x) print(y) plt.plot(x,y) plt.title("Mejor elemento global: desarrollo de volúmen del cilindo por cada generació por cada generación") plt.xlabel("Generación") plt.ylabel("Volúmen de la lata") plt.show() if __name__ == "__main__": lata = Lata() lata.gen_population() lata.calc() lata.drawCyl() lata.draw()

 Referencias

[1] K. Deb and D. Kalyanmoy, Multi-Objective Optimization Using Evolutionary Algorithms. New York, NY, USA: John Wiley & Sons, Inc., 2001.

Share This