Использование генетического алгоритма для поиска решения уравнения

Программирование
Как-то приходилось писать программу для составление инвестиционного портфеля, анализу и прогнозированию операций по этому портфелю. Саму программу по некоторым причинам опубликовать не могу, а вот рассмотреть один из основных участков — использование генетического алгоритма — можно.
Итак, полностью с описанием генетического алгоритма можно ознакомиться здесь.
Прочитав информацию по ссылке, Вы узнали, что есть такие понятия как популяция, особь (хромосом), ген, приспособленность особи.
Популяция состоит из особей, особи в свою очередь состоят из генов.

Популяция — множество решений некого уравнения. Т.е. допустим, имеем уравнение x*x + y*y + z*z = 10. Набор различных комбинаций переменных этого уравнения (x, y, z), причем далеко не оптимальных и есть популяция. Пример популяции: { (1, 1, 1), (2, 3, 1), (2, 0, 4), (1, 2, 3), (0, 0, 2)}. Эти значения выбраны абсолютно случайным образом. Популяция состоит из 5 особей.
Особь 1: (1, 1, 1)
Особь 2: (2, 3, 1)
Особь 3: (2, 0, 4)
Особь 4: (1, 2, 3)
Особь 5: (0, 0, 2)

Каждая особь состоит из трех генов (неизвестные в уравнение).
Гены для особи 1 — 1, 1, 1, для особи 2 — 2, 3, 1, и так далее.

Приспособленность особи показывает — как конкретная особь далека от решения уравнения.
Перепишем наше уравнение x*x + y*y + z*z = 10 к виду x*x + y*y + z*z — 10 = 0.
Тогда получается, если в это уравнение подставить каждую особь — мы получим некоторое значение, которое будет показывать — как далеко мы находимся от решения. Так, как нас интересует как далеко мы находимся от решения, то приспособленность не имеет знака.
Посчитаем приспособленность для каждой особи:
f(1, 1, 1) = abs(1*1 + 1*1 + 1*1 - 10) = 7
f(2, 3, 1) = abs(2*2 + 3*3 + 1*1 - 10) = 4
f(2, 0, 4) = abs(2*2 + 0*0 + 4*4 - 10) = 10
f(1, 2, 3) = abs(1*1 + 2*2 + 3*3 - 10) = 4
f(0, 0, 2) = abs(0*0 + 0*0 + 2*2 - 10) = 6

Видно, что некоторые из особей (особь 2 и особь 4) находятся ближе к решения, чем остальные.
У нас получилась первая популяция и мы узнали приспособленность каждой особи в этой популяции. Именно от приспособленности будет зависеть следующая популяция.

Теперь рассчитываем следующую популяцию. Каждая особь следующей популяции есть скрещивание двух особей от предыдущего поколения. Причем, для скрещивания с большей вероятностью надо отбирать те особи, приспособленность которых лучше. Случайным образом с учетом вышесказанного выбираем пять пар особей из предыдущей популяции. Например:
особь 2 - особь 3
особь 4 - особь 1
особь 2 - особь 4
особь 1 - особь 2
особь 4 - особь 5

При скрещивании появляются такие операции как кроссовер и мутация.
Кроссовер — выбора части генов из одной особи и остальной части генов из другой.
Берем первую пару, особь 4 — особь 1:
(2, 3, 1)
(2, 0, 4)

и случайным образом разбиваем ее на две части:
(2, 3,| 1)
(2, 0,| 4)

затем склеиваем из этих 4 кусочков одну новую особь:
(2, 3, 4)

т.е. из первой особи взяли гены 2, 3, а из второй особи ген 4.

Таким же образом поступаем с остальным четырьмя парами:
особь 4 - особь 1
(1, 2,| 3)
(1, 1,| 1)
(1, 2, 1)
особь 2 - особь 4
(2, |3, 1)
(1, |2, 3)
(2, 2, 3)
особь 1 - особь 2
(1, |1, 1)
(2, |3, 1)
(1, 3, 1)
особь 4 - особь 5
(1, 2,| 3)
(0, 0,| 2)
(1, 2, 2)


В результате у нас получилась вторая популяция, состоящая из особей:
Особь 1: (2, 3, 4)
Особь 2: (1, 2, 1)
Особь 3: (2, 2, 3)
Особь 4: (1, 3, 1)
Особь 5: (1, 2, 2)


Операция мутации изменяет значения генов в особи на случайную небольшую величину, например, мы будем менять в пределах 10%. Берем первую особь и изменяем значение каждого в соответствии с вышесказанным:
Особь 1: (2 + 5%, 3 - 3%, 4 + 2%) = (2.10, 2.91, 4.08)

Аналогично поступаем с остальными особями:
Особь 2: (1 + 3%, 2 + 2%, 1 - 4%) = (1.03, 2.04, 0.96)
Особь 3: (2 + 3%, 2 - 4%, 3 + 1%) = (2.06, 1.92, 3.03)
Особь 4: (1 - 1%, 3 + 3%, 1 - 2%) = (0.99, 3.09, 0.98)
Особь 5: (1 + 4%, 2 - 2%, 2 - 5%) = (1.04, 1.96, 1.90)

После операций кроссовера и мутации у нас получилось такое поколение:
Особь 1: (2.10, 2.91, 4.08)
Особь 2: (1.03, 2.04, 0.96)
Особь 3: (2.06, 1.92, 3.03)
Особь 4: (0.99, 3.09, 0.98)
Особь 5: (1.04, 1.96, 1.90)

Рассчитываем приспособленность особей в новом поколении:
f(2.10, 2.91, 4.08) = abs(2.10*2.10 + 2.91*2.91 + 4.08*4.08 - 10) = 19.5245
f(1.03, 2.04, 0.96) = abs(1.03*1.03 + 2.04*2.04 + 0.96*0.96 - 10) = 3.8559
f(2.06, 1.92, 3.03) = abs(2.06*2.06 + 1.92*1.92 + 3.03*3.03 - 10) = 7.1109
f(0.99, 3.09, 0.98) = abs(0.99*0.99 + 3.09*3.09 + 0.98*0.98 - 10) = 1.4886
f(1.04, 1.96, 1.90) = abs(1.04*1.04 + 1.96*1.96 + 1.90*1.90 - 10) = 1.4668
Получилось 2 особи, которые еще ближе подобрались к решению уравнения.

Таким же образом формируются последующие популяции. С каждой новой популяцией приспособленность будет все лучше и лучше.
Прекратить поиск можно по нескольким критериям: количество популяция, приспособленность одной из особей попадает в заданную точность, изменения приспособленности не изменяется и так далее, вариантов очень много. Я реализовал в приведенной ниже программе — по количеству популяций, либо пока приспособленность не попадет в заданную точность.

Решение уравнения x*x + y*y + z*z = 10(решение).
[dk@moon ~]$ ./gaTest
Populations: 226
x = { -1.9082 -2.1811 -1.2654 }
f(x) = -0.000018


Решение уравнения x + x*x = 5 (решение):
[dk@moon ~]$ ./gaTest
Populations: 541
x = { -2.7913 }
f(x) = -0.000010


Листинг gaTest.c:
#include <stdio.h>
#include <stdlib.h>
#include "gen.h"

/* Целевая функция */
float f(float *x)
{
    return x[0] + x[0]*x[0] - 5;
}

int main(int argc, char *argv[])
{
    sranddev();
    float x[3];
    int iterations = gaSearch(30, 1, -10, 10, 1000, 0.00001, f, x);
    // Количество итераций
    printf("Populations: %d\n", iterations);
    // Найденое решение
    printf("x = { ");
    for (int i = 0; i < 1; i ++) {
        printf("%0.4f ", x[i]);
    }
    printf("}\n");
    // Значение функции
    printf("f(x) = %f\n", f(x));
    return(0);
}


Листинг gen.h:
/* Особь */
typedef struct {
    float               *gens;                  // Список генов
    float               fitness;                // Приспособленность особи
    int                 length;                 // Размер списка генов
} INDIVIDUAL;


/* Популяция */
typedef struct {
    INDIVIDUAL          **individuals;          // Список особей
    int                 length;                 // Размер популяции
} POPULATION;


/* Список отбора особей */
typedef struct {
    INDIVIDUAL          **individuals;          // Список особей
    int                 length;                 // Размер списка особей
} SELECTION_ARRAY;


/* Поиск решения уравнения с помощью генетического алгоритма
При выходе из функции возвращается число затраченных популяций на поиск решения
*/
int gaSearch(
    int                 populationLength,       // Число особей в популяции
    int                 individualLength,       // Число генов в особи
    float               min,                    // Минимальное и максимальное значение для формирования
    float               max,                    // генов в первой популяции
    int                 maxPopulations,         // Максимальное число популяций
    float               delta,                  // Заданная точность
    float               f(float *),             // Целевая функция
    float               *x                      // Указатель на массив, куда поместится решение уравнения
);


Листинг gen.c:
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include "gen.h"

/* Функция расчета приспособленности особи f(X1 .. Xn) */
float fitness(INDIVIDUAL *individual, float f(float *))
{
    if (f) {
	return fabsf(f(individual->gens));
    }
    else return 0;
}


/* Создание новой популяции */
POPULATION *newPopulation(int populationLength, int individualLength, float min, float max)
{
    POPULATION *population = NULL;
    // Проверяеем входные параметры
    if (populationLength > 0 && individualLength > 0 && min < max){
	// Выделяем память под популяцию
	population = (POPULATION *) malloc(sizeof(POPULATION));
	// Выделяем память под список особей
	population->individuals = (INDIVIDUAL **) malloc(sizeof(INDIVIDUAL *) * populationLength);
	population->length = populationLength;
	// Проходимся в цикле по всем особям популяции
	for(int i = 0; i < population->length; i ++) {
	    // Выделяем память под особь
	    INDIVIDUAL *individual = (INDIVIDUAL *) malloc(sizeof(INDIVIDUAL));
	    individual->gens = (float *) malloc(sizeof(float) * individualLength);
	    individual->length = individualLength;
	    // Устанавливаем значения генов особи случайным образом в соотв. с параметрами min и max
	    for (int n = 0; n < individual->length; n ++) {
		individual->gens[n] = (float) rand() / RAND_MAX * (max - min) + min;
	    }
	    // Приспособленность особи
	    individual->fitness = 0;
	    // Добавляем особь в популяцию
	    population->individuals[i] = individual;
	}
    }
    // Возвращаем популяцию
    return population;
}

/* Удаление популяции */
void deletePopulation(POPULATION *population)
{
    // Проверяем входные параметры
    if (population) {
	// Освобождаем память, занятую популяцией
	for (int i = 0; i < population->length; i ++) {
	    free(population->individuals[i]->gens);
	    free(population->individuals[i]);
	}
	free(population->individuals);
	free(population);
    }
}


/* Разчет приспособленности каждой особи в популяции */
void calcPopulationFitness(POPULATION *population, float f(float *))
{
    // Проверяем входные параметры
    if (population && population->length > 0) {
	// Рассчитывем приспособленность каждой особи в популяции
	for (int i = 0; i < population->length; i ++) {
	    population->individuals[i]->fitness = fitness(population->individuals[i], f);
	}
    }
}


SELECTION_ARRAY *buildSelectionArray(POPULATION *population)
{
    SELECTION_ARRAY *selectionArray = NULL;
    // Проверяем входные параметры
    if (population && population->length > 0) {
	// Ищем максимальную и минимальную приспособленность особи
	float maxFitness = population->individuals[0]->fitness;
	float minFitness = population->individuals[0]->fitness;
	for (int i = 1; i < population->length; i ++) {
	    if (population->individuals[i]->fitness > maxFitness)
		maxFitness = population->individuals[i]->fitness;
	    if (population->individuals[i]->fitness < minFitness)
		minFitness = population->individuals[i]->fitness;
	}
	// Определяем веса особей для участия в отборе
	int *weight = (int *) malloc(sizeof(int) * population->length);
	int length = 0;
	for (int i = 0; i < population->length; i ++) {
	    // Вес особи
	    weight[i] = (int) exp((1.0 - (population->individuals[i]->fitness - minFitness) / (maxFitness - minFitness)) * 3.0);
	    length += weight[i];
	}
	// Создаем список особей и заполняем в соответствии с весом
	if (length > 0) {
	    selectionArray = (SELECTION_ARRAY *) malloc(sizeof(SELECTION_ARRAY));
	    selectionArray->individuals = (INDIVIDUAL **) malloc(sizeof(INDIVIDUAL *) * length);
	    selectionArray->length = 0;
	    for (int i = 0; i < population->length; i ++) {
		for (int n = 0; n < weight[i]; n ++) {
		    selectionArray->individuals[selectionArray->length] = population->individuals[i];
		    selectionArray->length ++;
		}
	    }
	}
    }
    // Возвращем веса отбора особей
    return selectionArray;
}


void deleteSelectionArray(SELECTION_ARRAY *selectionArray)
{
    // Проверяем входные параметры
    if (selectionArray) {
	// Освобождаем память
        free(selectionArray->individuals);
	free(selectionArray);
    }
}


/* Выбор особи для отбора */
INDIVIDUAL *selection(SELECTION_ARRAY *selectionArray)
{
    INDIVIDUAL *individual = NULL;
    // Проверяем входные параметры
    if (selectionArray) {
	int index = (int) ((float) rand() / RAND_MAX * selectionArray->length);
	individual = selectionArray->individuals[index];
    }
    // Возвращаем особь
    return individual;
}


/* Кроссинговер */
INDIVIDUAL *cross(INDIVIDUAL *i1, INDIVIDUAL *i2)
{
    INDIVIDUAL *individual = NULL;
    // Проверяем входящие параметры
    if (i1 && i2 && i1->length == i2->length) {
	// Создаем новую особь
	individual = (INDIVIDUAL *) malloc(sizeof(INDIVIDUAL));
	individual->gens = (float *) malloc(sizeof(float) * i1->length);
	individual->length = i1->length;
	individual->fitness = 0;
	// В хромосоме только 1 ген - выбираем произвольно одну из особей
	if (individual->length == 1) {
	    int rnd = (int) roundf((float) rand() / RAND_MAX);
	    if (rnd)
		individual->gens[0] = i1->gens[0];
	    else
		individual->gens[0] = i2->gens[0];
	}
	// В хромосоме 2 гена - выбираем по гену из каждой особи
	else if (individual->length == 2) {
	    int rnd = (int) roundf((float) rand() / RAND_MAX);
	    if (rnd) {
		individual->gens[0] = i1->gens[0];
		individual->gens[1] = i2->gens[1];
	    }
	    else {
		individual->gens[0] = i2->gens[0];
		individual->gens[1] = i1->gens[1];
	    }
	}
	// В хромосоме больше 2 ген
	else {
	    // Определяем точку пересечения генов
	    int crossPoint = (int) floorf((float) rand() / RAND_MAX * (individual->length - 1)) + 1;
	    // Выбираем гены из первой особи
	    for (int i = 0; i < crossPoint; i ++) {
		individual->gens[i] = i1->gens[i];
	    }
	    // Выбираем гены из второй особи
	    for (int i = crossPoint; i < individual->length; i ++) {
		individual->gens[i] = i2->gens[i];
	    }
	}
    }
    return individual;
}


/* Мутация */
void mutation(INDIVIDUAL *individual)
{
    // Проверка входных парметров
    if (individual) {
	// Случайным образом изменяем значения каждого гена в особи ( в пределах 10%)
	for (int i = 0; i < individual->length; i ++) {
	    individual->gens[i] += ((float) rand() / (float) RAND_MAX - 0.5) * individual->gens[i] / 10.0;
	}
    }
}

/* Генерация следующей популяции */
POPULATION *nextPopulation(POPULATION *population)
{
    POPULATION *next = NULL;
    // Проверяем входные параметры
    if (population) {
	// Выделяем память под новую популяцию
	next = (POPULATION *) malloc(sizeof(POPULATION));
	next->length = population->length;
	next->individuals = (INDIVIDUAL **) malloc(sizeof(INDIVIDUAL *) * next->length);
	// Строим список отбора особей
	SELECTION_ARRAY *selectionArray = buildSelectionArray(population);
	// Заполняем популяцию
	for (int i = 0; i < next->length; i ++) {
	    // Выбираем две разных особи
	    INDIVIDUAL *i1 = selection(selectionArray);
	    INDIVIDUAL *i2 = selection(selectionArray);
	    while(i1 == i2) {
		i2 = selection(selectionArray);
	    }
	    // Кроссинговер
	    INDIVIDUAL *i3 = cross(i1, i2);
	    // Мутация
	    mutation(i3);
	    // Добавляем особь в популяцию
	    next->individuals[i] = i3;
	}
	// Список отбора больше не понадобится, освобождаем ресурсы
	deleteSelectionArray(selectionArray);
    }
    // Вовзращяем популяцию
    return next;
}


void bestIndividual(POPULATION *population, INDIVIDUAL *individual)
{
    // Проверяем входные параметры
    if (population && population->length > 0 &&individual) {
	// Выбираем лучшую особь
	individual->fitness = population->individuals[0]->fitness;
	memcpy(individual->gens, population->individuals[0]->gens, sizeof(float) * population->individuals[0]->length);
	for (int i = 1; i < population->length; i ++) {
	    if (population->individuals[i]->fitness < individual->fitness) {
		individual->fitness = population->individuals[i]->fitness;
		memcpy(individual->gens, population->individuals[i]->gens, sizeof(float) * population->individuals[i]->length);
	    }
	}
    }
}


int gaSearch(int populationLength, int individualLength, float min, float max, int maxPopulations, float delta, float f(float *), float *x)
{

    // Проверяем входные параметры
    if (maxPopulations < 0 || !x)
	return -1;
    // формируем первую популяцию
    POPULATION		*p = newPopulation(populationLength, individualLength, min, max);
    if (!p) return -1;
    // Рассчитываем приспособленность особей в популяции
    calcPopulationFitness(p, f);
    // Выбираем лучшую особь
    INDIVIDUAL		best;
    float		bestFitness;
    best.gens = (float *) malloc(sizeof(float) * individualLength);
    best.length = individualLength;
    bestIndividual(p, &best);
    bestFitness = best.fitness;
    memcpy(x, best.gens, sizeof(float) * best.length);
    // Ищем рещение, пока не достигнем указанной точности или не превысим количество популяций
    int iterations = 0;
    for (int i = 0; i < maxPopulations; i ++) {
	// Формируем след. популяцию
	POPULATION *p2 = nextPopulation(p);
	// Рассчитываем приспособленность особей в популяции
	calcPopulationFitness(p2, f);
	// Удаляем старую популяцию
	deletePopulation(p);
	p = p2;
	// Выбираем лучшую особь
	bestIndividual(p, &best);
	if (best.fitness < bestFitness) {
	    bestFitness = best.fitness;
	    memcpy(x, best.gens, sizeof(float) * best.length);
	}
	// Проверяем точность
	if (bestFitness < delta) {
	    break;
	}
	// Увеличиваем число итераций
	iterations ++;
    }
    // Освобождаем ресурсы
    deletePopulation(p);
    free(best.gens);
    // Возвращаем число итераций,  затраченное на поиск решения
    return iterations;
}

2 комментария

avatar
попробуй отредактировать, увеличил до 35000
avatar
Совсем другое дело :)
Только зарегистрированные и авторизованные пользователи могут оставлять комментарии.