Estrategia de Evolución de Adaptación de la Matriz de Covarianza — Covariance Matrix Adaptation Evolution Strategy (CMA-ES)
Contenido
Introducción
En el mundo de la optimización existen muchos algoritmos, pero buscamos los más potentes para resolver los problemas de optimización de nuestros robots comerciales. El CMA-ES (Covariance Matrix Adaptation Evolution Strategy) es uno de esos raros ejemplos en los que el rigor matemático se combina con la intuición biológica, creando un algoritmo que no solo resuelve problemas de optimización, sino que también aprende a comprender su estructura.
La historia del CMA-ES comenzó a finales de la década de 1990 en laboratorios de investigación en Alemania, donde Nikolaus Hansen y Andreas Ostermeier plantearon una pregunta fundamental: ¿es posible crear un algoritmo de optimización que no solo busque una solución, sino que se adapte a la geometría del problema? Los algoritmos evolutivos tradicionales generaban descendencia en regiones esféricas alrededor de los individuos progenitores, lo cual funcionaba bien para funciones simples, pero resultaba ineficaz para problemas complejos y mal condicionados. Así pues, vamos a echar un vistazo a este interesante algoritmo.
Implementación del algoritmo
Imagine buscar un tesoro en una isla de forma irregular. El método habitual consiste en buscar en todas las direcciones por igual, como si la isla fuera redonda. El CMA-ES estudia gradualmente la forma de la isla y comienza a buscar principalmente en aquellas direcciones donde la probabilidad de encontrar tesoros es mayor. Además, recuerda las rutas exitosas y usa esta memoria para planificar búsquedas futuras.
El CMA-ES se basa en una ecuación engañosamente simple: x_k ~ N(m, σ²C), Pero tras esta simplicidad se esconde una profunda estructura matemática. Cada símbolo aquí contiene información importante sobre el estado de la búsqueda: m es la mejor estimación actual sobre la ubicación del óptimo, σ es una medida de cuánto estamos dispuestos a alejarnos de los conocido, y C es una matriz de covarianza que codifica nuestra comprensión de la geometría de la función. El único cambio que podemos justificar es reemplazar la distribución normal por una distribución de potencia, lo cual significa que la implementación seguirá una fórmula modificada: x_k ~ PowerDist(m, σ²C). Esta modificación cambia la naturaleza de la exploración espacial (saltos más amplios), pero conserva la naturaleza adaptativa fundamental del algoritmo.
La matriz de covarianza C es el verdadero corazón del algoritmo. Esta comienza su vida como una sencilla matriz identidad que representa una distribución esférica, pero con cada iteración evoluciona, expandiéndose en áreas de rápida mejora y reduciéndose donde el progreso es lento. Gradualmente, la esfera se transforma en una elipse, y luego en un elipsoide alargado, idealmente orientado a lo largo de los contornos de la función que se está optimizando.
La principal innovación de CMA-ES es el concepto de trayectorias evolutivas. Se trata de una especie de memoria genética del algoritmo que recuerda no solo dónde se produjeron los aciertos, sino también cómo llegó el algoritmo a ellos. La trayectoria evolutiva acumula información sobre los sucesivos pasos exitosos, creando un vector direccional que apunta a las áreas de búsqueda más prometedoras. La segunda vía evolutiva cumple una función más sutil, a saber, controla el tamaño del paso. Si los pasos sucesivos del algoritmo están correlacionados, es decir, cada paso subsiguiente continúa la dirección del anterior, entonces el tamaño del paso aumenta: el algoritmo "siente" que se está moviendo en la dirección correcta. Si los pasos son aleatorios y no están correlacionados, el tamaño del paso disminuye; quizá el algoritmo ya esté cerca del óptimo y necesite una búsqueda más cuidadosa.
Detrás de la metáfora biológica de la evolución en CMA-ES se esconde un principio matemático riguroso: la maximización de la verosimilitud. El algoritmo se pregunta constantemente: ¿qué parámetros de la distribución hacen que los puntos de éxito observados sean los más probables? Esto transforma la optimización evolutiva de un método heurístico en un método basado en estadísticas.
Cada actualización de la matriz de covarianza C consta de dos componentes: la actualización rank-one basada en la trayectoria evolutiva y la actualización rank-μ que utiliza información de toda la población muestreada. La actualización rank-one proporciona estabilidad y la capacidad de capturar tendencias a largo plazo, mientras que la actualización rank-μ permite una rápida adaptación a la nueva información de la generación actual.
En el contexto del CMA-ES, se utiliza una función de Heaviside modificada que desempeña un papel importante en el mecanismo de detección de estancamiento del algoritmo. La función compara la longitud de la trayectoria evolutiva con la longitud esperada; si la trayectoria es demasiado corta, es señal de "estancamiento". Condiciones de desactivación (hsig = 0): el algoritmo "vaga" aleatoriamente y los pasos se cancelan entre sí, entonces el tamaño del paso es probablemente demasiado grande. Durante el estancamiento, la influencia de las actualizaciones rank-one disminuye y dependemos más de la información de toda la población. Condiciones de activación (hsig = 1): el algoritmo avanza en una dirección determinada, los pasos sucesivos están correlacionados y, por lo tanto, el tamaño del paso es adecuado a la situación actual.
Una de las propiedades más importantes del CMA-ES es su invariancia ante las transformaciones del espacio de búsqueda. El algoritmo resuelve con igual eficiencia la función f(x) y la función f(Ax + b), donde A es cualquier matriz no singular, mientras que b es cualquier vector de desplazamiento. Esto significa que CMA-ES es independiente de la elección del sistema de coordenadas. Esta invariancia no es accidental; supone una consecuencia directa del uso del principio de máxima verosimilitud y de la adaptación de la matriz de covarianza. El algoritmo detecta automáticamente un sistema de coordenadas natural para el problema donde los ejes coinciden con las direcciones principales de variación de la función.
La belleza de la teoría debe combinarse con la aplicabilidad práctica. El CMA-ES requiere una memoria de O(n²) para almacenar la matriz de covarianza y un cálculo de O(n³) para su propia descomposición, lo cual hace que el algoritmo sea aplicable a problemas con dimensionalidades de hasta varios cientos de variables. Para dimensionalidades grandes, se han desarrollado modificaciones especializadas: sep-CMA-ES se limita a matrices de covarianza diagonales, VkD-CMA utiliza dimensiones variables y LM-CMA aplica los principios de memoria limitada. La implementación de estos métodos actualmente escapa al alcance de este artículo, pero es posible que se retomen en artículos posteriores.

Figura 1. Ilustración del algoritmo CMA-ES
La ilustración muestra la evolución a través de las generaciones: tres instantáneas del estado del algoritmo (generaciones 1, 5 y 15) que muestran cómo la población converge gradualmente al óptimo.
Adaptación de la matriz de covarianza: las elipsis azules muestran cómo cambia la forma de la distribución conforme avanza la generación: 1 - redonda (matriz identidad), 5 - alargada y rotada, 15 - ajustada con precisión a la geometría local.
Componentes del algoritmo: los puntos rojos representan la población total (λ descendientes), los puntos verdes son las mejores soluciones (μ padres), el punto negro m es la media de la población y los círculos punteados son las isolíneas de la función objetivo.
La fórmula clave se encuentra a continuación, con una nota sobre la modificación de la ley de potencias. La ilustración demuestra claramente cómo el algoritmo adapta su estrategia de búsqueda al entorno de la función optimizada, reduciendo gradualmente el espacio de búsqueda y acercándose al óptimo.
Pasemos ahora a escribir el pseudocódigo del algoritmo.
Inicialización
- Establecemos los parámetros del algoritmo:
- Tamaño de la población (lambda) = 50
- Número de padres (mu) = 25
- Tasa de aprendizaje para la actualización de rango 1 (c1) = 0,01
- Tasa de aprendizaje para la actualización de rango-μ (cμ) = 0,8
- Amortiguación del tamaño del paso = 0,6
- Tamaño del paso inicial (sigma) = 0,3
- Calcular los pesos de recombinación:
- Para cada uno de los μ padres, calculamos el peso como log(μ + 0.5) - log(i + 1).
- Normalizamos los pesos de modo que su suma sea igual a 1.
- Calculamos la masa efectiva de elección μ_eff
- Inicializamos los parámetros de la estrategia:
- Calculamos las tasas de aprendizaje para las trayectorias evolutivas (cs, cc)
- Determinamos el valor esperado de la distribución normal estándar (chiN).
- Determinamos el intervalo para la descomposición en valores propios.
- Inicializamos las estructuras de datos:
- Matriz de covarianza C = matriz de identidad
- Matriz de autovectores B = matriz de identidad
- Vector de valores propios D = vector unitario
- Trayectorias de evolución pc y ps = vectores cero
- Media poblacional m = punto aleatorio en el espacio de búsqueda
- Creamos una población inicial:
- Generamos puntos aleatorios lambda en el espacio de búsqueda
El ciclo básico de la evolución
Repetimos el procedimiento hasta que se cumpla el criterio de parada:
Etapa 1: Generación de descendientes
- Para cada uno de los descendientes de lambda:
- Generamos un vector aleatorio z a partir de la distribución normal estándar.
- Aplicamos la transformación: y = B × D × z
- Creamos un descendiente: x_k = m + σ × y
- Aplicamos las restricciones de espacio de búsqueda a x_k
Etapa 2: Evaluación y actualización
- Evaluación de la aptitud de todos los descendientes.
- Clasificación de la población:
- Clasificamos los descendientes en orden descendente de aptitud
- Actualizamos la mejor solución encontrada si es necesario.
- Actualización de la media poblacional:
- Guardamos la media antigua m_old
- Calculamos la nueva media como la suma ponderada de los μ mejores descendientes: m_new = Σ(w_i × x_i) para i desde 1 hasta μ
- Actualización de las rutas de evolución:
- Calculamos el desplazamiento medio: y_w = (m_new - m_old) / σ
- Calculamos C^(-1/2) × y_w mediante descomposición en valores propios.
- Actualizamos la ruta de evolución para el tamaño del paso: ps = (1 - cs) × ps + √(cs × (2 - cs) × μ_eff) × C^(-1/2) × y_w
- Comprobamos la condición de estancamiento (función de Heaviside).
- Actualizamos la trayectoria de evolución para la matriz de covarianza: pc = (1 - cc) × pc + indicador_de_estancamiento × √(cc × (2 - cc) × μ_eff) × y_w
- Actualización de la matriz de covarianza:
- Preparamos la actualización de rango 1 a partir del producto exterior de pc
- Preparamos una actualización de rango-μ a partir de la suma ponderada de las desviaciones de los mejores descendientes.
- Actualizamos C: C = (1 - c1 - cμ) × C + c1 × (pc × pc^T) + cμ × Σ(w_i × y_i × y_i^T)
- Aseguramos la simetría de la matriz
- Actualización del tamaño del paso:
- Calculamos el coeficiente de adaptación según la longitud de la ruta ps.
- Actualizamos de σ: σ = σ × exp((cs/damps) × (||ps||/chiN - 1))
- Limitamos σ dentro de límites razonables para la estabilidad numérica.
- Descomposición en valores propios (si es necesario):
- Si han transcurrido suficientes iteraciones desde la última descomposición:
- Realizamos la expansión de Jacobi de C = B × D² × B^T
- Clasificamos los valores propios y los vectores propios en orden descendente.
- Aseguramos la positividad definida de la matriz
- Si han transcurrido suficientes iteraciones desde la última descomposición:
- Comprobación y corrección de la matriz de covarianza:
- Comprobamos periódicamente la positividad definida de C.
- Si es necesario, añadimos un pequeño valor a la diagonal.
- Aseguramos la simetría de la matriz
Ahora vamos a escribir una clase " C_AO_CMAES", derivada de la clase " C_AO", que implementará el algoritmo de optimización CMA-ES.
Campos públicos, método:
- SetParams () - establece los valores de los parámetros CMA-ES del array "params".
- Init () - inicializa el algoritmo. Admite valores mínimos y máximos, así como el tamaño del paso para cada variable y el número de épocas.
- Moving () - implementa el ciclo principal del algoritmo, que es el responsable de generar nuevas soluciones.
- Revision () - evalúa nuevas soluciones y actualiza el estado del algoritmo.
Campos privados:
- Parámetros CMA-ES - variables que contiene parámetros específicos del algoritmo CMA-ES: número de padres (mu), tamaño del paso (sigma), tasas de aprendizaje (learningRateC1, learningRateCMu), factor de amortiguación (stepSizeDamping) y otros necesarios para que el algoritmo funcione.
- Estructuras de datos CMA-ES - arrays utilizados para almacenar los pesos de recombinación (weights), la matriz de covarianza (covMatrix), los vectores propios (B), los valores propios (D), las rutas de evolución (pc, ps), así como otros datos auxiliares. Las variables "mu_eff", "counteval", "eigeneval", "chiN", "eigenInterval" se utilizan para controlar y gestionar el progreso del algoritmo.
- Variables de optimización del rendimiento: están diseñadas para acelerar los cálculos durante la ejecución del algoritmo, como las tasas de aprendizaje y la amortiguación.
- Arrays auxiliares: se utilizan como almacenamiento temporal para vectores y matrices al realizar diversas operaciones.
- Métodos auxiliares: realizan pasos individuales del algoritmo CMA-ES, como la inicialización de la distribución, la actualización de la distribución, el cálculo de la factorización de valores propios, la clasificación de la población, el cálculo de los pesos, la actualización de la media, la comprobación de la positividad definida y la garantía de la positividad definida.
En general, la clase " C_AO_CMAES" encapsula la lógica del algoritmo CMA-ES, ofreciendo una interfaz para la inicialización, el ajuste de parámetros y la realización de la optimización. Contiene los parámetros, los arrays de datos y los métodos necesarios para implementar el algoritmo. La separación de campos públicos y privados ofrece control de acceso a los componentes internos de la clase.
//———————————————————————————————————————————————————————————————————— class C_AO_CMAES : public C_AO { public: //---------------------------------------------------------- ~C_AO_CMAES () { } C_AO_CMAES () { ao_name = "CMAES"; ao_desc = "Covariance Matrix Adaptation Evolution Strategy"; ao_link = "https://www.mql5.com/en/articles/18227"; // Default parameters popSize = 50; // Default population size (lambda) mu = 25; // Number of parents (half of the population) learningRateC1 = 0.01; // Learning rate for rank-1 update learningRateCMu = 0.8; // Learning rate for rank-μ update stepSizeDamping = 0.6; // Damping for step size // Create and initialize the parameters array ArrayResize (params, 5); params [0].name = "popSize"; params [0].val = popSize; params [1].name = "mu"; params [1].val = mu; params [2].name = "learningRateC1"; params [2].val = learningRateC1; params [3].name = "learningRateCMu"; params [3].val = learningRateCMu; params [4].name = "stepSizeDamping"; params [4].val = stepSizeDamping; } void SetParams () { popSize = (int)params [0].val; mu = (int)params [1].val; learningRateC1 = params [2].val; learningRateCMu = params [3].val; stepSizeDamping = params [4].val; } bool Init (const double &rangeMinP [], // minimum values const double &rangeMaxP [], // maximum values const double &rangeStepP [], // step size const int epochsP = 0); void Moving (); void Revision (); //------------------------------------------------------------------ private: //--------------------------------------------------------- // CMA-ES specific parameters int mu; // Number of parents (selected points) double sigma; // Step size double learningRateC1; // Learning rate for rank-1 update double learningRateCMu; // Learning rate for rank-μ update double stepSizeDamping; // Damping factor for step size update // CMA-ES specific data structures double weights []; // Recombination weights double covMatrix []; // Covariance matrix (stored as a one-dimensional array) double B []; // C eigenvectors double D []; // C eigenvalues (square roots) double pc []; // Evolution path for C double ps []; // Evolution path for step-size control double mu_eff; // Effective selection mass int counteval; // Counter of function evaluations since the last decomposition int eigeneval; // Generation counter when decomposition was performed double chiN; // Expected norm N(0,I) int eigenInterval; // Interval for eigenvalue decomposition // Variables for performance optimization double cs; // Learning rate for the sigma path double cc; // Learning rate for rank-1 path double damps; // Damping for sigma double hsig_threshold; // Threshold for the Heaviside function // Auxiliary arrays double y_vec []; // Mutation vector double arindex []; // Array of indices for sorting double arfitness []; // Fitness array for sorting double temp_vec []; // Temporary vector for matrix operations double invsqrtC_times_yw []; // Temporary storage for C^(-1/2) * y_w // Caching variables for Box-Muller double cached_normal; bool has_cached; // Auxiliary methods void InitDistribution (); void UpdateDistribution (); void ComputeEigendecomposition (); double GetChiN (); void SortPopulation (); void ComputeWeights (); void UpdateMean (); bool CheckPositiveDefinite (); void EnforcePositiveDefinite (); }; //————————————————————————————————————————————————————————————————————
El método "Init" de la clase "C_AO_CMAES" inicializa el algoritmo CMA-ES. Recibe como entrada los valores mínimo y máximo de los parámetros, el tamaño del paso y el número de épocas. En primer lugar, llamamos al método "StandardInit" para realizar la inicialización estándar. El indicador "has_cached" se inicializa en "false" para indicar que no hay ningún valor guardado en la caché para la generación de números aleatorios de Box-Muller. A continuación, el tamaño del paso "sigma" se inicializa con un valor inicial de 0,3 (30% del rango de búsqueda). El método ComputeWeights se llama para calcular los pesos de recombinación, mientras que el método GetChiN se llama para calcular la norma esperada N(0, I).
Los parámetros "cs", "cc", "damps" y "hsig_threshold" se calculan según "mu_eff" (masa efectiva) y el número de coordenadas (coords). Estos parámetros se usan para adaptar la estrategia CMA-ES durante el funcionamiento. El parámetro "eigenInterval" se calcula y se establece para determinar con qué frecuencia se realiza la descomposición de valores propios. Esta configuración está diseñada para optimizar el rendimiento.
Se asigna memoria para los arrays específicos de CMA-ES, como la matriz de covarianza "covMatrix", los vectores propios "B", los valores propios "D", las rutas de evolución "pc" y "ps", los vectores temporales "y_vec", "arindex", "arfitness", "temp_vec" e "invsqrtC_times_yw". Los arrays se recrean únicamente si es necesario cambiar la dimensión del problema.
Los arrays de trayectoria de evolución "pc" y "ps" se inicializan a cero, la matriz de covarianza "covMatrix" y la matriz de vectores propios "B" se inicializan a la matriz identidad, y el array de valores propios "D" se inicializa a unidades. El método "InitDistribution" se llama para inicializar la distribución inicial. A continuación, se reinician los contadores de evaluación de las funciones "counteval" y "eigeneval". El indicador "revision" se establece en "true" para indicar que los valores de aptitud deben recalcularse. Se retornará "true" si la inicialización es exitosa.
Finalmente, el método " Init" realiza toda la inicialización necesaria para que funcione el algoritmo CMA-ES, incluyendo la configuración de parámetros, la asignación de memoria para los arrays y la inicialización de los valores iniciales.
//———————————————————————————————————————————————————————————————————— bool C_AO_CMAES::Init (const double &rangeMinP [], // minimum values const double &rangeMaxP [], // maximum values const double &rangeStepP [], // step size const int epochsP = 0) // number of epochs { if (!StandardInit (rangeMinP, rangeMaxP, rangeStepP)) return false; //------------------------------------------------------------------ // Initialize Box-Muller caching has_cached = false; // Initialize CMA-ES specific variables sigma = 0.3; // Initial step size (30% of search range) // Calculate the effective mass of the variance selection ComputeWeights (); // Expected norm N(0,I) chiN = GetChiN (); // Calculate and save strategy parameters cs = (mu_eff + 2.0) / (coords + mu_eff + 5.0); cc = (4.0 + mu_eff / coords) / (coords + 4.0 + 2.0 * mu_eff / coords); damps = 1.0 + 2.0 * MathMax (0.0, MathSqrt ((mu_eff - 1.0) / (coords + 1.0)) - 1.0) + cs; hsig_threshold = 1.4 + 2.0 / (coords + 1.0); // Set the eigenvalue decomposition interval - tuning for performance eigenInterval = (int)(coords / (10.0 * MathSqrt (learningRateC1 + learningRateCMu))); eigenInterval = MathMax (1, eigenInterval); // Allocate arrays only once ArrayResize (covMatrix, coords * coords); ArrayResize (B, coords * coords); ArrayResize (D, coords); ArrayResize (pc, coords); ArrayResize (ps, coords); ArrayResize (y_vec, coords); ArrayResize (arindex, popSize); ArrayResize (arfitness, popSize); ArrayResize (temp_vec, coords); ArrayResize (invsqrtC_times_yw, coords); // Initialize evolution paths with zeros ArrayInitialize (pc, 0); ArrayInitialize (ps, 0); // Fast initialization of the identity covariance matrix and decomposition ArrayInitialize (covMatrix, 0.0); ArrayInitialize (B, 0.0); for (int i = 0; i < coords; i++) { covMatrix [i * coords + i] = 1.0; B [i * coords + i] = 1.0; D [i] = 1.0; } // Initialize initial distribution InitDistribution (); // Reset calculation counters counteval = 0; eigeneval = 0; // Forced fitness recalculation revision = true; return true; } //————————————————————————————————————————————————————————————————————
El método "Moving" de la clase "C_AO_CMAES" es responsable de generar nuevos descendientes (soluciones) en la población. Transformación y = B * D * z: dentro del ciclo sobre la población hay un ciclo anidado que itera sobre todas las "coordenadas" de cada individuo. En cada coordenada "i" se calcula el valor "y_vec[i]". Esto se realiza usando la multiplicación de matrices. De hecho, supone la multiplicación de un vector de números aleatorios "z" por una matriz "B * D", donde "B" son los vectores propios de la matriz de covarianza, "D" son los valores propios de la matriz de covarianza y "z" son números aleatorios generados mediante la función "u.PowerDistribution".
Generación de descendientes: para cada coordenada "i" de cada descendiente "k", se calcula el valor (a[k].c[i]). Esto se hace sumando a la media actual "cB[i]" (el centro de la distribución actual) el vector escalado "y_vec[i]". El escalado se logra multiplicando "y_vec[i]" por el tamaño del paso "sigma". Por consiguiente, a[k].c[i] = cB[i] + sigma * y_vec[i].
La función "u.SeInDiSp" garantiza que los valores de coordenadas de los descendientes permanezcan dentro de los rangos válidos definidos por "rangeMin[i]", "rangeMax[i]" y "rangeStep[i]". Esta puede recortar valores, hacerlos rebotar en los límites y cuantificarlos al valor de paso válido más cercano. Una vez generados todos los descendientes, el indicador "revision" se establece en "true" y se deben recalcular los valores de la función de aptitud (calidad) para los descendientes generados.
Como resultado, el método genera una nueva población de soluciones basada en la distribución actual (media " cB" y la matriz de covarianza representada por " B". y " D") y el tamaño del paso " sigma". Este usa números aleatorios para crear pequeñas variaciones alrededor de la media actual y así explorar el espacio de búsqueda. La función " SeInDiSp" garantiza que las nuevas decisiones se mantengan dentro de límites aceptables y la bandera de "revision" garantiza que se evaluará la calidad de estas decisiones.
//———————————————————————————————————————————————————————————————————— void C_AO_CMAES::Moving () { // Generate new lambda descendants for (int k = 0; k < popSize; k++) { // Apply the transformation y = B*D*z for (int i = 0; i < coords; i++) { y_vec [i] = 0.0; for (int j = 0; j < coords; j++) { y_vec [i] += B [i * coords + j] * D [j] * u.PowerDistribution (0.0, -8, 8, 20); } // Generate a descendant: x_k = m + σ * y a [k].c [i] = cB [i] + sigma * y_vec [i]; a [k].c [i] = u.SeInDiSp (a [k].c [i], rangeMin [i], rangeMax [i], rangeStep [i]); } } // Fitness recalculation mark revision = true; } //————————————————————————————————————————————————————————————————————
El método "revision" de la clase " C_AO_CMAES" es responsable de actualizar los parámetros del algoritmo después de evaluar la función de aptitud para una nueva población. En primer lugar, el método comprueba el valor del indicador "revision". A continuación, se llama al método "SortPopulation()", que clasifica la población actual de manera que los mejores individuos se coloquen al principio. El método "UpdateDistribution()" actualiza los parámetros de la distribución CMA-ES (la media "cB" del vector de estrategia, la matriz de covarianza, el tamaño del paso "sigma" y las rutas evolutivas) basándose en información sobre los mejores individuos de la población. La actualización de la distribución supone un paso clave de CMA-ES, ya que permite que el algoritmo adapte su estrategia de búsqueda en función de los resultados obtenidos.
Como resultado, el método " Revision" representa el ciclo principal de adaptación de CMA-ES. Este clasifica la población, actualiza los parámetros de distribución en función de los mejores individuos e incrementa el contador de cálculos.
//———————————————————————————————————————————————————————————————————— void C_AO_CMAES::Revision () { if (!revision) return; revision = false; // Sort the population by fitness SortPopulation (); // Update distribution parameters based on selected individuals UpdateDistribution (); // Update the calculation counter counteval++; } //————————————————————————————————————————————————————————————————————
El método "InitDistribution" de la clase "C_AO_CMAES" se encarga de inicializar la distribución inicial de la búsqueda de soluciones. El método determina un valor medio inicial (cB) para cada coordenada en el espacio de búsqueda. Para cada coordenada "i", el valor de "cB[i]" se establece aleatoriamente en el rango entre "rangeMin[i]" y "rangeMax[i]" usando la función "u.RNDfromCI()", con el valor medio inicial colocado en el centro del rango válido para cada variable. El ciclo exterior itera sobre toda la población (popSize individuos). Dentro de este ciclo hay un ciclo interno que itera sobre todas las coordenadas de cada individuo. Para cada descendiente, se generan sus coordenadas. Primero se eligen aleatoriamente del rango dado "RNDfromCI" y luego se redondean al paso de muestreo usando "SeInDiSp".
La función " u.RNDfromCI" genera un número aleatorio con distribución uniforme en un intervalo dado. La función " u.SeInDiSp" garantiza que las coordenadas generadas para cada individuo estén dentro de los límites aceptables definidos rangoMin[i], rangoMax[i] y rangeStep[i].
El método inicializa la población inicial con soluciones aleatorias distribuidas uniformemente en el espacio de búsqueda. El centro de distribución (cB) también se selecciona aleatoriamente. Esto ofrece un punto de partida para las iteraciones posteriores de CMA-ES, en las que el algoritmo adaptará la estrategia de búsqueda para encontrar la solución óptima.
//+------------------------------------------------------------------+ //| Initialize search distribution | //+------------------------------------------------------------------+ void C_AO_CMAES::InitDistribution () { // Set the initial mean to the center of the search space for (int i = 0; i < coords; i++) { cB [i] = u.RNDfromCI (rangeMin [i], rangeMax [i]); } for (int k = 0; k < popSize; k++) { for (int i = 0; i < coords; i++) { // Generate a uniformly distributed point a [k].c [i] = u.RNDfromCI (rangeMin [i], rangeMax [i]); a [k].c [i] = u.SeInDiSp (a [k].c [i], rangeMin [i], rangeMax [i], rangeStep [i]); } } } //————————————————————————————————————————————————————————————————————
El método "GetChiN" calcula la esperanza matemática de la norma de la distribución normal estándar N (0,I) para un número dado de coordenadas (coords). El resultado se aproxima según la dimensionalidad del espacio de búsqueda.
//+------------------------------------------------------------------+ //| Calculate the expected norm N(0,I) | //+------------------------------------------------------------------+ double C_AO_CMAES::GetChiN () { double n = (double)coords; return MathSqrt (n) * (1.0 - 1.0 / (4.0 * n) + 1.0 / (21.0 * n * n)); } //————————————————————————————————————————————————————————————————————
Método "SortPopulation" está diseñado para clasificar una población de soluciones según los valores de su función objetivo (aptitud) con el fin de determinar la mejor solución. En primer lugar, el método copia los valores de aptitud de cada individuo de la población en un array auxiliar llamado "arfitness". También se crea un array llamado "arindex", que inicialmente contiene los índices de los individuos. Esto es necesario para que, después de clasificar según la aptitud, sea posible restablecer la correspondencia entre la aptitud clasificada y el individuo original en la población.
Se utiliza el algoritmo de clasificación por inserción (insertion sort). Este itera a través del array "arfitness" y, para cada elemento, lo inserta en el lugar correcto de la parte clasificada del array, desplazando los elementos más grandes hacia la derecha. Es importante considerar que la clasificación se realiza en orden descendente (arfitness[j] < tempFitness en el ciclo while), el objetivo es maximizar la función objetivo (fitness). Junto con "arfitness", "arindex" también se clasifica de forma sincrónica para rastrear la posición de los índices originales de los individuos.
Tras la clasificación, el método comprueba si la aptitud del mejor individuo es mejor que la de la mejor solución actual y, en caso afirmativo, se actualiza "fB" y la mejor solución se sustituye por las coordenadas del mejor individuo de la población. De este modo, el método no solo clasifica la población según su aptitud, sino que también realiza un seguimiento de la mejor solución encontrada hasta el momento.
//+------------------------------------------------------------------+ //| Sort the population by fitness | //+------------------------------------------------------------------+ void C_AO_CMAES::SortPopulation () { // Copy fitness values and indices for (int i = 0; i < popSize; i++) { arindex [i] = i; arfitness [i] = a [i].f; } for (int i = 1; i < popSize; i++) { double tempFitness = arfitness [i]; double tempIndex = arindex [i]; int j = i - 1; // Sort in descending order (for maximization) while (j >= 0 && arfitness [j] < tempFitness) { arfitness [j + 1] = arfitness [j]; arindex [j + 1] = arindex [j]; j--; } arfitness [j + 1] = tempFitness; arindex [j + 1] = tempIndex; } // Update the best solution if necessary if (arfitness [0] > fB) { fB = arfitness [0]; int bestIdx = (int)arindex [0]; ArrayCopy (cB, a [bestIdx].c, 0, 0, coords); } } //————————————————————————————————————————————————————————————————————
El método "UpdateMean" actualiza el valor medio de la población (cB) mediante la recombinación ponderada de los mejores individuos. Para cada coordenada, se calcula una nueva media como la suma de los valores de coordenadas ponderadas de los mejores individuos, donde los pesos se especifican en el array " weights" .
//+------------------------------------------------------------------+ //| Update the mean using weighted recombination | //+------------------------------------------------------------------+ void C_AO_CMAES::UpdateMean () { // Weighted recombination: m^(g+1) = Σ w_i * x_{i:λ}^(g+1) for (int j = 0; j < coords; j++) { double meanSum = 0.0; for (int i = 0; i < mu; i++) { int idx = (int)arindex [i]; meanSum += weights [i] * a [idx].c [j]; } cB [j] = meanSum; } } //————————————————————————————————————————————————————————————————————
El método ComputeWeights calcula los pesos para la recombinación ponderada, asigna un array "weights" de tamaño "mu" y calcula "log(mu + 0.5)" para usarlo en la fórmula de pesos. Para cada "i" desde "0" hasta "mu-1", asigna un peso como la diferencia entre log(mu + 0.5) y log(i+1); sumando estos pesos, el método divide todos los pesos por la suma para que la suma sea igual a 1, y calcula la suma de los cuadrados de los pesos, calcula "mu_eff", la eficiencia del uso de los pesos, que es importante para ajustar la velocidad y la varianza en el algoritmo.
Este método ofrece pesos óptimos para la recombinación, considerando la contribución de las mejores soluciones.
//+------------------------------------------------------------------+ //| Calculate weighted recombination weights | //+------------------------------------------------------------------+ void C_AO_CMAES::ComputeWeights () { // Allocate the array of weights ArrayResize (weights, mu); // Pre-calculate log(mu + 0.5) double log_mu_plus_half = MathLog (mu + 0.5); // Calculate positive weights double sum = 0.0; for (int i = 0; i < mu; i++) { weights [i] = log_mu_plus_half - MathLog (i + 1); sum += weights [i]; } // Normalize weights double sum_weights = 0.0; double sum_squares = 0.0; for (int i = 0; i < mu; i++) { weights [i] /= sum; sum_weights += weights [i]; sum_squares += weights [i] * weights [i]; } // Calculate the effective mass of the variance selection mu_eff = sum_weights * sum_weights / sum_squares; } //————————————————————————————————————————————————————————————————————
El método "UpdateDistribution" actualiza los parámetros de distribución, a saber, la matriz de covarianza (covMatrix) y el tamaño del paso (sigma) en el algoritmo CMA-ES. Si han transcurrido suficientes generaciones (counteval - eigeneval > eigenInterval), entonces se calculará la descomposición de valores propios y la media actual (cB) se almacenará en el array temporal "oldMean". El método llama a "UpdateMean" para calcular una nueva media, calcula la diferencia entre la media nueva y la antigua dividida por "sigma", y también realiza la multiplicación de "y_w" por C^(-1/2), dividida en varios pasos usando la descomposición "B" y "D", y almacena el resultado en "invsqrtC_times_yw".
El método actualiza la ruta de evolución para el tamaño de paso "ps" usando "invsqrtC_times_yw", determina si el progreso se considera "estancado" en función de la norma "ps" y la longitud esperada, y luego actualiza la ruta de evolución para la matriz de covarianza "pc" usando "y_w" y el indicador de progreso "hsig". El método primero calcula una matriz de actualización de rango 1, c1a, luego una matriz de actualización de rango μ, cmu, basada en la diferencia entre los individuos y la media anterior dividida por sigma usando los pesos, actualiza la matriz de covarianza covMatrix usando las actualizaciones de rango 1 y rango μ, y ajusta c1 cuando el progreso se estanca. La función "EnforcePositiveDefinite" se llama periódicamente para garantizar que la matriz de covarianza sea definida de forma positiva. El tamaño del paso "sigma" se actualiza según la norma "ps" y está restringido entre "min_sigma" y "max_sigma" para garantizar la estabilidad numérica.
En general, " UpdateDistribution" ajusta los parámetros de distribución (matriz de covarianza y tamaño del paso) según la historia de búsqueda (trayectoria evolutiva) y los datos de la población para adaptarse al entorno de la función objetivo.
//+------------------------------------------------------------------+ //| Update distribution parameters | //+------------------------------------------------------------------+ void C_AO_CMAES::UpdateDistribution () { // Check the necessity of eigenvalue decomposition if (counteval - eigeneval > eigenInterval) { ComputeEigendecomposition (); eigeneval = counteval; } // Preserve the old average double oldMean []; ArrayResize (oldMean, coords); ArrayCopy (oldMean, cB, 0, 0, coords); // Update average UpdateMean (); // Calculate the mean shift double y_w []; ArrayResize (y_w, coords); for (int j = 0; j < coords; j++) { y_w [j] = (cB [j] - oldMean [j]) / sigma; } // Calculate C^(-1/2) * y_w // Step 1: B^T * y_w ArrayInitialize (temp_vec, 0.0); for (int i = 0; i < coords; i++) { for (int j = 0; j < coords; j++) { temp_vec [i] += B [j * coords + i] * y_w [j]; } } // Step 2: D^(-1) * (B^T * y_w) for (int i = 0; i < coords; i++) { temp_vec [i] /= D [i]; } // Step 3: B * D^(-1) * B^T * y_w ArrayInitialize (invsqrtC_times_yw, 0.0); for (int i = 0; i < coords; i++) { for (int j = 0; j < coords; j++) { invsqrtC_times_yw [i] += B [i * coords + j] * temp_vec [j]; } } // Update the evolution path for sigma double norm_ps_squared = 0.0; for (int i = 0; i < coords; i++) { ps [i] = (1.0 - cs) * ps [i] + MathSqrt (cs * (2.0 - cs) * mu_eff) * invsqrtC_times_yw [i]; norm_ps_squared += ps [i] * ps [i]; } // Heaviside function double norm_ps = MathSqrt (norm_ps_squared); double expected_length = MathSqrt (1.0 - MathPow (1.0 - cs, 2.0 * counteval)) * chiN; bool hsig = norm_ps / expected_length < hsig_threshold; // Update the evolution path for C double delta_hsig = hsig ? 1.0 : 0.0; for (int i = 0; i < coords; i++) { pc [i] = (1.0 - cc) * pc [i] + delta_hsig * MathSqrt (cc * (2.0 - cc) * mu_eff) * y_w [i]; } // Prepare rank-1 update double c1a []; ArrayResize (c1a, coords * coords); for (int i = 0; i < coords; i++) { for (int j = 0; j <= i; j++) { c1a [i * coords + j] = c1a [j * coords + i] = pc [i] * pc [j]; } } // Prepare rank-μ update double cmu []; ArrayResize (cmu, coords * coords); ArrayInitialize (cmu, 0.0); for (int k = 0; k < mu; k++) { int idx = (int)arindex [k]; // Calculate y_i = (x_i - m_old) / sigma for (int i = 0; i < coords; i++) { temp_vec [i] = (a [idx].c [i] - oldMean [i]) / sigma; } // Add the weighted outer product for (int i = 0; i < coords; i++) { for (int j = 0; j <= i; j++) { double update = weights [k] * temp_vec [i] * temp_vec [j]; cmu [i * coords + j] += update; if (i != j) cmu [j * coords + i] += update; } } } // Update C the covariance matrix double c1 = learningRateC1; double cmu_rate = learningRateCMu; // Adjust c1 if hsig is 'false' (progress stalled) if (!hsig) { c1 *= (1.0 - (1.0 - delta_hsig) * cc * (2.0 - cc)); } double one_minus_c1_cmu = 1.0 - c1 - cmu_rate; // Update C with rank-1 and rank-μ updates for (int i = 0; i < coords; i++) { for (int j = 0; j <= i; j++) { covMatrix [i * coords + j] = one_minus_c1_cmu * covMatrix [i * coords + j] + c1 * c1a [i * coords + j] + cmu_rate * cmu [i * coords + j]; // Maintain symmetry if (i != j) { covMatrix [j * coords + i] = covMatrix [i * coords + j]; } } } // Ensure positive definiteness if (counteval % (10 * eigenInterval) == 0) { EnforcePositiveDefinite (); } // Update the sigma step size double exponent = (cs / damps) * (norm_ps / chiN - 1.0); sigma *= MathExp (exponent); // Limit sigma for numerical stability double min_sigma = 1e-16; double max_eigenvalue = D [0]; // D sorted in descending order double max_sigma = 1e4 * MathMax (1.0, MathSqrt (max_eigenvalue)); if (sigma < min_sigma) sigma = min_sigma; else if (sigma > max_sigma) sigma = max_sigma; } //————————————————————————————————————————————————————————————————————
El método "ComputeEigendecomposition" realiza la descomposición de la matriz de covarianza simétrica " covMatrix" en valores propios y vectores propios utilizando el método Jacobi mejorado; el método crea una copia de la "covMatrix" para no modificar los datos originales. Inicialmente, "B" se define como la matriz identidad que representa la base inicial. Realiza el número máximo de iteraciones o hasta que los elementos distintos de cero en la diagonal principal sean menores que "tolerance".
Búsqueda del elemento máximo fuera de la diagonal: encuentra el elemento con el valor absoluto máximo fuera de la diagonal para seleccionarlo para la rotación. A continuación, si el valor máximo es inferior a la "tolerance", el ciclo se interrumpe. El método calcula "phi" para realizar la rotación jacobiana, modifica los elementos de "C_copy" aplicando la rotación y actualiza las columnas de la matriz "B" (vectores propios) de acuerdo con la rotación. Tras varias iteraciones, el método halla la raíz cuadrada de los elementos diagonales (garantizando una positividad mínima de 1e-14) y clasifica los valores propios y los vectores propios correspondientes en orden descendente para su posterior uso. De este modo, el método permite calcular con precisión los valores propios y los vectores propios de la matriz de covarianza.
//+------------------------------------------------------------------+ //| Calculate the eigenvalue decomposition using the Jacobi method | //+------------------------------------------------------------------+ void C_AO_CMAES::ComputeEigendecomposition () { // Create a copy of the covariance matrix for decomposition double C_copy []; ArrayResize (C_copy, coords * coords); ArrayCopy (C_copy, covMatrix); // Initialize B as the identity matrix for (int i = 0; i < coords; i++) { for (int j = 0; j < coords; j++) { B [i * coords + j] = (i == j) ? 1.0 : 0.0; } } // Improved Jacobi eigenvalue decomposition int max_iterations = 10; //50 * coords; double tolerance = 0.01; //1e-14 * coords * coords; for (int iter = 0; iter < max_iterations; iter++) { // Find the largest off-diagonal element double max_val = 0.0; int p = 0, q = 1; for (int i = 0; i < coords - 1; i++) { for (int j = i + 1; j < coords; j++) { double val = MathAbs (C_copy [i * coords + j]); if (val > max_val) { max_val = val; p = i; q = j; } } } // Check convergence if (max_val < tolerance) break; // Calculate the rotation angle double app = C_copy [p * coords + p]; double aqq = C_copy [q * coords + q]; double apq = C_copy [p * coords + q]; double phi = 0.5 * MathArctan (2.0 * apq / (aqq - app + 1e-14)); double c = MathCos (phi); double s = MathSin (phi); // Update the matrix elements double app_new = c * c * app - 2 * c * s * apq + s * s * aqq; double aqq_new = s * s * app + 2 * c * s * apq + c * c * aqq; C_copy [p * coords + p] = app_new; C_copy [q * coords + q] = aqq_new; C_copy [p * coords + q] = C_copy [q * coords + p] = 0.0; // Update other elements in p and q rows/columns for (int i = 0; i < coords; i++) { if (i != p && i != q) { double aip = C_copy [i * coords + p]; double aiq = C_copy [i * coords + q]; C_copy [i * coords + p] = C_copy [p * coords + i] = c * aip - s * aiq; C_copy [i * coords + q] = C_copy [q * coords + i] = s * aip + c * aiq; } } // Update eigenvectors for (int i = 0; i < coords; i++) { double bip = B [i * coords + p]; double biq = B [i * coords + q]; B [i * coords + p] = c * bip - s * biq; B [i * coords + q] = s * bip + c * biq; } } // Extract eigenvalues and ensure positivity double min_eigenvalue = 1e-14; for (int i = 0; i < coords; i++) { D [i] = MathSqrt (MathMax (min_eigenvalue, C_copy [i * coords + i])); } // Sort eigenvalues and vectors in descending order for (int i = 0; i < coords - 1; i++) { int max_idx = i; for (int j = i + 1; j < coords; j++) { if (D [j] > D [max_idx]) max_idx = j; } if (max_idx != i) { // Exchange eigenvalues double temp = D [i]; D [i] = D [max_idx]; D [max_idx] = temp; // Exchange eigenvectors for (int k = 0; k < coords; k++) { temp = B [k * coords + i]; B [k * coords + i] = B [k * coords + max_idx]; B [k * coords + max_idx] = temp; } } } } //————————————————————————————————————————————————————————————————————
El método " CheckPositiveDefinite" comprueba si la matriz de covarianza es definida positiva. Comprueba rápidamente la positividad de los elementos diagonales y compara el valor propio mínimo (del array " D" , clasificado en orden descendente) con un número positivo pequeño 1e-14 . Si ambas comprobaciones se superan, el método devuelve "true".
//+------------------------------------------------------------------+ //| Check the positive definiteness of the covariance matrix | //+------------------------------------------------------------------+ bool C_AO_CMAES::CheckPositiveDefinite () { // Quick check: all diagonal elements should be positive for (int i = 0; i < coords; i++) { if (covMatrix [i * coords + i] <= 0) return false; } // Check if the minimum eigenvalue is positive double min_eigenvalue = D [coords - 1]; // D is sorted in descending order return min_eigenvalue > 1e-14; } //————————————————————————————————————————————————————————————————————
El método " EnforcePositiveDefinite" garantiza la positividad definida de la matriz de covarianza " covMatrix". Realiza varios pasos: encuentra el elemento mínimo en la diagonal y agrega el valor de "correction" faltante a los elementos diagonales para que el mínimo sea igual a 1e-10, maximiza la simetría de la matriz promediando cada elemento fuera de la diagonal con su transpuesta.
Si la matriz todavía no es definida positiva, se realiza una descomposición en valores propios utilizando "ComputeEigendecomposition", reemplazando todos los valores propios menores que √1e-10 con √1e-10 para garantizar la positividad. El método luego recalcula "covMatrix" usando los vectores propios de "B" y los valores propios corregidos de "D" para obtener una matriz definida positiva corregida. Este enfoque garantiza que la matriz de covarianza sea definida positiva y adecuada para cálculos posteriores en el algoritmo CMA-ES.
//+------------------------------------------------------------------+ //| Ensuring positive definiteness of the covariance matrix | //+------------------------------------------------------------------+ void C_AO_CMAES::EnforcePositiveDefinite () { // Method 1: Adding a small value to the diagonal double min_diag = 1e308; // A very large number for (int i = 0; i < coords; i++) { if (covMatrix [i * coords + i] < min_diag) { min_diag = covMatrix [i * coords + i]; } } if (min_diag < 1e-10) { double correction = 1e-10 - min_diag; for (int i = 0; i < coords; i++) { covMatrix [i * coords + i] += correction; } } // Method 2: Ensuring symmetry for (int i = 0; i < coords; i++) { for (int j = i + 1; j < coords; j++) { double avg = (covMatrix [i * coords + j] + covMatrix [j * coords + i]) * 0.5; covMatrix [i * coords + j] = covMatrix [j * coords + i] = avg; } } // If still not positive definite, perform the factorization and fix if (!CheckPositiveDefinite ()) { ComputeEigendecomposition (); double min_eigenvalue = 1e-10; for (int i = 0; i < coords; i++) { if (D [i] < MathSqrt (min_eigenvalue)) { D [i] = MathSqrt (min_eigenvalue); } } // Reconstruct C = B * D^2 * B^T ArrayInitialize (covMatrix, 0.0); for (int i = 0; i < coords; i++) { for (int j = 0; j <= i; j++) { double sum = 0.0; for (int k = 0; k < coords; k++) { sum += B [i * coords + k] * D [k] * D [k] * B [j * coords + k]; } covMatrix [i * coords + j] = covMatrix [j * coords + i] = sum; } } } } //————————————————————————————————————————————————————————————————————
Resultados de las pruebas
Ahora por fin podemos ver los resultados. Como podemos observar de inmediato, el algoritmo funciona bien y con rapidez con problemas de tamaño mediano y, en ciertas configuraciones, con problemas de baja dimensionalidad, pero los problemas multidimensionales más difíciles requieren una cantidad desproporcionada de tiempo para ejecutarse, por lo que fueron excluidos de los resultados. ¿Con qué está relacionado esto? Tradicionalmente, los cálculos matriciales en MQL5 presentan un grave problema computacional. Al implementar operaciones matriciales manualmente, el algoritmo muestra un rendimiento insatisfactorio en problemas de alta dimensionalidad, lo cual limita significativamente su aplicabilidad práctica. Sin embargo, con la llegada de las clases integradas para trabajar con matrices, la situación cambia drásticamente. Ahora es fundamental usar la implementación nativa de la plataforma para las operaciones matriciales en todas las tareas que requieren una gran capacidad de cálculo.
=============================
5 Hilly's; Func runs: 10000; result: 0,7625797883550075
25 Hilly's; Func runs: 10000; result: 0,7208874560706138
500 Hilly's; Func runs: 10000; result: 0.0
=============================
5 Forest's; Func runs: 10000; result: 0,8205636421348295
25 Forest's; Func runs: 10000; result: 0,7961602627346933
500 Forest's; Func runs: 10000; result: 0.0
=============================
5 Megacity's; Func runs: 10000; result: 0,7584615384615383
25 Megacity's; Func runs: 10000; result: 0,49076923076923074
500 Megacity's; Func runs: 10000; result: 0.0
=============================
All score: 4,34942 (48,33%)
En la visualización, podemos ver que el algoritmo funciona tanto con "grandes saltos" como centrándose en los extremos de la función, examinándolos cuidadosamente.

CMA-ES en la función de prueba Hilly

CMA-ES en la función de prueba Forest

CMA-ES en la función de prueba Megacity
Según los resultados del trabajo, el algoritmo CMA-ES ocupa el puesto 38 en la clasificación de los algoritmos de optimización basados en poblaciones más potentes.
# | AO | Description | Hilly | Hilly Final | Forest | Forest Final | Megacity (discrete) | Megacity Final | Final Result | % de MAX | ||||||
| 10 p (5 F) | 50 p (25 F) | 1000 p (500 F) | 10 p (5 F) | 50 p (25 F) | 1000 p (500 F) | 10 p (5 F) | 50 p (25 F) | 1000 p (500 F) | ||||||||
| 1 | ANS | across neighbourhood search | 0.94948 | 0.84776 | 0.43857 | 2.23581 | 1.00000 | 0.92334 | 0.39988 | 2.32323 | 0.70923 | 0.63477 | 0.23091 | 1.57491 | 6.134 | 68.15 |
| 2 | CLA | code lock algorithm (joo) | 0.95345 | 0.87107 | 0.37590 | 2.20042 | 0.98942 | 0.91709 | 0.31642 | 2.22294 | 0.79692 | 0.69385 | 0.19303 | 1.68380 | 6.107 | 67.86 |
| 3 | AMOm | animal migration ptimization M | 0.90358 | 0.84317 | 0.46284 | 2.20959 | 0.99001 | 0.92436 | 0.46598 | 2.38034 | 0.56769 | 0.59132 | 0.23773 | 1.39675 | 5.987 | 66.52 |
| 4 | (P+O)ES | (P+O) evolution strategies | 0.92256 | 0.88101 | 0.40021 | 2.20379 | 0.97750 | 0.87490 | 0.31945 | 2.17185 | 0.67385 | 0.62985 | 0.18634 | 1.49003 | 5.866 | 65.17 |
| 5 | CTA | comet tail algorithm (joo) | 0.95346 | 0.86319 | 0.27770 | 2.09435 | 0.99794 | 0.85740 | 0.33949 | 2.19484 | 0.88769 | 0.56431 | 0.10512 | 1.55712 | 5.846 | 64.96 |
| 6 | TETA | time evolution travel algorithm (joo) | 0.91362 | 0.82349 | 0.31990 | 2.05701 | 0.97096 | 0.89532 | 0.29324 | 2.15952 | 0.73462 | 0.68569 | 0.16021 | 1.58052 | 5.797 | 64.41 |
| 7 | SDSm | stochastic diffusion search M | 0.93066 | 0.85445 | 0.39476 | 2.17988 | 0.99983 | 0.89244 | 0.19619 | 2.08846 | 0.72333 | 0.61100 | 0.10670 | 1.44103 | 5.709 | 63.44 |
| 8 | BOAm | billiards optimization algorithm M | 0.95757 | 0.82599 | 0.25235 | 2.03590 | 1.00000 | 0.90036 | 0.30502 | 2.20538 | 0.73538 | 0.52523 | 0.09563 | 1.35625 | 5.598 | 62.19 |
| 9 | AAm | archery algorithm M | 0.91744 | 0.70876 | 0.42160 | 2.04780 | 0.92527 | 0.75802 | 0.35328 | 2.03657 | 0.67385 | 0.55200 | 0.23738 | 1.46323 | 5.548 | 61.64 |
| 10 | ESG | evolution of social groups (joo) | 0.99906 | 0.79654 | 0.35056 | 2.14616 | 1.00000 | 0.82863 | 0.13102 | 1.95965 | 0.82333 | 0.55300 | 0.04725 | 1.42358 | 5.529 | 61.44 |
| 11 | SIA | simulated isotropic annealing (joo) | 0.95784 | 0.84264 | 0.41465 | 2.21513 | 0.98239 | 0.79586 | 0.20507 | 1.98332 | 0.68667 | 0.49300 | 0.09053 | 1.27020 | 5.469 | 60.76 |
| 12 | BBO | biogeography based optimization | 0.94912 | 0.69456 | 0.35031 | 1.99399 | 0.93820 | 0.67365 | 0.25682 | 1.86867 | 0.74615 | 0.48277 | 0.17369 | 1.40261 | 5.265 | 58.50 |
| 13 | ACS | artificial cooperative search | 0.75547 | 0.74744 | 0.30407 | 1.80698 | 1.00000 | 0.88861 | 0.22413 | 2.11274 | 0.69077 | 0.48185 | 0.13322 | 1.30583 | 5.226 | 58.06 |
| 14 | DA | dialectical algorithm | 0.86183 | 0.70033 | 0.33724 | 1.89940 | 0.98163 | 0.72772 | 0.28718 | 1.99653 | 0.70308 | 0.45292 | 0.16367 | 1.31967 | 5.216 | 57.95 |
| 15 | BHAm | black hole algorithm M | 0.75236 | 0.76675 | 0.34583 | 1.86493 | 0.93593 | 0.80152 | 0.27177 | 2.00923 | 0.65077 | 0.51646 | 0.15472 | 1.32195 | 5.196 | 57.73 |
| 16 | ASO | anarchy society optimization | 0.84872 | 0.74646 | 0.31465 | 1.90983 | 0.96148 | 0.79150 | 0.23803 | 1.99101 | 0.57077 | 0.54062 | 0.16614 | 1.27752 | 5.178 | 57.54 |
| 17 | RFO | royal flush optimization (joo) | 0.83361 | 0.73742 | 0.34629 | 1.91733 | 0.89424 | 0.73824 | 0.24098 | 1.87346 | 0.63154 | 0.50292 | 0.16421 | 1.29867 | 5.089 | 56.55 |
| 18 | AOSm | atomic orbital search M | 0.80232 | 0.70449 | 0.31021 | 1.81702 | 0.85660 | 0.69451 | 0.21996 | 1.77107 | 0.74615 | 0.52862 | 0.14358 | 1.41835 | 5.006 | 55.63 |
| 19 | TSEA | turtle shell evolution algorithm (joo) | 0.96798 | 0.64480 | 0.29672 | 1.90949 | 0.99449 | 0.61981 | 0.22708 | 1.84139 | 0.69077 | 0.42646 | 0.13598 | 1.25322 | 5.004 | 55.60 |
| 20 | DE | differential evolution | 0.95044 | 0.61674 | 0.30308 | 1.87026 | 0.95317 | 0.78896 | 0.16652 | 1.90865 | 0.78667 | 0.36033 | 0.02953 | 1.17653 | 4.955 | 55.06 |
| 21 | SRA | successful restaurateur algorithm (joo) | 0.96883 | 0.63455 | 0.29217 | 1.89555 | 0.94637 | 0.55506 | 0.19124 | 1.69267 | 0.74923 | 0.44031 | 0.12526 | 1.31480 | 4.903 | 54.48 |
| 22 | CRO | chemical reaction optimization | 0.94629 | 0.66112 | 0.29853 | 1.90593 | 0.87906 | 0.58422 | 0.21146 | 1.67473 | 0.75846 | 0.42646 | 0.12686 | 1.31178 | 4.892 | 54.36 |
| 23 | BIO | blood inheritance optimization (joo) | 0.81568 | 0.65336 | 0.30877 | 1.77781 | 0.89937 | 0.65319 | 0.21760 | 1.77016 | 0.67846 | 0.47631 | 0.13902 | 1.29378 | 4.842 | 53.80 |
| 24 | BSA | bird swarm algorithm | 0.89306 | 0.64900 | 0.26250 | 1.80455 | 0.92420 | 0.71121 | 0.24939 | 1.88479 | 0.69385 | 0.32615 | 0.10012 | 1.12012 | 4.809 | 53.44 |
| 25 | HS | harmony search | 0.86509 | 0.68782 | 0.32527 | 1.87818 | 0.99999 | 0.68002 | 0.09590 | 1.77592 | 0.62000 | 0.42267 | 0.05458 | 1.09725 | 4.751 | 52.79 |
| 26 | SSG | saplings sowing and growing | 0.77839 | 0.64925 | 0.39543 | 1.82308 | 0.85973 | 0.62467 | 0.17429 | 1.65869 | 0.64667 | 0.44133 | 0.10598 | 1.19398 | 4.676 | 51.95 |
| 27 | BCOm | bacterial chemotaxis optimization M | 0.75953 | 0.62268 | 0.31483 | 1.69704 | 0.89378 | 0.61339 | 0.22542 | 1.73259 | 0.65385 | 0.42092 | 0.14435 | 1.21912 | 4.649 | 51.65 |
| 28 | ABO | african buffalo optimization | 0.83337 | 0.62247 | 0.29964 | 1.75548 | 0.92170 | 0.58618 | 0.19723 | 1.70511 | 0.61000 | 0.43154 | 0.13225 | 1.17378 | 4.634 | 51.49 |
| 29 | (PO)ES | (PO) evolution strategies | 0.79025 | 0.62647 | 0.42935 | 1.84606 | 0.87616 | 0.60943 | 0.19591 | 1.68151 | 0.59000 | 0.37933 | 0.11322 | 1.08255 | 4.610 | 51.22 |
| 30 | FBA | Fractal-Based Algorithm | 0.79000 | 0.65134 | 0.28965 | 1.73099 | 0.87158 | 0.56823 | 0.18877 | 1.62858 | 0.61077 | 0.46062 | 0.12398 | 1.19537 | 4.555 | 50.61 |
| 31 | TSm | tabu search M | 0.87795 | 0.61431 | 0.29104 | 1.78330 | 0.92885 | 0.51844 | 0.19054 | 1.63783 | 0.61077 | 0.38215 | 0.12157 | 1.11449 | 4.536 | 50.40 |
| 32 | BSO | brain storm optimization | 0.93736 | 0.57616 | 0.29688 | 1.81041 | 0.93131 | 0.55866 | 0.23537 | 1.72534 | 0.55231 | 0.29077 | 0.11914 | 0.96222 | 4.498 | 49.98 |
| 33 | WOAm | wale optimization algorithm M | 0.84521 | 0.56298 | 0.26263 | 1.67081 | 0.93100 | 0.52278 | 0.16365 | 1.61743 | 0.66308 | 0.41138 | 0.11357 | 1.18803 | 4.476 | 49.74 |
| 34 | AEFA | artificial electric field algorithm | 0.87700 | 0.61753 | 0.25235 | 1.74688 | 0.92729 | 0.72698 | 0.18064 | 1.83490 | 0.66615 | 0.11631 | 0.09508 | 0.87754 | 4.459 | 49.55 |
| 35 | AEO | artificial ecosystem-based optimization algorithm | 0.91380 | 0.46713 | 0.26470 | 1.64563 | 0.90223 | 0.43705 | 0.21400 | 1.55327 | 0.66154 | 0.30800 | 0.28563 | 1.25517 | 4.454 | 49.49 |
| 36 | CAm | camel algorithm M | 0.78684 | 0.56042 | 0.35133 | 1.69859 | 0.82772 | 0.56041 | 0.24336 | 1.63149 | 0.64846 | 0.33092 | 0.13418 | 1.11356 | 4.444 | 49.37 |
| 37 | ACOm | ant colony optimization M | 0.88190 | 0.66127 | 0.30377 | 1.84693 | 0.85873 | 0.58680 | 0.15051 | 1.59604 | 0.59667 | 0.37333 | 0.02472 | 0.99472 | 4.438 | 49.31 |
| 38 | CMAES | covariance_matrix_adaptation_evolution_strategy | 0.76258 | 0.72089 | 0.00000 | 1.48347 | 0.82056 | 0.79616 | 0.00000 | 1.61672 | 0.75846 | 0.49077 | 0.00000 | 1.24923 | 4.349 | 48.33 |
| 39 | BFO-GA | bacterial foraging optimization - ga | 0.89150 | 0.55111 | 0.31529 | 1.75790 | 0.96982 | 0.39612 | 0.06305 | 1.42899 | 0.72667 | 0.27500 | 0.03525 | 1.03692 | 4.224 | 46.93 |
| 40 | SOA | simple optimization algorithm | 0.91520 | 0.46976 | 0.27089 | 1.65585 | 0.89675 | 0.37401 | 0.16984 | 1.44060 | 0.69538 | 0.28031 | 0.10852 | 1.08422 | 4.181 | 46.45 |
| 41 | ABHA | artificial bee hive algorithm | 0.84131 | 0.54227 | 0.26304 | 1.64663 | 0.87858 | 0.47779 | 0.17181 | 1.52818 | 0.50923 | 0.33877 | 0.10397 | 0.95197 | 4.127 | 45.85 |
| 42 | ACMO | atmospheric cloud model optimization | 0.90321 | 0.48546 | 0.30403 | 1.69270 | 0.80268 | 0.37857 | 0.19178 | 1.37303 | 0.62308 | 0.24400 | 0.10795 | 0.97503 | 4.041 | 44.90 |
| 43 | ADAMm | adaptive moment estimation M | 0.88635 | 0.44766 | 0.26613 | 1.60014 | 0.84497 | 0.38493 | 0.16889 | 1.39880 | 0.66154 | 0.27046 | 0.10594 | 1.03794 | 4.037 | 44.85 |
| 44 | CGO | chaos game optimization | 0.57256 | 0.37158 | 0.32018 | 1.26432 | 0.61176 | 0.61931 | 0.62161 | 1.85267 | 0.37538 | 0.21923 | 0.19028 | 0.78490 | 3.902 | 43.35 |
| 45 | CROm | coral reefs optimization M | 0.78512 | 0.46032 | 0.25958 | 1.50502 | 0.86688 | 0.35297 | 0.16267 | 1.38252 | 0.63231 | 0.26738 | 0.10734 | 1.00703 | 3.895 | 43.27 |
| RW | random walk | 0.48754 | 0.32159 | 0.25781 | 1.06694 | 0.37554 | 0.21944 | 0.15877 | 0.75375 | 0.27969 | 0.14917 | 0.09847 | 0.52734 | 2.348 | 26.09 | |
Conclusiones
La capacidad del CMA-ES para adaptar la estrategia de búsqueda a la geometría local de la función objetivo la hace prácticamente indispensable para problemas complejos y mal condicionados con estructura desconocida. Las colas pesadas de la distribución de ley de potencias permiten que el algoritmo realice "saltos largos", lo que resulta fundamental para escapar de los óptimos locales. La complejidad computacional de O(n²) en la memoria y O(n³) en el tiempo limita severamente la aplicabilidad del algoritmo. Para problemas de alta dimensionalidad (n > 100), la intensidad de los recursos se vuelve desproporcionada con respecto a los beneficios obtenidos. El tiempo de ejecución en funciones multidimensionales crece tan rápidamente que hace que el algoritmo resulte prácticamente inaplicable para valores grandes de n.
El CMA-ES funciona con funciones ruidosas, maneja funciones objetivo discontinuas, es eficaz en entornos con múltiples extremos, y el secreto de esta versatilidad reside en la filosofía fundamental del algoritmo: en lugar de hacer suposiciones sobre la estructura de la función, la explora cuidadosamente, adaptando su estrategia de búsqueda según la experiencia acumulada. El CMA-ES transforma el panorama de la computación evolutiva al demostrar que los algoritmos inspirados en la biología pueden tener fundamentos matemáticos rigurosos y que la adaptación de algoritmos no es simplemente un ajuste de parámetros, sino un principio fundamental para aprender la estructura de un problema.
El algoritmo es idóneo para problemas de tamaño mediano, pues supone una herramienta especializada que parece una alternativa valiosa en su nicho. El algoritmo combina elegancia matemática con eficiencia práctica, pero requiere una aplicación cuidadosa debido a sus limitaciones computacionales. En estos casos, los desarrolladores del lenguaje MQL5 han introducido recientemente operaciones rápidas con matrices para calcularlas; en esta implementación del algoritmo, se usaron cálculos secuenciales convencionales.
Figura 2. Clasificación por colores de los algoritmos según las pruebas respectivas

Figura 3. Histograma de los resultados de las pruebas de algoritmos (en una escala de 0 a 100, cuanto más mejor, donde 100 es el máximo resultado teórico posible, el script para calcular la tabla de puntuación está en el archivo)
Ventajas y desventajas del algoritmo CMAES:
Ventajas:
- Buena convergencia en funciones de dimensionalidad media.
Desventajas:
- Dispersión de resultados en funciones de baja dimensionalidad.
- Consumo crítico de recursos en funciones de alta dimensionalidad.
Adjuntamos al artículo un archivo con las versiones actuales de los códigos de los algoritmos. El autor de este artículo no se responsabiliza de la exactitud absoluta de la descripción de los algoritmos canónicos, muchos de ellos han sido modificados para mejorar las capacidades de búsqueda. Las conclusiones y juicios presentados en los artículos se basan en los resultados de los experimentos realizados.
Programas usados en el artículo
| # | Nombre | Tipo | Descripción |
|---|---|---|---|
| 1 | #C_AO.mqh | Archivo de inclusión | Clase padre de algoritmos de optimización basados en la población |
| 2 | #C_AO_enum.mqh | Archivo de inclusión | Enumeración de los algoritmos de optimización basados en la población |
| 3 | TestFunctions.mqh | Archivo de inclusión | Biblioteca de funciones de prueba |
| 4 | TestStandFunctions.mqh | Archivo de inclusión | Biblioteca de funciones del banco de pruebas |
| 5 | Utilities.mqh | Archivo de inclusión | Biblioteca de funciones auxiliares |
| 6 | CalculationTestResults.mqh | Archivo de inclusión | Script para calcular los resultados en una tabla comparativa |
| 7 | Testing AOs.mq5 | Script | Banco de pruebas único para todos los algoritmos de optimización basados en la población |
| 8 | Simple use of population optimization algorithms.mq5 | Script | Ejemplo sencillo de utilización de algoritmos de optimización basados en la población sin visualización |
| 9 | Test_AO_CMAES.mq5 | Script | Banco de pruebas para CMAES |
Traducción del ruso hecha por MetaQuotes Ltd.
Artículo original: https://www.mql5.com/ru/articles/18227
Advertencia: todos los derechos de estos materiales pertenecen a MetaQuotes Ltd. Queda totalmente prohibido el copiado total o parcial.
Este artículo ha sido escrito por un usuario del sitio web y refleja su punto de vista personal. MetaQuotes Ltd. no se responsabiliza de la exactitud de la información ofrecida, ni de las posibles consecuencias del uso de las soluciones, estrategias o recomendaciones descritas.
Aprendizaje automático y Data Science (Parte 41): Detección de patrones en los mercados de divisas y de valores mediante YOLOv8
Características del Wizard MQL5 que debe conocer (Parte 65): Uso de los patrones FrAMA y Force Index
Redes neuronales en el trading: Pipeline de pronóstico inteligente (Time-MoE)
Redes neuronales en el trading: Pipeline inteligente de previsiones (Mezcla dispersa de expertos)
- Aplicaciones de trading gratuitas
- 8 000+ señales para copiar
- Noticias económicas para analizar los mercados financieros
Usted acepta la política del sitio web y las condiciones de uso
