ТАУ-Дарвинизм: реализация на Ruby

    Предисловие


    Послушайте, ворона, а может быть собака,
    А может быть корова, но тоже хорошо!
    У вас такие перья, у вас рога такие,
    Копыта очень стройные и добрая душа.

    Мультфильм «Пластилиновая ворона».

    В этой статье представляю Вашему вниманию реализацию эволюционного подхода к идентификации динамической системы. Программа написана на языке Ruby версии 1.9.2 (gems: NArray, GNUPlot). Заглянув под кат найдете пример вещественного кодирования генной информации и подходящего для него алгоритма скрещивания («flat crossover»).

    Введение


    Начну с описания задачи «черного ящика». Предположим у нас есть некий объект, передаточную функцию которого мы бы хотели узнать (упрощенную, близкую или эквивалентную к реальной). Все, что мы можем сделать, это подать на вход объекту некоторое воздействие и измерить некоторый выходной параметр (например, выходное напряжение). Как нам подглядеть, что у ящика внутри? Иными словами: как определить примерную математическую модель исследуемого объекта?
    В теории автоматического управления (ТАУ) данная задача называется идентификацией. Существуют разные способы идентификации, основанные на серьезной математике и сложные в реализации. Однако есть возможность получить приближенную математическую модель сравнительно простым способом, а именно используя стохастический подход с применение эволюционных алгоритмов.
    Подробности Вы можете прочитать в первой статье, посвященной данному вопросу: ТАУ-Дарвинизм, а здесь я остановлюсь на программной реализации.

    Реализация


    При написании алгоритмов использована библиотека «gga4r» (rubygems.org). Данная библиотека была переписана для приведение в работоспособное состояние для совместимости с ruby 1.9. Автор этой библиотеки (спасибо ему большое!) реализовал подход эволюционной обертки над сущностями, реализующими представления различных по содержанию, но одинаковых по форме эволюционирующих особей. Предлагается реализовывать свой собственный класс, имплементирующий генетические операторы. У меня получилось это так:
    class TFIndivid
        attr_reader :filter, :fitness, :dna, :num_len, :denom_len
        def initialize(num,denom,k=1.0,ts=@@ts)
          @num_len,@denom_len,@filter = num.length,denom.length,SSF::tf2ssd(num, denom, k, ts)
          @fitness = 1/(step_response_deviation(@filter)+@@eps)
          @fitness = 0.0 if @fitness.infinite?
          @dna = num + denom + [k] + [ts]
          self
        end
        def integrity_ok?(dna)
          return dna.length == @num_len + @denom_len + 2
        end
        def TFIndivid.flat_crossover(par1, par2)
          r = Random.new
          child = []
          par1.size.times { |i|
            if par1[i]>par2[i]
              child << r.rand(par2[i]..par1[i])
            else
              child << r.rand(par1[i]..par2[i])
            end
          }
          child
        end
        def recombine(other)
          child = TFIndivid::flat_crossover(self.dna, other.dna)
          res = self.integrity_ok?(child) ? TFIndivid.new(child[0..@num_len-1], child[@num_len..-3], child[-2], child[-1]) : []
          [res]
        end
        def mutate
          rnd = Random.new
          mutate_point = rnd.rand(@dna.size-1)
          @dna[mutate_point] *= rnd.rand(1.0)
        end
    end

    В примере на сайте «gga4r» класс особи наследовался от Array. Я решил поменять немного логику и завел дополнительное свойство "@dna". Сделано это было для удобства и использования конструктора, в котором происходит создание цифрового фильтра для данных коэффициента усиления, числителя и знаменателя передаточной функции (в непрерывном времени), а также периода квантования по времени. Также в конструкторе переданные параметры сливаются в один одномерный массив для работы с ними посредством генетических алгоритмов. Фитнесс особи вычисляется сразу в конструкторе и хранится в соответствующем свойстве.
    Функция «flat_crossover» реализует упомянутый в самом начале статьи вещественного скрещивания. Это, наверное, самая простая его реализация. Существует еще несколько, в которых детей реализуется несколько (здесь генерируется один потомок) и которые мне наверняка предстоит реализовать. В данном же случае порождаемое дитя формируется путем выбора для каждого гена случайного числа в интервале между значениями генов родителей.
    Оператор мутации реализован в виде умножения значения гена на случайное число в диапазоне [0.0..1.0]. Таким образом некоторые из генов «притягиваются» к нулю. Выбранный оператор скрещивания формирует в основном особи со значениями генов в средней части допустимого диапазона значений. Выбранный способ мутирования играет роль противовеса.
    Сама обертка эволюционного алгоритма работает как описано в предыдущей статье ТАУ-Дарвинизм: формируется популяция, проводится селекция, рекомбинация (скрещивание), мутация. Внутренности сможете посмотреть в исходниках программы (ссылка в конце статьи). Ее пришлось немного переписать, как написано выше, для возможности использования с ruby 1.9. Но на самом деле библиотека у меня отказалась работать и с версией 1.8 на kubuntu 10.10 (не выяснял почему). Как бы то ни было, можете сравнить с оригиналом «gga4r» — изменений немного. Кое что еще предстоит изменить.
    Фитнесс-функция реализована как среднеквадратическое отклонение между выходным сигналом тестового фильтра и сгенерированного по генам особи:
    
    def step_response_deviation(filter)
      x, test_y,y = step_response(filter)
      av_err = 0.0
      100.times {|i|
        err = test_y[i] - y[i]
        av_err += err*err
      }
      Math.sqrt(av_err / 99)
    end
    def step_response(filter)
      @@test_filter.reset
      filter.reset
      x = []; 100.times { |i| x << i*@@ts }
      test_y,y = [],[]
      100.times {
        test_y << @@test_filter.step(@@test_input)
        y << filter.step(@@test_input)
      }
      return x,test_y,y
    end
    

    Это немного отличается от постановки задачи («дана запись выхода и входа черного ящика; требуется...»), но так удобнее, а суть та же.

    Результаты


    Играясь с новой бирюлькой В ходе предварительных исследований было выяснено, что выгоднее всего задавать большой размер популяции при меньшем числе итераций. Подходящее решение находится довольно быстро при достаточном числе особей в популяции. Задачу поставил так:
    • Имеется эталонная (тестовая) передаточная функция — числитель: [1.5, 1.0], знаменатель: [0.75, 0.35, 1.1]
    • Коэффициент усиления тестового фильтра выбрал равным единице (коэффициент передаточной функции при построении дискретного фильтра к нему приводится в функции SSF::tf2ssd — модуль «ssf.rb» в исходниках)
    • Период квантования — 0.1 сек., для выбранного тестового фильтра вполне достаточно малый период
    • Размер популяции — 350 в начале и 500 максимально
    • Число итераций — 50

    Оказалось, что не смотря на заданный максимум размер популяции превышал заданное значение. После запуска получил следующий результат:
    $ ruby1.9.1 main.rb
    Iteration: 1
    Population size: 387
    best Individ: [ 1.552605923140199 1.8418549861379885 4.4154048910142265 5.556332886121151 5.320061971918478 1.6334223860651258 0.1 ]
    best fitness: 2.6508867715774027
    mean fitness: 1.1722382322077871
    mean fitness recalc 1.1722382322077871
    ******************************
    Iteration: 2
    Population size: 421
    best Individ: [ 0.37594352736713244 0.41988779230454 6.092873656120998 3.575618643442387 8.786369744495486 1.3186424599689655 0.1 ]
    best fitness: 2.9780762386062767
    mean fitness: 1.381001802845182
    mean fitness recalc 1.381001802845182
    ******************************
    Iteration: 3
    Population size: 456
    best Individ: [ 5.355702236617921 2.9021864520499134 3.979973441887163 1.1587534831116912 5.818934614234973 0.6934218845420625 0.1 ]
    best fitness: 4.008150118088937
    mean fitness: 1.6192001539010286
    mean fitness recalc 1.6192001539010286
    ******************************
    Iteration: 4
    Population size: 501
    best Individ: [ 6.981596537433821 7.0034200556235495 6.592795290155172 3.1179936580946577 7.236199584387684 1.4664747534746492 0.1 ]
    best fitness: 4.629268673426835
    mean fitness: 1.8244067682429481
    mean fitness recalc 1.8244067682429481
    ******************************
    Iteration: 5
    Population size: 537
    best Individ: [ 5.917012124258587 6.182284010197045 5.214809836688884 3.9333556599339055 8.349166334410114 1.544147649767568 0.1 ]
    best fitness: 4.900537478511452
    mean fitness: 2.024810983956715
    mean fitness recalc 2.024810983956715
    ******************************
    Iteration: 6
    Population size: 551
    best Individ: [ 5.917012124258587 6.182284010197045 5.214809836688884 3.9333556599339055 8.349166334410114 1.544147649767568 0.1 ]
    best fitness: 4.900537478511452
    mean fitness: 2.3489882909253152
    mean fitness recalc 2.3489882909253152
    ******************************
    Iteration: 7
    Population size: 540
    best Individ: [ 7.667944581517669 7.685145812466274 3.6039889820553768 2.4390348703480536 5.087422184655714 1.4689450977845646 0.1 ]
    best fitness: 8.166839319133336
    mean fitness: 2.652950357867813
    mean fitness recalc 2.652950357867813
    ******************************
    Iteration: 8
    Population size: 552
    best Individ: [ 7.667944581517669 7.685145812466274 3.6039889820553768 2.4390348703480536 5.087422184655714 1.4689450977845646 0.1 ]
    best fitness: 8.166839319133336
    mean fitness: 2.9467084852969725
    mean fitness recalc 2.9467084852969725
    ******************************
    Iteration: 9
    Population size: 553
    best Individ: [ 5.7806312313223485 4.681216450534634 5.384196754050039 3.0234885021495357 7.238140602489246 1.1815749271908358 0.1 ]
    best fitness: 9.966075154360498
    mean fitness: 3.335389533280396
    mean fitness recalc 3.335389533280396
    ******************************
    Iteration: 10
    Population size: 543
    best Individ: [ 5.7806312313223485 4.681216450534634 5.384196754050039 3.0234885021495357 7.238140602489246 1.1815749271908358 0.1 ]
    best fitness: 9.966075154360498
    mean fitness: 3.818330258524996
    mean fitness recalc 3.818330258524996
    ******************************
    Iteration: 11
    Population size: 555
    best Individ: [ 5.7806312313223485 4.681216450534634 5.384196754050039 3.0234885021495357 7.238140602489246 1.1815749271908358 0.1 ]
    best fitness: 9.966075154360498
    mean fitness: 4.241020711174743
    mean fitness recalc 4.241020711174743
    ******************************
    Iteration: 12
    Population size: 563
    best Individ: [ 7.114438733870453 5.7051128886012386 3.8495059720675098 1.8164317233081309 5.543227442940334 1.1352374472380862 0.1 ]
    best fitness: 11.644557896655645
    mean fitness: 4.980294253355436
    mean fitness recalc 4.980294253355436
    ******************************
    .......
    ******************************
    Iteration: 23
    Population size: 542
    best Individ: [ 5.949029964618817 4.681516090696933 4.8007731550065325 2.562344672524173 6.969591955243219 1.1684490266387202 0.1 ]
    best fitness: 20.468430876331446
    mean fitness: 11.25191007878279
    mean fitness recalc 11.25191007878279
    ******************************
    ......
    ******************************
    Iteration: 36
    Population size: 550
    best Individ: [ 6.304563716617802 4.937028188768713 5.01636959411469 2.5019167234335344 7.143946485486306 1.1690443373254094 0.1 ]
    best fitness: 30.361894930102864
    mean fitness: 22.65217005076974
    mean fitness recalc 22.65217005076974
    ******************************
    ........
    ******************************
    Iteration: 46
    Population size: 542
    best Individ: [ 6.548169493382168 5.060607979367818 4.739775170940725 2.3444848651883947 6.741642980854007 1.1564890916088506 0.1 ]
    best fitness: 32.826565607333514
    mean fitness: 25.606090666108418
    mean fitness recalc 25.606090666108418
    ******************************
    Iteration: 47
    Population size: 550
    best Individ: [ 6.548169493382168 5.060607979367818 4.739775170940725 2.3444848651883947 6.741642980854007 1.1564890916088506 0.1 ]
    best fitness: 32.826565607333514
    mean fitness: 25.489911664932283
    mean fitness recalc 25.489911664932283
    ******************************
    Iteration: 48
    Population size: 561
    best Individ: [ 6.548169493382168 5.060607979367818 4.739775170940725 2.3444848651883947 6.741642980854007 1.1564890916088506 0.1 ]
    best fitness: 32.826565607333514
    mean fitness: 25.44926190395142
    mean fitness recalc 25.44926190395142
    ******************************
    Iteration: 49
    Population size: 550
    best Individ: [ 3.1439340250602816 5.060607979367818 4.739775170940725 2.3444848651883947 6.741642980854007 1.1564890916088506 0.1 ]
    best fitness: 32.826565607333514
    mean fitness: 26.108921277427932
    mean fitness recalc 26.108921277427932
    ******************************
    Iteration: 50
    Population size: 559
    best Individ: [ 3.1439340250602816 5.060607979367818 4.739775170940725 2.3444848651883947 6.741642980854007 1.1564890916088506 0.1 ]
    best fitness: 32.826565607333514
    mean fitness: 26.23635084771665
    mean fitness recalc 26.23635084771665
    ******************************

    Первые два числа в «best Individ» — числитель передаточной функции, следующие три — знаменатель, далее — масштабный коэффициент и период квантования. Наилучший фитнесс, равный 32 «попугаям» — достаточно хороший результат. В подтверждение следующий рисунок:

    На нем изображены графики переходных процессов тестового и полученного в ходе эволюции фильтров при подаче им на вход единичного ступенчатого воздействия. «Кто» из графиков «Кто» пока не могу сказать — недопилил функциональность по работе с gnuplot, но это и не имеет значение в данном случае, т.к. невооруженным глазом видно, что результаты очень близки.

    Заключение


    Программа реализована в ходе исследовательских работ по программе «У.М.Н.И.К». Код сырой, но работающий (по крайней мере у меня). В планах оснастить программу неким простеньким GUI для удобства. Ниже ссылка на код.
    UPD1: исходники размещены на GitHub.
    Репозиторий на GitHub
    Поделиться публикацией

    Похожие публикации

    Комментарии 8

      +2
      Интересно. Я бы порекомендовал вам выкладывать код на github. Кто-нибудь мог бы помочь избавиться от «сырости» кода :)
        0
        Спасибо за рекомендацию. Поправлено.
        0
        > Существуют разные способы идентификации, основанные на серьезной математике и сложные в реализации.

        > Фитнесс-функция реализована как среднеквадратическое отклонение между выходным сигналом тестового фильтра и сгенерированного по генам особи

        Вот уж не сказал бы, что банальный метод наименьших квадратов сложнее в реализации генетического алгоритма.

        Как еще одно применение генетического алгоритма — не оправдано и не интересно.
          0
          Поправлюсь — мне не интересно, а тем кто только начинает изучать тему, вполне может быть интересно как простой примерю
            0
            Здесь МНК не противопоставляется ГА, а наоборот используется как критерий оценки приспособленности. Генетика использована для минимизации функционала (в данном случае среднеквадратическое отклонение, а может быть и другой критерий). И почему использование генетики неоправдано? Я этой программой автоматизирую поиск оптимальной мат. модели. Реализация проста для понимания. Работает эффективно. Следующий шаг — натравить генетику на фильтр Калмана с ковариационными матрицами и его большой размерностью. Думаю как минимум интересно, но может и результат будет положительным.
          0
          x = []; 100.times { |i| x << i*@@ts } => x = Array.new(100) {|i| i*@@ts} и аналогично с test_y и т.п.
            0
            Красивше… согласен. А с точки зрения производительности есть выигрыш?
              +1
              В первом сразу же выделяется нужное количесво памяти. Да и заполнение быстрее будет.

          Только полноправные пользователи могут оставлять комментарии. Войдите, пожалуйста.

          Самое читаемое