Generalización de algoritmos de búsqueda estocástica local por medio de métodos de representación conjunta de soluciones

Trabajo de Fin de Grado Curso 2023–2024 Generalización de algoritmos de búsqueda estocástica local por medio de métodos de representación conjunta de soluciones Generalization of stochastic local search algorithms via methods of common representation of solutions Autora Sara Vicente Arroyo Directores Ismael Rodríguez Laguna Fernando Rubio Diez Doble Grado en Ingeniería Informática y Matemáticas Departamento de Sistemas Informáticos y Computación Facultad de Informática Universidad Complutense de Madrid Madrid, 13 de septiembre de 2024 A mis padres, por aguantar mis llantos durante diez periodos de exámenes consecutivos, y a mis compañeros, por hacerlos mucho más llevaderos iii Agradecimientos Al Dr. Damar Wicaksono y al Prof. Dr. Michael Hecht, por su tiempo y sus aclaraciones sobre el fundamento teórico y los detalles de implementación de la librería Minterpy. Las conversaciones mantenidas con ellos mediante correo electrónico y videoconferencia han resultado fundamentales a la hora de encauzar este Trabajo de Fin de Grado. A mis directores, Ismael y Fernando, por todo lo que he aprendido de ellos. v Resumen Generalización de algoritmos de búsqueda estocástica local por medio de métodos de representación conjunta de soluciones Durante las últimas décadas, se ha producido un aluvión de contribuciones al campo de la computación evolutiva. Se han desarrollado algoritmos inspirados en todo tipo de especies naturales, además de híbridos entre varios de ellos. No obstante, el comportamiento intrínseco de muchos es esencialmente el mismo. Esto ha permitido dividirlos en tres grandes grupos: los que se parecen a los algoritmos genéticos (GA), a los de enjambre de partículas (PSO) y a los de colonia de hormigas (ACO), respectivamente. En este trabajo se da un paso más y se diseña un algoritmo generalizado al que pueden reducirse todos. Además, se propone una particularización del mismo (AEMP) que centra su atención únicamente en la distribución de probabilidad de las soluciones consideradas en cada momento. Este nuevo algoritmo permite crear variantes e híbridos entre ellas de forma completamente directa, sin necesidad de tomar inspiración en la naturaleza. Los experimentos realizados en este trabajo muestran, además, que no existen diferencias significativas entre los resultados obtenidos por ciertas variantes de AEMP y métodos clásicos como PSO al enfrentarse a problemas de optimización NP-difíciles. Palabras clave Computación evolutiva, algoritmos genéticos, algoritmos de enjambre de partículas, algoritmos de colonia de hormigas, optimización combinatoria, generalización, distribuciones de probabilidad, interpolación numérica, problemas NP-difíciles, pruebas no paramétricas. vii Abstract Generalization of stochastic local search algorithms via methods of common representation of solutions Over the last decades, there has been a flood of contributions to the field of evolutionary computation. Algorithms inspired in all kinds of natural species have been developed, as well as hybrids between several of them. However, the intrinsic behavior of many is essentially the same. This has allowed authors to sort them into three large groups: those similar to genetic algorithms (GA), to particle swarm optimization (PSO) algorithms and to ant colony optimization (ACO) algorithms, respectively. In this work, we go one step forward and design an algorithm that generalizes all of them. Moreover, we propose a particularization of it (AEMP) that only focuses on the probability distribution of the solutions considered at any moment. This novel algorithm allows us to create variants and hybrids between them in a way that is completely direct, without the need for natural inspiration. Experiments conducted in this work also show that no significant difference exists between results obtained by certain AEMP variants and classical methods such as PSO when facing NP-hard optimization problems. Keywords Evolutionary computation, genetic algorithms, particle swarm optimization algorithms, ant colony optimization algorithms, combinatorial optimization, generalization, probability distributions, numerical interpolation, NP-hard problems, nonparametric tests. ix Capítulo 1 Introducción Durante las últimas décadas, cientos de algoritmos evolutivos han invadido la literatura, inspirados en todo tipo de especies, poblaciones o fenómenos naturales (o no tan naturales). Gatos [6], delfines [49, 24], galaxias [37], bacterias [36, 7, 32], abejas [2, 22, 54], tumores [53], futbolistas [43, 44, 1] o incluso zombis [39] son tan solo algunas de las fuentes de inspiración o metáforas seleccionadas para alargar la lista de integrantes de este paradigma. Además, otros autores han desarrollado algoritmos híbridos, que consisten simplemente en la aplicación consecutiva de los operadores de unos y otros [23, 56]. En los últimos años han surgido voces críticas con la mencionada explosión de contribuciones al campo de la computación evolutiva. En 2015, Sörensen [52] insta a analizar la dinámica subyacente de todos estos algoritmos, sin distraerse con las particularidades de métodos supuestamente innovadores. Por su parte, Molina et al. [33] proponen en 2020 una taxonomía que clasifica más de trescientos de ellos en únicamente tres grandes grupos, según su comportamiento real (movimiento de vector diferencial,combinación de soluciones yestigmergia). Argumentan, además, que dicho comportamiento es más relevante que el fenómeno biológico, social o natural en el que se inspiran. El objetivo de este trabajo es dar un paso más y definir un algoritmo evolutivo generalizado al que puedan reducirse todos ellos. De esta forma, crear ligeras variantes e incluso híbridos entre distintas particularizaciones será trivial y no requerirá del desarrollo de todo un escenario metafórico nuevo para cada pequeño cambio en su mecánica subyacente. Asimismo, se pretende diseñar e implementar nuevas particularizaciones de dicho algoritmo capaces de enfrentarse a problemas de optimización en los que resulte completamente prohibitivo el uso de una simple búsqueda exhaustiva. Más concretamente, problemas NP-difíciles de dimensión superior a 30. En base a los objetivos antedichos, se ha seguido el plan de trabajo expuesto a continuación, que consta de una parte teórica y otra práctica. En primer lugar, se ha estudiado en profundidad la literatura relativa a tres algoritmos bioinspirados clásicos: los genéticos [21], los de enjambre de partículas [26] y los de colonia de hormigas [10]. Se trata de los representantes indiscutibles de cada uno de los grupos en los que Molina et al. [33] clasifican los algoritmos 1 2Capítulo 1. Introducción evolutivos. En el Capítulo 3 se introduce brevemente al lector a su funcionamiento. A continuación, se han explorado trabajos previos de generalización o, al menos, de alejamiento de metáforas y otras particularidades. Los algoritmos de estimación de la distribución (EDA, por sus siglas en inglés) [28] construyen en cada iteración un modelo probabilístico que, como sugiere su nombre, estima la distribución de las mejores soluciones del problema de optimización correspondiente y permite obtener nuevas muestras según dicha distribución. Por otro lado, el algoritmo CMA-ES (estrategia de evolución por adaptación de la matriz de covarianza) [19] también utiliza el muestreo aleatorio para dirigir la búsqueda hacia las mejores soluciones. Asume que estas siguen una distribución normal, cuyos parámetros (media y matriz de covarianza) se van actualizando iterativamente. En el Capítulo 4 se aporta una introducción más detallada a los EDA y al algoritmo CMA-ES. Al comienzo del Capítulo 5 se propone, por fin, un algoritmo propio capaz de generalizar todos los algoritmos evolutivos. Se basa en la existencia de una variable aleatoria subyacente que privilegia unas soluciones frente a otras (es decir, les asigna una mayor o menor probabilidad) y que se actualiza tras cada iteración. La parte práctica de este trabajo ha consistido, en primer lugar, en la investigación de distintas técnicas de interpolación y muestreo para su uso en nuevas particularizaciones del algoritmo generalizado, tal y como se justifica en la Sección 5.2 del Capítulo 5. Para ello, se ha contactado con expertos del ámbito de la optimización numérica del Helmholtz-Zentrum Dresden-Rossendorf (HZDR) en Dresde, Alemania, cuyos consejos han resultado fundamentales para encauzar las ideas de diseño hacia algoritmos verdaderamente implementables. Finalmente, se ha desarrollado una particularización del algoritmo generalizado, el algoritmo evolutivo por medias ponderadas (AEMP), que pone el foco en la distribución subyacente de cada conjunto de soluciones consideradas. Se han diseñado e implementado cuatro variantes del mismo de forma completamente directa, sin necesidad de definir operadores o escenarios metafóricos. Los resultados obtenidos tras enfrentar dichas variantes, así como los algoritmos clásicos, a problemas de optimización NP-difíciles se presentan y se discuten en el Capítulo 6. Capítulo 2 Preliminares En este capítulo se introducen definiciones que serán de utilidad en los capítulos posteriores, así como la notación utilizada. 2.1. Teoría de la probabilidad Definición 2.1.1. Una variable aleatoria sobre un espacio de probabilidad Ωes una función X: Ω →Rnmedible. 2.2. Optimización Definición 2.2.1. Dados un conjunto Θ⊂Rny una función f:Rn→R, definimos el problema de minimización como minimizar f(x) sujeto a x∈Θ Decimos que Θes la región factible, el espacio de búsqueda o el espacio de soluciones del problema, y que fes la función objetivo. Definición 2.2.2. Dados un conjunto Θ⊂Rny una función f:Rn→R, decimos que x∈Rnes solución factible del problema de minimización si x∈Θ. Definición 2.2.3. Dados un conjunto Θ⊂Rny una función f:Rn→R, decimos que x∗∈Rnes solución óptima del problema de minimización si es solución factible y f(x∗)≤f(x)∀x∈Θ Observación 2.2.4. El problema de maximización se define de forma análoga. Definición 2.2.5. Un problema de optimización combinatoria es un problema de optimización (es decir, de minimización o de maximización) cuya región factible Θes un conjunto finito. 3 Capítulo 3 Algoritmos evolutivos clásicos En este capítulo se resumen las ideas fundamentales sobre tres de los algoritmos evolutivos inspirados en la naturaleza más importantes de la literatura, como son los algoritmos genéticos (véase la Sección 3.1), los de enjambre de partículas (véase la Sección 3.2) y los de colonia de hormigas (véase la Sección 3.3). Se han elegido estos tres por ser, de alguna manera, los representantes canónicos de cada una de las categorías en las que Molina et al. [33] dividen los algoritmos evolutivos: combinación de soluciones,movimiento de vector diferencial yestigmergia, respectivamente. 3.1. Algoritmos genéticos Los algoritmos genéticos (genetic algorithms, GA) [21] toman su inspiración en la teoría clásica de la evolución y la selección natural [8]. Se emplean principalmente para resolver de forma aproximada un problema de optimización en Θ⊂Rn optimizar f(x) sujeto a x∈Θ donde el espacio de soluciones o población Θpuede ser continuo o discreto. Cada solución factible x∈Θes un individuo y la función objetivo f:Rn→Ro función de fitness describe cuán adaptado al medio se encuentra cada uno de ellos. Si el problema es de maximización, los individuos xcon un alto valor de fitness f(x) serán las mejores soluciones. Por el contrario, si el problema es de minimización, buscaremos individuos con niveles bajos de fitness. Observamos que x= (x1, . . . , xn) es un vector de ncomponentes. A menudo este se interpreta como un cromosoma donde cada componente es un gen. El esquema general de un algoritmo genético consta de tres fases que se repiten en cada iteración. Dada una generación de individuos, es decir, un subconjunto de Θ, se seleccionan aquellos mejor adaptados para dar lugar a la siguiente. Los individuos seleccionados se cruzan dos a dos y así se origina su descendencia, que habrá heredado características de ambos. Finalmente, al igual que en la naturaleza, 5 6Capítulo 3. Algoritmos evolutivos clásicos pueden producirse mutaciones, es decir, pequeños cambios en el genoma de los nuevos individuos. Los operadores de selección, cruce y mutación se definirán de forma distinta según el problema particular al que nos enfrentemos, y su aplicación se repetirá de forma iterativa hasta que se cumpla una cierta condición de parada. También será preciso definir la manera en la que se elegirá la generación inicial de individuos, habitualmente por muestreo aleatorio. 3.1.1. Operadores A lo largo de los años se han desarrollado multitud de operadores de selección, cruce y mutación distintos. En esta sección se enumeran los más relevantes, aunque no se profundiza en ellos. Sus respectivas definiciones se pueden consultar en gran cantidad de fuentes [p. ej. 13, Capítulos 4-5]. Selección Un conjunto de individuos es seleccionado para transmitir sus genes a la siguiente generación, y los demás son descartados. Por norma general, este operador se define de forma que asigne una menor, pero no nula, probabilidad de ser escogidos a aquellos con peor fitness. Algunos métodos vastamente utilizados son la denominada selección por torneo, la selección por ruleta y la basada en ranking. Cruce Los cromosomas seleccionados se consideran ahora progenitores, y se aparean dos a dos. Cada pareja da lugar a dos nuevos cromosomas hijos, definidos como una combinación de ambos. Los operadores de cruce más conocidos son los llamados cruce de un punto y cruce de dos puntos. No obstante, también cabe destacar el cruce uniforme [51], que da lugar a una mayor variabilidad genética. Mutación Por último, cada cromosoma es susceptible de sufrir un pequeño cambio según una probabilidad de mutación elegida. De esta forma, nos aseguramos poder explorar regiones del espacio de soluciones inaccesibles mediante aplicaciones de los otros dos operadores a la generación inicial de individuos. Mutaciones típicas son el cambio en el valor o en la posición de un gen, o la permutación de dos. 3.1.2. Variantes La definición de algoritmo genético y la de sus componentes se pueden adaptar a otros problemas de optimización, en los que el espacio de búsqueda no está compuesto por vectores numéricos sino por recorridos sobre grafos o permutaciones de elementos [16]. También existen variantes que eliminan la intuición de reproducción sexual. En ellas los nuevos individuos no se generan a partir de dos progenitores, 3.2. Algoritmos de enjambre de partículas 7 sino potencialmente a partir de todo el material genético presente en la generación previa [35]. 3.2. Algoritmos de enjambre de partículas Los algoritmos de enjambre de partículas (particle swarm optimization, PSO) [26] tienen un punto en común fundamental con los algoritmos genéticos: en cada iteración definen explícitamente un conjunto de soluciones factibles del problema optimizar f(x) sujeto a x∈Θ Se alejan, sin embargo, de ellos en la forma de obtener un nuevo conjunto de soluciones a partir del anterior. Cada solución se visualiza en este caso como una partícula en movimiento a través del espacio de búsqueda, cuya posición se modifica según su propia velocidad einercia. Además, las partículas forman parte de un enjambre y reciben información sobre las mejores soluciones halladas por las demás, de forma que es el enjambre en su conjunto el que resuelve el problema de optimización. Estos algoritmos se inspiran en la naturaleza cooperativa observada en poblaciones de determinadas especies del reino animal, como enjambres de abejas, bancos de peces o bandadas de aves. La notación utilizada a lo largo de esta sección es una adaptación de la empleada por Dorigo et al. [11]. Definición 3.2.1. Una partícula es un par (x, v)con x, v ∈Rn. Llamamos a x posición y a vvelocidad de la partícula p= (x, v). Definición 3.2.2. Llamamos enjambre a una tupla ordenada de partículas (p1= (x1, v1), . . . , pm= (xm, vm)) y denominamos topología del enjambre a cualquier grafo G= (V={1, . . . , m}, E). Definición 3.2.3. Dado un enjambre de partículas con topología G, llamamos entorno ovecindario de la partícula pi,N(i), al conjunto de vértices adyacentes de ien G. Observación 3.2.4. Las definiciones de topología y de vecindario dependen exclusivamente de la posición que ocupa cada partícula en la tupla. Intuitivamente, esto nos permitirá seguir la pista de cada partícula individual a lo largo de su trayectoria, aunque sus coordenadas cambien. Nos encontramos en posición de formalizar el algoritmo básico de enjambre de partículas [25]. En la iteración tconsideramos un enjambre Pt= (pt 1= (xt 1, vt 1), . . . , pt m= (xt m, vt m)) = (pt i= (xt i, vt i))m i=1 8Capítulo 3. Algoritmos evolutivos clásicos donde cada posición xies una solución factible del problema, y una cierta topología Gsobre Pt. Consideramos también la tupla Mt= (bt 1, . . . , bt m) de las mejores posiciones para cada partícula hasta el momento, es decir, Mtverifica que para todo i∈ {1, . . . , m}existe t∗∈ {0, . . . , t}tal que bt i=xt∗ iyf(xt∗ i)≤f(xs i)∀s∈ {0, . . . , t} si el problema es de minimización. Finalmente, consideramos una tercera tupla Vt= (lt 1, . . . , lt m) de las mejores posiciones de los vecinos de cada partícula. Esto es, para un problema de minimización, Vtverifica que para todo i∈ {1, . . . , m}existe j∈ N (i)tal que lt i=bt jyf(lt i)≤f(bt k)∀k∈ N(i) En la iteración t+ 1 actualizamos el enjambre de la siguiente manera Pt+1 = (pt+1 i= (xt+1 i, vt+1 i))m i=1 donde vt+1 i=wvt i+α1Ut 1(bt i−xt i) + α2Ut 2(lt i−xt i) xt+1 i=xt i+vt+1 i y la topología Gqueda fija. El parámetro w∈Rse denomina inercia y los parámetros α1yα2,coeficientes de aceleración. Las matrices Ut 1yUt 2son diagonales y sus elementos de la diagonal se generan aleatoriamente en cada iteración según una distribución uniforme. Podemos observar que el desplazamiento de una partícula depende en esencia de tres factores. El término wvt iestá relacionado con la inercia o el momento lineal, la tendencia de la partícula a continuar desplazándose en la misma dirección y sentido. El término α1Ut 1(bt i−xt i)es la componente cognitiva omemorística, que lleva a la partícula de vuelta hacia las mejores posiciones que ha visitado. Por último, la componente α2Ut 2(lt i−xt i)es la social, que guía a la partícula hacia las mejores soluciones halladas por su entorno. 3.2.1. Variantes El primer problema que salta a la vista con el algoritmo clásico es la falta de garantía de permanecer en la región factible al avanzar una iteración. Por ello, existen variantes en las que se adapta la velocidad de las partículas en pos de mantener la factibilidad de las soluciones [14]. 3.3. Algoritmos de colonia de hormigas 9 3.3. Algoritmos de colonia de hormigas Los algoritmos de colonia de hormigas (ant colony optimization, ACO) [10] se inspiran en el comportamiento colaborativo, y a la vez distribuido, de estos insectos en su búsqueda de alimento. Cada una de ellas deja un rastro de sustancias químicas denominadas feromonas por el camino seguido hasta la fuente de sustento y de vuelta al nido, que posteriormente servirá de guía para otras hormigas de su colonia. Esta forma de cooperación y comunicación indirecta a través de modificaciones en el medio físico se conoce como estigmergia. En este paradigma, a diferencia de los dos anteriores, cada agente u hormiga artificial no representa una solución de forma explícita, sino que la construye iterativamente al recorrer una cierta estructura de datos. En este caso nos reducimos a problemas de optimización combinatoria, es decir, optimizar f(x) sujeto a x∈Θ donde Θes un conjunto finito. Más concretamente, si x= (x1, . . . , xn)∈Θ, cada variable de decisión xipodrá tomar valores en un conjunto finito Vi={vi 1, vi 2, . . . , vi Ni}. Observación 3.3.1. El recíproco no es cierto. Es decir, un vector x∈ V1× · · ·× Vn no tiene por qué ser factible, puede haber restricciones adicionales. Un ejemplo son los problemas de permutaciones, donde V1=· · · =Vn={v1, v2, . . . , vn}y no se puede dar xi=xjpara i=j. Las siguientes definiciones pueden resultar poco intuitivas, sobre todo si se conoce previamente el concepto de algoritmo de colonia de hormigas como método de optimización en grafos. Se incluyen porque son las que se utilizan típicamente en la literatura para introducir estos algoritmos [12], para dotarlos de aplicabilidad en la resolución de cualquier problema de optimización combinatoria. Al final de esta sección, no obstante, se da una interpretación gráfica de estas definiciones y se propone una alternativa más natural para problemas de grafos, que fueron la verdadera motivación de estos algoritmos. Definición 3.3.2. Dado un problema de optimización combinatoria, denominamos componente a una instancia de una variable de decisión, xi=vi j. Definición 3.3.3. Un grafo de construcción GC= (C, E)es un grafo completo tal que Ces el conjunto de todas las componentes del problema. A los elementos de Elos llamamos conexiones. Definición 3.3.4. Un estado es una secuencia finita de componentes. Además, decimos que es un estado factible si se puede completar a una solución factible. Definición 3.3.5. Si aes una conexión y s= (c0, c1, . . . , ck)es un estado, decimos que aforma parte de s, y lo denotamos por a⊑s, si a= (ci, ci+1)para algún i∈ {0, . . . , k −1}. Definición 3.3.6. Un rastro de feromonas es una función τ:E→R+. 16 Capítulo 4. Algunos intentos de generalización x1x2x3x4x5x6f(x) 1 1 1 1 1 1 1 6 2 1 0 1 0 1 1 4 3 1 1 1 1 1 0 5 4 0 1 0 1 1 1 4 5 1 1 1 1 0 1 5 6 1 0 0 1 1 1 4 7 0 1 0 1 1 0 3 8 1 1 1 0 1 0 4 9 1 1 1 0 0 1 4 10 1 0 0 1 1 1 4 11 1 1 0 0 1 1 4 12 1 0 1 1 1 0 4 13 0 1 1 0 1 1 4 14 0 1 1 1 1 0 4 15 0 1 1 1 1 1 5 16 0 1 1 0 1 1 4 17 1 1 1 1 1 0 5 18 0 1 0 0 1 0 2 19 0 0 1 1 0 1 3 20 1 1 0 1 1 1 5 Tabla 4.3: Ejemplo de conjunto A1tomado de [28, Tabla 3.3]. desarrollado algoritmos de estimación de la distribución más complejos, capaces de capturar relaciones de dependencia entre pares de variables aleatorias. Destacan en el ámbito de la optimización discreta los EDA que emplean redes bayesianas, esto es, grafos dirigidos acíclicos, para construir un modelo probabilístico que estime p(x|B)y obtener nuevas muestras a partir de él. Un representante importante de este paradigma es el BOA (Bayesian Optimization Algorithm) [42]. Para optimización continua, por su parte, destacan algoritmos como el EGNA (Estimation of Gaussian Networks Algorithm) [29], que utilizan redes gaussianas para capturar y modelar la distribución de las soluciones seleccionadas. Una fuente útil para profundizar en el funcionamiento de las redes bayesianas y gaussianas es [27]. 4.2. Adaptación de la matriz de covarianza El algoritmo CMA-ES (Covariance Matrix Adaptation - Evolution Strategy), presentado por Hansen y Ostermeier [19] en 2001, se encuadra dentro del conjunto de métodos denominados estrategias de evolución. Las estrategias de evolución surgen en la década de 1960 con el trabajo de Rechenberg [45], con el objetivo de ser utilizadas en problemas de optimización continua. Se caracterizan por la aplicación iterativa de operadores de selección, cruce y mutación, como en el caso de los algoritmos genéticos. La manera específica y conceptualmente diferente de definir y aplicar dichos operadores, sin embargo, ha 4.2. Adaptación de la matriz de covarianza 17 provocado que tradicionalmente se los considere clases distintas de algoritmos evolutivos. Una gran diferencia entre ellos es la importancia que le otorgan a la mutación. En un algoritmo genético clásico es habitual que este operador provoque pequeños cambios en algunos individuos, con el objetivo de alcanzar regiones inexploradas del espacio de soluciones. En una estrategia de evolución, por el contrario, la mutación dirige completamente el camino de los individuos hacia el óptimo. De hecho, en las versiones más antiguas, ni siquiera existe operador de cruce (algo impensable en un algoritmo genético), y los hijos se obtienen únicamente aplicando mutaciones a los progenitores. El operador de mutación canónico de esta clase de algoritmos consiste en sumar a cada individuo (x1, . . . , xn)∈Rnun vector obtenido como muestra de una distribución normal multivariante con media cero. Veamos un primer ejemplo sencillo, en el que además definimos el operador de cruce de la manera más habitual: como la media aritmética de los progenitores. N(µ, σ2)denota la distribución normal univariable de media µy desviación típica σ. Ejemplo 4.2.1. Seleccionamos los progenitores p1= (x1 1, x1 2, x1 3),p2= (x2 1, x2 2, x2 3), p3= (x3 1, x3 2, x3 3)∈R3. Tras aplicar el operador de cruce, obtenemos el hijo h=1 3 3 X j=1 pj, donde la suma vectorial se realiza componente a componente. A continuación, obtenemos una muestra {y= (y1, y2, y3)} ∼ X= (X1, X2, X3), con Xi∼ N(0, σ2) para todo i∈ {1,2,3}, dado un cierto σ > 0. Finalmente, añadimos a la nueva generación de individuos el hijo mutado, esto es, h+y. El Ejemplo 4.2.1 nos invita a formularnos algunas preguntas. ¿Cómo se elige el parámetro σ? ¿Se mantiene constante? ¿Tiene sentido utilizar el mismo para todas las componentes? Una característica importante de la mayoría de estrategias de evolución es, precisamente, la autoadaptación dinámica de los parámetros que rigen la mutación, para facilitar la convergencia. Podemos sustituir la expresión h+y por h+sy, donde s > 0es el paso, que aumentará cuando el algoritmo se esté dirigiendo de manera consistente hacia mejores soluciones, y viceversa. Por otro lado, efectivamente, la naturaleza de numerosos problemas de optimización requiere tratar de forma distinta cada dimensión. Algunas estrategias de evolución definen Xi∼ N(0, σ2 i), con una desviación típica σi>0distinta para cada i∈ {1, . . . , n}. No obstante, podemos dar un paso más y considerar también las dependencias entre variables. Para ello, obtenemos muestras de X∼ N(0,C), que ahora denota la distribución normal multivariante de media el origen y matriz de covarianza C, que captura dichas dependencias. Recordemos que la posición (C)i,j de la matriz representa la covarianza de las variables XiyXj, que es siempre nula cuando estas son independientes, pero no lo será en general. En el Algoritmo 2 se muestra el esquema general de la estrategia CMA-ES. El lector interesado en conocer todos los detalles matemáticos puede consultar [18]. 18 Capítulo 4. Algunos intentos de generalización Algoritmo 2 CMA-ES. Es necesario definir la condición de parada. Obtener una muestra A={y1, . . . , ym}∼N(0,I). mientras no se cumpla la condición de parada hacer µ←1 mPm j=1 yj Calcular el nuevo paso s. Calcular la nueva matriz de covarianza C. Obtener una muestra A={y1, . . . , ym} ∼ s· N (µ, C)∼ N(µ, s2C). fin mientras Capítulo 5 Algoritmo generalizado La idea central de este trabajo surje de observar que el comportamiento de todos los algoritmos evolutivos se puede describir de la manera siguiente: Al comienzo de cada iteración, está definido un cierto conjunto A⊂Θ(finito o infinito) de soluciones al problema de optimización considerado. En el caso de GA o PSO, A es el propio conjunto de individuos seleccionados o de partículas generadas. En ACO, por el contrario, Aes el conjunto de todos los posibles recorridos sobre el grafo de construcción con probabilidad positiva de ser escogidos por las hormigas. A continuación, se actualiza el conjunto Ay comienza la siguiente iteración. Este esquema nos lleva a imaginar una variable aleatoria subyacente en cada iteración, que privilegia unas soluciones frente a otras (es decir, les asigna una mayor probabilidad de formar parte del nuevo conjunto A). El Algoritmo 3 es una abstracción de todos los algoritmos evolutivos en estos términos. Xdenota el conjunto de todas las variables aleatorias sobre Θ, y Pes el conjunto potencia o partes de. Algoritmo 3 Algoritmo evolutivo generalizado. Es necesario definir la región factible Θ, la variable aleatoria inicial X0: Θ →R, la función de actualización act:X × P (Θ) → X y la condición de parada. X←X0 Obtener una muestra A⊂Θsegún la distribución inducida por X. mientras no se cumpla la condición de parada hacer X←act (X, A) Obtener una muestra A⊂Θsegún la distribución inducida por X. fin mientras Hasta ahora no ha aparecido la función fa optimizar, que a menudo llamaremos función de fitness por influencia de GA. Esto se debe a que, en principio, el Algoritmo 3no tiene por qué conducirnos a las mejores soluciones. Esto es importante. En su expresión más general, este solo busca capturar las actualizaciones sucesivas de su variable aleatoria intrínseca a partir de sí misma, que pueden ir, a priori, en cualquier sentido. Unas particularizaciones del algoritmo funcionarán mejor que otras, e intuitivamente estas serán precisamente las que utilicen la función fen la definición de act. 19 20 Capítulo 5. Algoritmo generalizado 5.1. Los algoritmos clásicos como particularizaciones del generalizado En cada iteración de GA, act (X, A)devuelve una nueva variable aleatoria que asigna probabilidad positiva e igual a todos los elementos del conjunto Bde los hijos de A, generado tras aplicar en este último los correspondientes operadores de selección, cruce y mutación, y nula al resto de elementos de Θ. De manera completamente análoga, en cada iteración de PSO la función act (X, A)devuelve la variable aleatoria que asigna probabilidad positiva e igual a todos los elementos del conjunto de nuevas partículas alcanzadas a partir de las de Asegún las reglas de actualización propias del algoritmo, y cero al resto. El caso de ACO es ligeramente distinto. Recordemos que en una cierta iteración de dicho algoritmo, Arepresenta el conjunto de todos los posibles recorridos sobre el grafo de construcción con probabilidad positiva, dada por la distribución de la variable X, de ser escogidos por las hormigas. Por tanto, act (X, A)devuelve la variable aleatoria que asigna probabilidad positiva (y posiblemente distinta) a cada recorrido sobre el grafo de construcción tras haber actualizado consecuentemente el rastro de feromonas. 5.2. Definición de otras particularizaciones A la hora de crear distintas particularizaciones del Algoritmo 3, parece razonable definir las correspondientes funciones de actualización de forma que nuevas variables aleatorias asignen probabilidades mayores a elementos con mayor fitness. Una manera de lograr esto es mediante el uso de funciones de interpolación. Imaginemos que en una determinada iteración los elementos del conjunto A= {y1, . . . , ym}alcanzan los valores de fitness {f(y1), . . . , f(ym)}. Nos proponemos definir la nueva variable aleatoria o, más concretamente, su función de probabilidad (de densidad en el caso continuo o de masa en el caso discreto) asociada pnueva. Decimos que una función ginterpola afen Asi g(yi) = f(yi)∀i∈ {1, . . . , m}. Es evidente que si pnueva es proporcional a alguna función gque interpola a fen A, tendremos lo que perseguimos: que los elementos de Acon mayor fitness cuenten con una mayor probabilidad de ser seleccionados para la siguiente iteración y viceversa. Además, si exigimos un cierto grado de continuidad ag, los puntos próximos a cada uno de ellos serán escogidos también con probabilidad similar. Esto será útil bajo la premisa de que, hasta cierto punto, el fitness de individuos parecidos es parecido, en la que se fundamentan todos los algoritmos evolutivos. En la Sección 5.2.1 se analizan las funciones de interpolación ideadas durante el desarrollo del trabajo, la motivación detrás de cada una y sus ventajas e inconvenientes. La Sección 5.2.2 está dedicada al otro punto fundamental del Algoritmo 3: precisamos obtener muestras según una distribución concreta arbitraria, con algún 5.2. Definición de otras particularizaciones 21 método implementable en la práctica. Se presentan y se analizan también distintas formas de lograrlo. 5.2.1. Funciones de interpolación La interpolación polinómica es, sin duda, la más estudiada en el campo del análisis numérico. También fue la primera considerada durante el desarrollo de este trabajo. En el caso unidimensional es ampliamente conocido el polinomio interpolador de Newton o de Lagrange [46, Capítulo 6]. Hecht et al. [20] han desarrollado todo un marco teórico para extender el objeto mencionado a dimensión arbitraria. Además, la librería Minterpy [3] implementa en lenguaje Python todos los algoritmos diseñados en el artículo citado y las estructuras de datos necesarias. Tratar con funciones de interpolación gpolinómicas supone grandes ventajas debido a su regularidad. En particular, al ser fácilmente integrables analíticamente, tenemos garantizado que la función de densidad proporcional p(x) = g(x) RΘg(x)dx puede ser calculada. Esta propiedad también ha permitido diseñar un método de muestreo específico que se detalla en la Sección 5.2.2. Lamentablemente, tras charlar con los desarrolladores de Minterpy, comprendimos que dicha librería no se ajusta adecuadamente a nuestro propósito. En primer lugar, no permite obtener un polinomio gque interpola a fen un conjunto fijado de puntos A. Por el contrario, el propio algoritmo selecciona los nodos de interpolación de entre todos los puntos del espacio multidimensional, pues estos deben cumplir la llamada propiedad de unisolvencia [20, Sección 2]. ¿Es por tanto completamente imposible utilizar Minterpy para obtener un polinomio que pase por un conjunto arbitrario de puntos? Lo cierto es que no, siempre y cuando se admita un determinado error de aproximación. La librería también ofrece una funcionalidad interesante: la regresión polinómica. En este caso, a partir de un conjunto A={y1, . . . , ym}de puntos, el algoritmo calcula los nodos unisolventes de forma que la función de regresión gobtenida (que será un polinomio de interpolación en dichos nodos) aproxima fen A. Es decir, ya no se cumple necesariamente que g(yi) = f(yi)∀i∈ {1, . . . , m}, pero esto no nos supone un inconveniente. Por la naturaleza de nuestro problema, únicamente necesitamos distinguir las regiones de mayor y menor fitness, sin preocuparnos por la exactitud. A pesar de todo, las primeras pruebas hicieron patente la necesidad de espacio en disco del orden de petabytes al considerar puntos en un espacio 30-dimensional, incluso restringiéndonos a polinomios de grado menor o igual que tres. De nuevo, los desarrolladores de la librería arrojaron algo de luz sobre lo que ocurría. La estructura de datos interna en la que se fundamentan prácticamente todos los algoritmos incluidos en Minterpy, los conjuntos de multiíndices (clase MultiIndexSet), aumenta 22 Capítulo 5. Algoritmo generalizado su tamaño exponencialmente con la dimensión del espacio. Por ello, la librería queda finalmente descartada para dimensión mayor que cinco. En realidad, tanto Minterpy como otras librerías de interpolación numérica persiguen un objetivo concreto: obtener una función que aproxime otra desconocida (como podría ser la de fitness) con el mayor grado de precisión posible en todo punto de su dominio. Esta tarea se ve afectada por el fenómeno denominado maldición de la dimensión, descrito por primera vez por Bellman [5] en 1957. La cantidad de puntos necesarios para aproximar adecuadamente una función crece exponencialmente con la dimensión del espacio. No obstante, recordemos que nuestro objetivo dista de conocer los valores de fitness asociados a todo punto de la región factible Θ. Por el contrario, pretendemos ir desplazándonos hacia mejores soluciones iterativamente, según el espíritu de cualquier algoritmo evolutivo. Podemos alejarnos completamente de la idea de la interpolación polinómica y emplear otra técnica, mucho más sencilla a priori, pero que ha resultado ser poderosa y versátil: la interpolación por medias ponderadas. Imaginemos el siguiente ejemplo en R2, ilustrado en la Figura 5.1. Ejemplo 5.2.1. Sean los puntos y1= (0,1), y2= (2,0), y3= (0,3), con fitness f(y1)=1, f(y2) = 4, f(y3) = 10. y1 y2 y3 y Figura 5.1: Representación gráfica del Ejemplo 5.2.1. Si definimos g(yi) = f(yi)para todo i∈ {1,2,3}, ¿qué valor de gasignamos, por ejemplo, al punto y= (0,0)? Una forma de conseguir que dicho valor sea parecido al de los puntos cercanos es emplear simplemente la media de {f(y1), f(y2), f(y3)} ponderada por el inverso de la distancia de ya cada uno de ellos. Es decir, si d denota la métrica euclídea (o Manhattan, que en este caso coinciden), g(y) = 1 d(y,y1)f(y1) + 1 d(y,y2)f(y2) + 1 d(y,y3)f(y3) 1 d(y,y1)+1 d(y,y2)+1 d(y,y3) = 1 1·1 + 1 2·4 + 1 3·10 1 1+1 2+1 3 =38 11. El ejemplo anterior nos conduce de manera natural a la definición de la fórmula general de función de interpolación por medias ponderadas. 5.2. Definición de otras particularizaciones 23 Definición 5.2.2. Sean el subconjunto A={y1, . . . , ym}de Rn, con valores de fitness {f(y1), . . . , f(ym)}, y una distancia d, la función de interpolación por medias ponderadas gMP :Rn→Rse define como gMP(y) =    f(y)si y∈A, Pm i=1 w(y, yi)f(yi)si y /∈A, donde w(y, z) = 1 d(y,z) Pm i=1 1 d(y,yi) . Una ventaja importante de gMP es que su definición se puede modificar fácilmente para obtener variantes que no emplean todos los puntos de A, sino únicamente los más interesantes en cada momento (los de mayor fitness, los más cercanos a y, etcétera). Esto, como veremos en el Capítulo 6, servirá para dar lugar a particularizaciones más similares a unos o a otros algoritmos evolutivos clásicos. 5.2.2. Técnicas de muestreo Una vez hemos definido una función de probabilidad pnueva : Θ →[0,1], nos preguntamos cómo obtener muestras según su distribución inducida, que es a priori arbitraria y desconocida. Un algoritmo ampliamente utilizado por su versatilidad es el llamado método de aceptación-rechazo [47, Sección 2.3], que es aplicable siempre y cuando Θesté acotado. Este aprovecha que es sencillo en la práctica obtener muestras siguiendo una distribución uniforme cualquiera. En primer lugar, lo presentamos para el caso univariable en el Algoritmo 4. Consideramos, por tanto, que Θ⊂(a, b)⊂R. Además, denotamos por U(a, b)la distribución uniforme en el intervalo (a, b). Algoritmo 4 Método de aceptación-rechazo univariable. Es preciso definir la región factible Θ, la función de probabilidad py el tamaño muestral deseado k∈N. El valor ysup >0es una cota superior del valor máximo de pen Θ. S← ∅ mientras |S|< k hacer Obtener una muestra {x} ∼ U(a, b). Obtener una muestra {y} ∼ U(0, ysup). si p(x)> y entonces S←S∪ {x} fin si fin mientras Si carecemos de más información, siempre podemos definir ysup = 1, aunque estimaciones más finas darán mejores resultados en términos de número necesario de iteraciones. La Figura 5.2 muestra una interpretación visual de la corrección del algoritmo. Los puntos (x, y)se distribuyen uniformemente en el rectángulo (a, b)× 24 Capítulo 5. Algoritmo generalizado (0, ysup). Únicamente aquellos que quedan debajo de la gráfica de p, es decir, verifican p(x)> y son aceptados y añadidos a S. Intuitivamente, en las regiones de Θcon mayor valor de pse concentrarán más puntos aceptados y al contrario. ab 0 ysup p(x) ysup Figura 5.2: Puntos cuya coordenada xes aceptada (círculos azules) y rechazada (cruces rojas). La extensión del método de aceptación-rechazo a dimensión arbitraria resulta muy natural, y se muestra en el Algoritmo 5. Consideramos en este caso que Θ⊂ (a1, b1)×· · ·×(an, bn)⊂Rn. Sabemos que si X= (X1, . . . , Xn)es un vector aleatorio definido en el rectángulo R= (a1, b1)× · · · × (an, bn)yX∼U(R), entonces Xi∼ U(ai, bi)para todo i∈ {1, . . . , n}y son independientes. Por tanto, podemos obtener muestras x= (x1, . . . , xn)componente a componente. Algoritmo 5 Método de aceptación-rechazo multivariable. Es preciso definir la región factible Θ, la función de probabilidad py el tamaño muestral deseado k∈N. El valor ysup >0es una cota superior del valor máximo de pen Θ. S← ∅ mientras |S|< k hacer para i∈ {1, . . . , n}hacer Obtener una muestra {xi} ∼ U(ai, bi). fin para Obtener una muestra {y} ∼ U(0, ysup). si p(x1, . . . , xn)> y entonces S←S∪ {(x1, . . . , xn)} fin si fin mientras El método de aceptación-rechazo es muy sencillo de implementar en la práctica, y se puede emplear con cualquier función de probabilidad y en dimensión arbitraria. Además, mantiene sus propiedades al sustituir ppor otra función proporcional, por lo que no es necesario que pnueva esté normalizada en [0,1] (siempre y cuando seamos capaces de estimar una cota superior, que ya no tendría por qué ser 1). Esto será 5.2. Definición de otras particularizaciones 25 especialmente útil al tratar con funciones de interpolación con integral desconocida o computacionalmente prohibitiva. Por supuesto, no todo son ventajas. El principal inconveniente de este algoritmo aparece cuando la función de interpolación considerada asigna valores muy próximos a cero a una región demasiado grande de su dominio. En este caso, el número de iteraciones del bucle externo crece considerablemente, debido a la alta cantidad de rechazos. En las situaciones más extremas el método de aceptación-rechazo es, por tanto, indistinguible de una búsqueda exhaustiva. En el caso concreto en que pnueva es una función polinómica, proponemos un algoritmo de muestreo totalmente distinto que no ve incrementado su coste en tiempo en presencia de grandes regiones de ceros. Describimos su funcionamiento por medio de un ejemplo en dimensión n= 2. Ejemplo 5.2.3. Sean Θ = [0,100] ×[0,100] y pnueva(x, y) = g(x, y) RΘg(x, y)d(x, y), donde g(x, y) = x3+ 2x2y+xy + 5xy2+ 2xy3+y4. Nos proponemos obterner una muestra según la distribución de probabilidad que tiene a pnueva por función de densidad. En primer lugar observamos que el denominador es fácilmente calculable: ZΘ g(x, y)d(x, y) = Z100 0Z100 0 [x3+ 2x2y+xy + 5xy2+ 2xy3+y4]dydx =Z100 0 [x3y+x2y2+1 2xy2+5 3xy3+1 2xy4+1 5y5] 100 y=0 dx =Z100 0 [100x3+ 104x2+ 5 ·103x+5·106 3x+ 5 ·107x+ 2 ·109]dx = [25x4+104 3x3+155015 ·103 6x2+ 2 ·109x] 100 x=0 =1392575000000 3. Si lo analizamos intuitivamente, únicamente hemos empleado operaciones sencillas con los exponentes y hemos recorrido en dos ocasiones todos los monomios. Llamemos D= 1392575000000/3. Como ya hemos mencionado, es factible en la práctica generar números según una distribución uniforme. Consideramos {r, s} ∼ U(0,1), y supongamos por un momento que conocemos el valor x0∈[0,100] tal que 1 DZx0 0Z100 0 g(x, y)dxdy=r, esto es, reutilizando los cálculos anteriores, aquel tal que 25x4 0+104 3x3 0+155015 ·103 6x2 0+ 2 ·109x0=rD. 32 Capítulo 6. Implementación práctica Algoritmo 6 AEMP generalizado. Es necesario definir la región factible Θ, el número de iteraciones k, el tamaño de las muestras m, la distancia dy la función de fitness fnormalizada en [0,1]. También se debe elegir el operador de selección sel(B, f, m), que devuelve msoluciones de Ben función de su valor de fitness. A← {s1, . . . , sm} ∼ U(Θ) para kiteraciones hacer S← ∅ mientras |S|< m hacer Obtener una muestra {x} ∼ U(Θ). Obtener una muestra {y} ∼ U(0,1). Elegir el subconjunto Bx⊂Ade soluciones privilegiadas. si gBx(x)> y entonces S←S∪ {x} fin si fin mientras A←sel(A∪S, f, m) fin para donde wB(y, z) = 1 d(y,z) Pb∈B 1 d(y,b) . Nótese que la media ponderada de los valores de fitness ahora se realiza únicamente sobre los elementos de un subconjunto de Aque depende de x. Esto es una diferencia fundamental con los EDA, donde recordemos que en cada iteración seleccionábamos un cierto B⊂Ade soluciones privilegiadas. Este determinaba también la probabilidad p(x|B)de obtener cada nuevo punto xpor muestreo, pero era el mismo para todos ellos. Si Bxestá formado por los dos puntos más próximos a x, entonces obtenemos una versión del AEMP más parecida a GA. Si, por el contrario, dicho subconjunto contiene el elemento de Amás cercano a xy el óptimo (el de mayor fitness), obtenemos otra versión más parecida a PSO. Lo ilustramos en la Figura 6.3 y en la Figura 6.4 respectivamente. Figura 6.3: Conjunto A(rojo) y punto x(azul). Los dos elementos de Amás próximos axforman Bx. 6.2. Algoritmos empleados 33 Figura 6.4: Conjunto A(elemento de mayor fitness ven verde, resto en rojo) y punto x(azul). El subconjunto Bxestá formado por vy el punto de Amás próximo a x. La intuición detrás de esta afirmación es sencilla. Al aplicar un algoritmo genético, la probabilidad de considerar una cierta solución en la iteración i+ 1 está fuertemente relacionada con la de haber seleccionado dos posibles progenitores (es decir, soluciones próximas a ella) en la iteración i. Por otro lado, en la iteración i+1 de PSO es muy probable que las partículas consideradas sean cercanas a las de la iteración iy que, además, se aproximen a la solución óptima del enjambre. Denominamos ambas versiones AEMP-GA y AEMP-PSO respectivamente. Es importante notar que no estamos tratando de imitar a la perfección los algoritmos clásicos, sino de definir y probar distintas variaciones de AEMP que capturen la esencia de cómo dichos algoritmos privilegian unas soluciones u otras, de manera fácil y directa. Además, hibridarlos ahora resulta trivial. Basta con, por ejemplo, definir Bxcomo el conjunto que contiene los dos elementos de Amás próximos a x y aquel con mayor valor de fitness. Llamamos a esta tercera variante AEMP-HIB en lo que resta de capítulo. Definir una particularización de AEMP similar a ACO ha supuesto un gran reto durante el desarrollo de este trabajo. El propio Marco Dorigo, padre de los algoritmos de colonia de hormigas, presenta la que es, según su criterio, la mejor forma de generalizar los mismos en [50]. Su objetivo es poder aplicar ACO a problemas de optimización continua, y aporta dos ideas fundamentales que podemos utilizar. En primer lugar, Socha y Dorigo [50] proponen descartar completamente el grafo de construcción y simplemente llevar un archivo, es decir, una tabla con las mejores soluciones y su valor de fitness. Esto es exactamente lo que se hace en AEMP, donde el conjunto Atoma el papel de archivo. Se basan en que, cuando una hormiga decide qué nueva componente añadir a su estado (solución parcial), simplemente está tomando una muestra a partir de una variable aleatoria discreta (que asigna distintas probabilidades positivas a cada arista), y esta variable podría sustituirse por otra totalmente arbitraria. Lo único verdaderamente importante según los autores es que cada solución se construya componente a componente, tal y como se hace en AEMP. La segunda modificación fundamental se centra en la actualización de feromonas. Como ya no existe el grafo de construcción, los autores proponen simplemente 34 Capítulo 6. Implementación práctica olvidar (es decir, eliminar del archivo) las Mpeores soluciones tras haber generado Mnuevas. Por tanto, para obtener una versión de AEMP más parecida a ACO, podemos definir el subconjunto Bx=Apara todo punto x. Llamamos a esta última particularización AEMP-ACO. 6.3. Resultados obtenidos En la Sección 6.3.1 y la Sección 6.3.2 se muestran y analizan los resultados obtenidos tras aplicar los seis algoritmos a diversas instancias del problema de la cobertura de vértices y de la mochila respectivamente. En todos los casos se ha empleado la distancia más natural en espacios binarios, la de Hamming, definida a continuación. Nótese que dHamming(x, y)simplemente mide el número de componentes de xque habría que invertir (cambiar 1por 0o viceversa) para llegar ay. Definición 6.3.1. Dados x, y ∈ {0,1}n, la distancia de Hamming dHamming entre xeyse define como dHamming(x, y) = n X i=1 |xi−yi|. Además, el operador de selección utilizado es el de truncamiento, es decir, en cada iteración se eligen las msoluciones con mayor fitness y se descartan las demás. 6.3.1. Problema de la cobertura de vértices Se han probado todos los algoritmos con distintas instancias del problema VERTEX COVER tomadas de [48]. En concreto, se han empleado las mostradas en la Tabla 6.1, por contar con un tamaño (número de vértices) ajustado a los objetivos de este trabajo. Algunas han sido preprocesadas para que todas representen un grafo no dirigido y no valorado con vértices numerados de 1an. Instancia n ENZYMES-G102 42 aves-sparrow-social-2010 40 aves-sparrowlyon-flock-season2 46 cat-mixed-species-brain-1 65 Tabla 6.1: Las cuatro instancias empleadas y su número de vértices n. Los resultados obtenidos en 20 ejecuciones independientes, con 1000 iteraciones y poblaciones de 300 individuos para los seis algoritmos considerados y la instancia ENZYMES-G102 se muestran en la Tabla 6.2. Cabe destacar que, aunque gran parte de las ejecuciones arrojan los mismos valores de fitness máximo, las soluciones en las que este se alcanza son en general distintas. Esto resalta la naturaleza no determinista de los algoritmos empleados. El valor de fitness máximo encontrado ha sido 17, y en la Figura 6.5 se muestra un ejemplo de solución con dicho valor, que tal y como hemos mencionado no es única. 6.3. Resultados obtenidos 35 1 2 3 4 5 6 7 8 9 10 1112 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 Figura 6.5: Cobertura de 42 −17 = 25 vértices para la instancia ENZYMES-G102. 36 Capítulo 6. Implementación práctica GA PSO AEMP-GA AEMP-PSO AEMP-HIB AEMP-ACO 1 17 16 16 16 16 16 2 17 16 16 16 16 17 3 17 16 16 16 16 16 4 17 16 16 16 16 16 5 17 16 16 16 16 16 6 17 16 16 16 16 16 7 17 16 15 16 16 16 8 17 16 16 17 16 16 9 17 16 16 16 16 17 10 17 16 16 16 16 16 11 17 16 16 16 16 16 12 17 17 16 15 16 16 13 17 16 16 16 16 16 14 17 16 17 16 16 16 15 17 17 16 16 16 16 16 17 17 16 16 16 16 17 17 16 16 16 16 16 18 17 16 16 16 16 16 19 17 17 16 16 16 16 20 17 17 15 16 16 16 Tabla 6.2: Resultados para ENZYMES-G102. Como se puede observar, los resultados son por lo menos aceptables. No forma parte de los objetivos de este trabajo diseñar algoritmos que superen en rendimiento a las implementaciones optimizadas de los algoritmos tradicionales. No obstante, es buena práctica emplear pruebas estadísticas para comparar de forma justificada los resultados. Tomamos como hipótesis nula H0que no existen diferencias significativas entre los seis algoritmos, es decir, que los 20 resultados de cada uno podrían ser muestras de la misma variable aleatoria subyacente, y realizamos un contraste estadístico. Debemos utilizar las llamadas pruebas no paramétricas, puesto que estas no asumen propiedades como normalidad, independencia u homogeneidad de varianzas. El lector interesado en profundizar en esta clase de pruebas y cómo aplicarlas puede encontrar toda la información en [9]. Tras emplear el test de Quade, una variación del más comúnmente conocido test de Friedman que ajusta los pesos a la hora de realizar la clasificación en función de la dificultad del problema, obtenemos que H0puede ser rechada. Es decir, sí existen diferencias significativas entre los seis algoritmos. Sin embargo, tras aplicar lo que se conoce como un método post-hoc para comparar los mismos dos a dos (en concreto, el de Holm, de los recomendados en [9]), resulta que únicamente existen diferencias significativas entre los resultados de GA y AEMP-GA, y no entre los otros 14 pares, con un nivel de significación α= 0,05. En la Tabla 6.3 y en la Tabla 6.4 se muestran los resultados obtenidos para las instancias aves-sparrow-social-2010 yaves-sparrowlyon-flock-season2 res- 6.3. Resultados obtenidos 37 pectivamente. Resulta llamativo que absolutamente todas han alcanzado el mismo valor de fitness máximo, 10 en el caso de la primera instancia y 13 en el de la segunda. Quizá incluso demasiado llamativo. No obstante, de nuevo, dichos valores se han alcanzado en soluciones diferentes. Es decir, sí se ha producido una adecuada exploración del espacio de búsqueda. Evidentemente, podemos afirmar que no existen diferencias significativas entre los seis algoritmos sin necesidad de emplear pruebas no paramétricas. GA PSO AEMP-GA AEMP-PSO AEMP-HIB AEMP-ACO 1 10 10 10 10 10 10 2 10 10 10 10 10 10 3 10 10 10 10 10 10 4 10 10 10 10 10 10 5 10 10 10 10 10 10 6 10 10 10 10 10 10 7 10 10 10 10 10 10 8 10 10 10 10 10 10 9 10 10 10 10 10 10 10 10 10 10 10 10 10 11 10 10 10 10 10 10 12 10 10 10 10 10 10 13 10 10 10 10 10 10 14 10 10 10 10 10 10 15 10 10 10 10 10 10 16 10 10 10 10 10 10 17 10 10 10 10 10 10 18 10 10 10 10 10 10 19 10 10 10 10 10 10 20 10 10 10 10 10 10 Tabla 6.3: Resultados para aves-sparrow-social-2010. Finalmente, en la Tabla 6.5 se presentan los resultados obtenidos para la instancia cat-mixed-species-brain-1. Si aplicamos de nuevo el test de Quade y el método post-hoc de Holm, obtenemos que únicamente existen diferencias significativas entre GA y AEMP-ACO, entre GA y AEMP-HIB y entre GA y AEMP-GA, con un nivel de significación 0,05. Resulta reseñable que las pruebas estadísticas no encuentran diferencias significativas entre GA y AEMP-PSO, ni entre PSO y ninguna de las variantes de AEMP. 6.3.2. Problema de la mochila Se han empleado cuatro instancias de KNAPSACK de dimensión 50 tomadas de [40]. Concretamente, las cuatro primeras del directorio 00Uncorrelated/n00050/R01000. Las nombramos kn0,kn1,kn2 ykn3 por simplicidad. En esta ocasión, se han necesitado 5000 iteraciones por ejecución para alcanzar la convergencia de los algoritmos. 38 Capítulo 6. Implementación práctica GA PSO AEMP-GA AEMP-PSO AEMP-HIB AEMP-ACO 1 13 13 13 13 13 13 2 13 13 13 13 13 13 3 13 13 13 13 13 13 4 13 13 13 13 13 13 5 13 13 13 13 13 13 6 13 13 13 13 13 13 7 13 13 13 13 13 13 8 13 13 13 13 13 13 9 13 13 13 13 13 13 10 13 13 13 13 13 13 11 13 13 13 13 13 13 12 13 13 13 13 13 13 13 13 13 13 13 13 13 14 13 13 13 13 13 13 15 13 13 13 13 13 13 16 13 13 13 13 13 13 17 13 13 13 13 13 13 18 13 13 13 13 13 13 19 13 13 13 13 13 13 20 13 13 13 13 13 13 Tabla 6.4: Resultados para aves-sparrowlyon-flock-season2. En la Tabla 6.6 se muestran los resultados obtenidos para la instancia kn0. El test de Quade y el método post-hoc de Holm únicamente hallan diferencias significativas entre GA y AEMP-GA, entre GA y AEMP-ACO y entre GA y AEMP-PSO, con un nivel de significación α= 0,05. Es decir, en este caso concreto AEMP-HIB supera en rendimiento a los dos algoritmos que hibrida y resulta indistinguible del mejor, que de nuevo es GA. En la Tabla 6.7 se presentan los resultados obtenidos para la instancia kn1. En esta ocasión sí hay diferencias significativas (con α= 0,05) entre GA y todas las variantes de AEMP, aunque nuevamente AEMP-HIB es la que mejor actúa de las cuatro. Los resultados para la instancia kn2 se muestran en la Tabla 6.8. Tras aplicar el test de Quade y el método de Holm, podemos observar algo curioso. A diferencia de los dos anteriores, en este caso no hay diferencias significativas entre AEMP-GA y GA (como sí las hay entre GA y las demás variantes de AEMP). Es decir, con un nivel de significación α= 0,05, el algoritmo genético clásico y la versión de AEMP más parecida a él actúan igual. Lo mismo se puede afirmar de AEMP-PSO y PSO. Por último, los resultados obtenidos para la instancia kn3 se presentan en la Tabla 6.9. En esta ocasión únicamente hay diferencias significativas (con α= 0,05) entre dos de los 15 pares de algoritmos: entre GA y AEMP-ACO y entre GA y AEMP-GA. 6.3. Resultados obtenidos 39 GA PSO AEMP-GA AEMP-PSO AEMP-HIB AEMP-ACO 1 12 11 11 11 11 11 2 12 11 11 11 11 10 3 12 11 10 11 10 11 4 12 11 10 11 11 11 5 12 11 10 10 11 10 6 12 11 10 10 10 11 7 12 10 10 10 11 10 8 12 11 10 11 10 10 9 12 11 10 11 11 10 10 12 12 10 11 11 10 11 12 11 11 10 11 10 12 12 11 11 11 10 10 13 12 10 10 11 10 11 14 12 10 11 11 10 10 15 12 11 11 10 10 10 16 12 10 11 11 10 10 17 12 11 10 11 10 11 18 12 10 11 11 10 10 19 12 10 11 11 10 10 20 12 11 10 10 10 10 Tabla 6.5: Resultados para cat-mixed-species-brain-1. GA PSO AEMP-GA AEMP-PSO AEMP-HIB AEMP-ACO 1 20783 19213 18931 19176 19076 18703 2 20716 19250 19042 19121 18888 19218 3 20776 19318 19204 18976 19003 19550 4 20806 19185 19110 19601 19220 18974 5 20877 19711 19098 19535 18993 19007 6 20741 19602 19025 19418 19693 18779 7 20936 19418 19277 19320 18798 19433 8 20817 19865 18819 19308 19239 18836 9 20884 19678 19079 18929 19213 19133 10 20842 19824 18846 19097 18923 18930 11 20768 19396 19490 18848 18896 19144 12 20723 19290 19106 19168 19291 18743 13 20612 19148 19107 18955 19113 19274 14 20653 20403 18927 18988 19147 19146 15 20984 19365 19186 18679 19202 19040 16 20864 19479 19609 18937 19001 19118 17 20583 19427 18762 18978 18897 19093 18 20842 20056 18884 19487 19485 19426 19 20831 19483 19085 19167 19262 18877 20 20723 19609 19452 19060 19018 18987 Tabla 6.6: Resultados para kn0. 40 Capítulo 6. Implementación práctica GA PSO AEMP-GA AEMP-PSO AEMP-HIB AEMP-ACO 1 19459 18062 17460 17020 17490 17154 2 19441 17951 17046 17219 17738 18273 3 19471 17556 17274 17449 17557 17444 4 19599 17682 17612 17039 17201 17494 5 19604 17862 17270 16860 17349 17718 6 19411 17850 17613 17169 17852 17328 7 19537 17991 16990 17226 16996 16897 8 19836 17604 17368 17114 17206 16914 9 19511 18464 17386 17257 17417 17273 10 19676 17537 17145 17490 17134 17250 11 19496 18098 17486 17197 17427 17395 12 19432 17940 17392 17543 17089 17278 13 19422 17637 17248 17596 17150 17197 14 19836 17508 17112 17387 17320 17914 15 19548 17613 17353 17535 17305 17340 16 19809 17703 17524 17361 18573 17329 17 19836 17555 17149 17840 17185 16953 18 19334 17769 17172 18293 17396 17311 19 19666 18041 17244 18387 17519 17282 20 19376 17717 17274 17295 17274 17796 Tabla 6.7: Resultados para kn1. GA PSO AEMP-GA AEMP-PSO AEMP-HIB AEMP-ACO 1 20370 19402 18960 18457 18555 18732 2 20339 19139 18640 18558 18611 18221 3 20360 19562 18959 18506 18814 18779 4 20243 19282 18821 18986 18795 18740 5 20382 19246 18553 18681 18608 18509 6 20127 18992 18720 18549 18483 18563 7 20285 19227 19003 18641 18549 18414 8 20390 19122 18780 18469 18621 18823 9 20471 19483 18599 18845 18496 18595 10 20419 19188 18743 18743 18937 18739 11 20371 19057 18637 18390 18822 18587 12 20436 18848 18980 18584 18778 19034 13 20410 18862 18503 18800 18755 18782 14 20361 19165 19001 18452 18819 18834 15 20553 19158 18903 18688 18448 18539 16 20317 18847 18816 18597 18505 18721 17 20373 19293 19093 18774 18563 18695 18 20401 18945 18842 18612 18778 18754 19 20507 18735 18586 18672 18603 18520 20 20507 19128 18781 19265 18672 19125 Tabla 6.8: Resultados para kn2. 6.3. Resultados obtenidos 41 GA PSO AEMP-GA AEMP-PSO AEMP-HIB AEMP-ACO 1 19450 18208 17751 18568 18190 17589 2 19547 18108 17839 17764 18595 17794 3 19651 18366 17957 17846 18039 17902 4 19691 18224 17847 18375 17925 17738 5 19513 17866 17626 17973 17828 17868 6 19609 18359 17707 17972 18088 17855 7 19640 18077 18171 17949 18247 17653 8 19691 17775 17900 18517 17565 17814 9 19605 18351 17891 18079 17991 17783 10 19684 18174 18013 18090 18674 17840 11 19736 18570 17881 17727 18261 17671 12 19596 18454 18021 17905 18065 18045 13 19571 18298 17987 18167 18127 17866 14 19530 18798 18057 17794 17869 17806 15 19574 17787 17896 17695 18009 17676 16 19487 17894 17808 18283 17624 17980 17 19523 18281 17455 17674 17904 17584 18 19677 18441 17795 18353 17763 17886 19 19581 17964 17779 17884 17937 18080 20 19709 18327 17700 17831 17799 17972 Tabla 6.9: Resultados para kn3. 48 Capítulo 7. Conclusiones y trabajo futuro proposed to handle polynomials and use it to solve continuous optimization problems in dimension five or lower. There are highly interesting problems in low-dimensional spaces, such as those in the field of so-called black-box optimization [41]. In these cases, it is a priority to minimize the number of evaluations for the fitness function, as these evaluations are computationally expensive or the function itself is unknown. Bibliografía [1] Ab. Rashid, M. Tiki-taka algorithm: a novel metaheuristic inspired by football playing style. Engineering Computations, vol. 38(1), páginas 313–343, 2020. [2] Abbass, H. MBO: marriage in honey bees optimization – a Haplometrosis polygynous swarming approach. En Proceedings of the 2001 Congress on Evolutionary Computation (IEEE Cat. No.01TH8546), vol. 1, páginas 207–214. 2001. [3] Acosta, U. H.,Veettil, S. K. T.,Wicaksono, D.,Schreiber, J.,Michelfeit, J.,Hoffman, N.,Schmerler, S. yChandrashekar, V. Repositorio de la librería Minterpy (Multivariate Interpolation in Python). 2021. Disponible en https://github.com/casus/minterpy (último acceso, 2024). [4] Baluja, S. Population-Based Incremental Learning: A Method for Integrating Genetic Search Based Function Optimization and Competitive Learning. Informe Técnico CMU-CS-94-163, Computer Science School, Carnegie Mellon University, 1994. [5] Bellman, R. Dynamic programming. Princeton University Press, Princeton, NJ, USA, 1957. [6] Chu, S.-C.,Tsai, P.-w. yPan, J.-S. Cat Swarm Optimization. En PRICAI 2006: Trends in Artificial Intelligence (editado por Q. Yang y G. Webb), páginas 854–858. Springer Berlin Heidelberg, 2006. [7] Chu, Y.,Mi, H.,Liao, H.,Ji, Z. yWu, Q. H. A Fast Bacterial Swarming Algorithm for high-dimensional function optimization. En 2008 IEEE Congress on Evolutionary Computation (IEEE World Congress on Computational Intelligence), páginas 3135–3140. 2008. [8] Darwin, C. On the Origin of Species. John Murray, 1859. [9] Derrac, J.,García, S.,Molina, D. yHerrera, F. A practical tutorial on the use of nonparametric statistical tests as a methodology for comparing evolutionary and swarm intelligence algorithms. Swarm and Evolutionary Computation, vol. 1(1), páginas 3–18, 2011. [10] Dorigo, M.,Maniezzo, V. yColorni, A. Positive Feedback as a Search Strategy. Informe Técnico 91-016, Dipartimento di Elettronica, Politecnico di Milano, Italia, 1999. 49 50 BIBLIOGRAFÍA [11] Dorigo, M.,Montes de Oca, M. A. yEngelbrecht, A. Particle swarm optimization. Scholarpedia, vol. 3(11), página 1486, 2008. Disponible en http: //www.scholarpedia.org/article/Particle_swarm_optimization (último acceso, 2024). [12] Dorigo, M. yStützle, T. Ant Colony Optimization. Bradford Company, 2004. [13] Eiben, A. E. ySmith, J. E. Introduction to Evolutionary Computing. Springer, 2015. [14] Engelbrecht, A. P. Fundamentals of Computational Swarm Intelligence. Wiley, 2005. [15] Gad, A. F. Pygad: An intuitive genetic algorithm python library. Multimedia Tools and Applications, páginas 1–14, 2023. [16] Goldberg, D. E. yLingle, R. Alleles, loci, and the traveling salesman problem. En Proceedings of the first international conference on genetic algorithms and their applications, páginas 154–159. Psychology Press, 1985. [17] Hagberg, A. A.,Schult, D. A. ySwart, P. J. Exploring network structure, dynamics, and function using NetworkX. En Proceedings of the 7th Python in Science Conference (SciPy2008) (editado por G. Varoquaux, T. Vaught y J. Millman), páginas 11–15. 2008. [18] Hansen, N. The CMA Evolution Strategy: A Tutorial. ArXiv, 2016. Disponible en https://arxiv.org/pdf/1604.00772 (último acceso, 2024). [19] Hansen, N. yOstermeier, A. Completely Derandomized Self-Adaptation in Evolution Strategies. Evolutionary Computation, vol. 9(2), páginas 159–195, 2001. [20] Hecht, M.,Gonciarz, K.,Michelfeit, J.,Sivkin, V. ySbalzarini, I. F. Multivariate Interpolation on Unisolvent Nodes – Lifting the Curse of Dimensionality. ArXiv, 2020. Disponible en https://arxiv.org/pdf/2010. 10824.pdf (último acceso, 2024). [21] Holland, J. H. Hidden Order, How Adaptation Builds Complexity. AddisonWesley Publishing Company, 1995. [22] Jung, S. H. Queen-bee evolution for genetic algorithms. Electronics Letters, vol. 39, páginas 575–576(1), 2003. [23] Kao, Y.-T. yZahara, E. A hybrid genetic algorithm and particle swarm optimization for multimodal functions. Applied Soft Computing, vol. 8(2), páginas 849–857, 2008. [24] Kaveh, A. yFarhoudi, N. A new optimization method: Dolphin echolocation. Advances in Engineering Software, vol. 59, páginas 53–70, 2013. BIBLIOGRAFÍA 51 [25] Kennedy, J. Swarm Intelligence. En Handbook of Nature-Inspired and Innovative Computing, Integrating Classical Models with Emerging Technologies (editado por A. Zomaya), páginas 187–219. 2006. [26] Kennedy, J. yEberhart, R. Particle swarm optimization. En Proceedings of ICNN’95 - International Conference on Neural Networks, vol. 4, páginas 1942–1948. 1995. [27] Larrañaga, P. An Introduction to Probabilistic Graphical Models. En Estimation of Distribution Algorithms: A New Tool for Evolutionary Computation (editado por P. Larrañaga y J. A. Lozano), páginas 27–56. Springer US, Boston, MA, 2002. [28] Larrañaga, P. A Review on Estimation of Distribution Algorithms. En Estimation of Distribution Algorithms: A New Tool for Evolutionary Computation (editado por P. Larrañaga y J. A. Lozano), páginas 57–100. Springer US, Boston, MA, 2002. [29] Larrañaga, P.,Etxeberria, R.,Lozano, J. yPeña, J. Optimization in continuous domains by learning and simulation of Gaussian networks. En Proceedings of the 2000 Genetic and Evolutionary Computation Conference Workshop Program, páginas 201–204. 2000. [30] López-Ibáñez, M.,Dubois-Lacoste, J.,Pérez Cáceres, L.,Stützle, T. yBirattari, M. The irace package: Iterated Racing for Automatic Algorithm Configuration. Operations Research Perspectives, vol. 3, páginas 43–58, 2016. [31] Miranda, L. J. V. PySwarms, a research-toolkit for Particle Swarm Optimization in Python. Journal of Open Source Software, vol. 3(21), 2018. [32] Mo, H. yXu, L. Magnetotactic bacteria optimization algorithm for multimodal optimization. En 2013 IEEE Symposium on Swarm Intelligence (SIS), páginas 240–247. 2013. [33] Molina, D.,Poyatos, J.,Del Ser, J.,García, S.,Hussain, A. yHerrera, F. Comprehensive Taxonomies of Natureand Bio-inspired Optimization: Inspiration versus Algorithmic Behavior, Critical Analysis Recommendations. Cognitive Computation, vol. 12, páginas 897––939, 2020. [34] Mühlenbein, H. yPaaß, G. From recombination of genes to the estimation of distributions I. Binary parameters. En Parallel Problem Solving from Nature — PPSN IV (editado por H.-M. Voigt, W. Ebeling, I. Rechenberg y H.- P. Schwefel), páginas 178–187. Springer Berlin Heidelberg, Berlin, Heidelberg, 1996. [35] Mühlenbein, H. yVoigt, H.-M. Gene Pool Recombination in Genetic Algorithms. En Meta-Heuristics: Theory and Applications (editado por I. H. Osman y J. P. Kelly), páginas 53–62. Springer, 1996. 52 BIBLIOGRAFÍA [36] Muller, S.,Marchetto, J.,Airaghi, S. yKournoutsakos, P. Optimization based on bacterial chemotaxis. IEEE Transactions on Evolutionary Computation, vol. 6(1), páginas 16–29, 2002. [37] Muthiah-Nakarajan, V. yNoel, M. M. Galactic Swarm Optimization: A new global optimization metaheuristic inspired by galactic motion. Applied Soft Computing, vol. 38, páginas 771–787, 2016. [38] Mühlenbein, H. The equation for response to selection and its use for prediction. Evolutionary Computation, vol. 5(3), páginas 303–346, 1997. [39] Nguyen, H. T. yBhanu, B. Zombie Survival Optimization: A swarm intelligence algorithm inspired by zombie foraging. En Proceedings of the 21st International Conference on Pattern Recognition (ICPR2012), páginas 987–990. 2012. [40] Onoue, Y. Repositorio de la librería kplib. 2013. Disponible en https: //github.com/likr/kplib (último acceso, 2024). [41] Pardalos, P. M.,Rasskazova, V. yVrahatis, M. N., editores. Black Box Optimization, Machine Learning, and No-Free Lunch Theorems. Springer, 2021. [42] Pelikan, M.,Goldberg, D. E. yCantú-Paz, E. BOA: the Bayesian optimization algorithm. En Proceedings of the 1st Annual Conference on Genetic and Evolutionary Computation, GECCO’99, páginas 525—-532. Morgan Kaufmann Publishers Inc., San Francisco, CA, USA, 1999. [43] Purnomo, H. D. yWee, H.-M. Soccer Game Optimization: An Innovative Integration of Evolutionary Algorithm and Swarm Intelligence Algorithm. En Meta-Heuristics Optimization Algorithms in Engineering, Business, Economics, and Finance (editado por P. M. Vasant), páginas 386–420. IGI Global, 2013. [44] Razmjooy, N.,Khalilpour, M. yRamezani, M. A New Meta-Heuristic Optimization Algorithm Inspired by FIFA World Cup Competitions: Theory and Its Application in PID Designing for AVR System. Journal of Control, Automation and Electrical Systems, vol. 27, página 419–440, 2016. [45] Rechenberg, I. Cybernetic Solution Path of an Experimental Problem. RAELT-1122. Royal Aircraft Establishment, 1965. [46] Infante del Río, J. A. yRey Cabezas, J. M. Métodos numéricos: Teoría, problemas y prácticas con MATLAB. Ciencia y Técnica. Ediciones Pirámide, sexta edición, 2022. [47] Robert, C. P. yCasella, G. Monte Carlo Statistical Methods. Springer Texts in Statistics. Springer-Verlag, segunda edición, 2004. BIBLIOGRAFÍA 53 [48] Rossi, R. A. yAhmed, N. K. The Network Data Repository with Interactive Graph Analytics and Visualization. En Proceedings of the Twenty-Ninth AAAI Conference on Artificial Intelligence. 2015. [49] Shiqin, Y.,Jianjun, J. yGuangxing, Y. A Dolphin Partner Optimization. En 2009 WRI Global Congress on Intelligent Systems, vol. 1, páginas 124–128. 2009. [50] Socha, K. yDorigo, M. Ant colony optimization for continuous domains. European Journal of Operational Research, vol. 185(3), páginas 1155–1173, 2008. [51] Sywerda, G. Uniform Crossover in Genetic Algorithms. En Proceedings of the Third International Conference on Genetic Algorithms (editado por J. D. Schaffer), páginas 2–9. Morgan Kaufmann Publishers Inc., 1989. [52] Sörensen, K. Metaheuristics—the metaphor exposed. International Transactions in Operational Research, vol. 22(1), páginas 3–18, 2015. [53] Tang, D.,Dong, S.,Jiang, Y.,Li, H. yHuang, Y. ITGO: Invasive tumor growth optimization algorithm. Applied Soft Computing, vol. 36, páginas 670– 698, 2015. [54] Teodorovic, D.,Lucic, P.,Markovic, G. yOrco, M. D. Bee Colony Optimization: Principles and Applications. En 2006 8th Seminar on Neural Network Applications in Electrical Engineering, páginas 151–156. 2006. [55] Vicente Arroyo, S. Repositorio de este trabajo. 2024. Disponible en https: //github.com/savicente2109/TFG-Infor (último acceso, 2024). [56] Zukhri, Z. yPaputungan, I. A Hybrid Optimization Algorithm based on Genetic Algorithm and Ant Colony Optimization. International Journal of Artificial Intelligence Applications, vol. 4(5), páginas 63–75, 2013.