ArrayFire机器学习—遗传算法

ArrayFire机器学习—遗传算法

 

1. 、遗传算法基本原理

遗传算法（Genetic Algorithm，GA）是借鉴生物界的进化规律（适者生存，优胜劣汰，遗传机制）演化而来的随机搜索方法。它是由美国的J.Holland教授1975年首先提出，其主要特点是直接对结构对象进行操作，不存在求导和函数连续性的限定；具有内在的隐并行性和更好的全局寻优能力；采用概率化的寻优方法，能自动获取和指导优化的搜索空间，自适应地调整搜索方向，不需要确定的规则。

遗传算法的基本运算过程如下：

a)初始化：设置进化代数计数器t=0，设置最大进化代数T，随机生成M个个体作为初始群体P(0)。

b)个体评价：计算群体P(t)中各个个体的适应度（fitness）。

c)选择运算:将选择算子作用于群体。选择的目的是把优化的个体直接遗传到下一代或通过配对交叉产生新的个体再遗传到下一代。选择操作是建立在群体中个体的适应度评估基础上的。

d)交叉运算：将交叉算子作用于群体。遗传算法中起核心作用的就是交叉算子。

e)变异运算：将变异算子作用于群体。即是对群体中的个体串的某些基因座上的基因值作变动。

群体P(t)经过选择、交叉、变异运算之后得到下一代群体P(t+1)。

f)终止条件判断:若t=T,则以进化过程中所得到的具有最大适应度个体作为最优解输出，终止计算。

上述中选择（selection）、交叉（crossover），变异（mutation）都是遗传操作。

选择操作的目的是把优化的个体直接遗传到下一代或通过配对产生新的个体再遗传到下一代。目前常用的选择算子有：适应度比例法、随机遍历抽样法、局部选择法，比如pi=fij=1nfi。

标准的比较和选择过程是生成新的点集和舍弃旧的点集，生成和舍弃都是基于现有种群的函数。这些GA算法中用已有的单个点集的组合推导出新的点集。

 

1. 、遗传算法的主要函数

#include <arrayfire.h>
#include <climits>
#include <cstdio>
#include <cstring>
#include <ctime>

using namespace af;
static const float DefaultTopFittest = 0.5;

array update(const array& searchSpace, const array& sampleX,
    const array& sampleY, const int n) {
    return searchSpace(sampleY * n + sampleX);
}

array selectFittest(const array& sampleZ, const int nSamples,
    const float topFit = DefaultTopFittest) {
    // pick top fittest
    array indices, values;
    sort(values, indices, sampleZ);
    int topFitElem = topFit * nSamples;
    int n = indices.elements();
    return (n > topFitElem) ? indices(seq(n - topFitElem, n - 1)) : indices;
}

void reproduce(array& searchSpace, array& sampleX, array& sampleY,
    array& sampleZ, const int nSamples, const int n) {
    // Get fittest parents
    array selection = selectFittest(sampleZ, nSamples);
    array parentsX = sampleX(selection);
    array parentsY = sampleY(selection);
    int bits = (int)log2(n);

    // Divide selection in two
    array parentsX1 = parentsX.rows(0, parentsX.elements() / 2 - 1);
    array parentsX2 =
        parentsX.rows(parentsX.elements() / 2, parentsX.elements() - 1);
    array parentsY1 = parentsY.rows(0, parentsY.elements() / 2 - 1);
    array parentsY2 =
        parentsY.rows(parentsY.elements() / 2, parentsY.elements() - 1);

    // Get crossover points (at which bit to crossover) and construct bit masks
    // from them
    array crossover = randu(nSamples / 4, u32) % bits;
    array lowermask = (1 << crossover) - 1;
    array uppermask = INT_MAX - lowermask;

    // Create children as the cross between two parents
    array childrenX1 = (parentsX1 & uppermask) + (parentsX2 & lowermask);
    array childrenY1 = (parentsY1 & uppermask) + (parentsY2 & lowermask);

    array childrenX2 = (parentsX2 & uppermask) + (parentsX1 & lowermask);
    array childrenY2 = (parentsY2 & uppermask) + (parentsY1 & lowermask);

    // Join two new sets
    sampleX = join(0, childrenX1, childrenX2);
    sampleY = join(0, childrenY1, childrenY2);

    // Create mutant children
    array mutantX = sampleX;
    array mutantY = sampleY;

    // Flip a random bit to vary the gene pool
    mutantX = mutantX ^ (1 << (randu(nSamples / 2, u32) % bits));
    mutantY = mutantY ^ (1 << (randu(nSamples / 2, u32) % bits));

    sampleX = join(0, sampleX, mutantX);
    sampleY = join(0, sampleY, mutantY);

    // Update the value of each sample with the new coordinates
    sampleZ = update(searchSpace, sampleX, sampleY, n);
}

//初始化样本
void initSamples(array& searchSpace, array& sampleX, array& sampleY,
    array& sampleZ, const int nSamples, const int n) {
    setSeed(time(NULL));
    sampleX = randu(nSamples, u32) % n;
    sampleY = randu(nSamples, u32) % n;
    sampleZ = update(searchSpace, sampleX, sampleY, n);
}

void init(array& searchSpace, array& searchSpaceXDisplay,
    array& searchSpaceYDisplay, array& sampleX, array& sampleY,
    array& sampleZ, const int nSamples, const int n) {
    // initialize space
    searchSpace = range(dim4(n / 2, n / 2), 0) + range(dim4(n / 2, n / 2), 1);
    searchSpace = join(0, searchSpace, flip(searchSpace, 0));
    searchSpace = join(1, searchSpace, flip(searchSpace, 1));

    // initialize display data
    searchSpaceXDisplay = iota(dim4(n, 1), dim4(1, n));
    searchSpaceYDisplay = iota(dim4(1, n), dim4(n, 1));

    // initalize searchers
    initSamples(searchSpace, sampleX, sampleY, sampleZ, nSamples, n);
}

//每一代的打印
void reproducePrint(float& currentMax, array& searchSpace, array& sampleX,
    array& sampleY, array& sampleZ, const float trueMax,
    const int nSamples, const int n) {
    if (currentMax < trueMax * 0.99) {
        float maximum = max<float>(sampleZ);
        array whereM = where(sampleZ == maximum);
        if (maximum < trueMax * 0.99) {
            printf("Current max at ");
        }
        else {
            printf("\nMax found at ");
        }
        printf("(%d,%d): %f (trueMax %f)\n",
            sampleX(whereM).scalar<unsigned int>(),
            sampleY(whereM).scalar<unsigned int>(), maximum, trueMax);
        currentMax = maximum;
        reproduce(searchSpace, sampleX, sampleY, sampleZ, nSamples, n);
    }
}

//遗传搜索--主控函数
void geneticSearch(bool console, const int nSamples, const int n) {
    array searchSpaceXDisplay = 0;
    array searchSpaceYDisplay = 0;
    array searchSpace;
    array sampleX;
    array sampleY;
    array sampleZ;

    init(searchSpace, searchSpaceXDisplay, searchSpaceYDisplay, sampleX,
        sampleY, sampleZ, nSamples, n);
    float trueMax = max<float>(searchSpace);
    float maximum = -trueMax;

    if (!console) {
        af::Window win(1600, 800, "Arrayfire Genetic Algorithm Search Demo");
        win.grid(1, 2);
        do {
            reproducePrint(maximum, searchSpace, sampleX, sampleY, sampleZ,
                trueMax, nSamples, n);
            win(0, 0).setAxesTitles("IdX", "IdY", "Search Space");
            win(0, 1).setAxesTitles("IdX", "IdY", "Search Space");
            win(0, 0).surface(searchSpaceXDisplay, searchSpaceYDisplay,
                searchSpace);
            win(0, 1).scatter(sampleX.as(f32), sampleY.as(f32), sampleZ.as(f32),
                AF_MARKER_CIRCLE);
            win.show();
        } while (!win.close());
    }
    else {
        do {
            reproducePrint(maximum, searchSpace, sampleX, sampleY, sampleZ,
                trueMax, nSamples, n);
        } while (maximum < trueMax * 0.99);
    }
}

1. 、main函数代码

int main(int argc, char** argv) {
    bool console = false;
    const int n = 32;
    const int nSamples = 16;
    if (argc > 2 || (argc == 2 && strcmp(argv[1], "-"))) {
        printf("usage: %s [-]\n", argv[0]);
        return -1;
    }
    else if (argc == 2 && argv[1][0] == '-') {
        console = true;
    }

    try {
        af::info();
        printf("** ArrayFire Genetic Algorithm Search Demo **\n\n");
        printf(
            "Search for trueMax in a search space where the objective function "
            "is defined as :\n\n");
        printf("SS(x ,y) = min(x, n - (x + 1)) + min(y, n - (y + 1))\n\n");
        printf("(x, y) belongs to RxR; R = [0, n); n = %d\n\n", n);
        if (!console) {
            printf("The left figure shows the objective function.\n");
            printf(
                "The figure on the right shows current generation's parameters "
                "and function values.\n\n");
        }
        geneticSearch(console, nSamples, n);
    }
    catch (af::exception& e) { fprintf(stderr, "%s\n", e.what()); }

    return 0;
}

1. 、示例结果

ArrayFire v3.8.0 (CUDA, 64-bit Windows, build d99887a)
Platform: CUDA Runtime 10.1, Driver: 11000
[0] GeForce GTX 1050 Ti, 4096 MB, CUDA Compute 6.1
** ArrayFire Genetic Algorithm Search Demo **

Search for trueMax in a search space where the objective function is defined as :

SS(x ,y) = min(x, n - (x + 1)) + min(y, n - (y + 1))

(x, y) belongs to RxR; R = [0, n); n = 32

The left figure shows the objective function.
The figure on the right shows current generation's parameters and function values.

Current max at (18,14): 27.000000 (trueMax 30.000000)
Current max at (18,14): 27.000000 (trueMax 30.000000)

Max found at (16,15): 30.000000 (trueMax 30.000000)