Atmosphere GM - Statistical package

Project site: 9il.github.io/atmosphere_gm

Features

1. Separating mixtures of probability distributions
2. Grid methods
3. Parameterized algorithms
4. Optimization over sliding window

Documentation

Documentation can be found here.

Installation

To use this package, put the following dependency into your project's dub.json into the dependencies section:

{
    ...
    "dependencies": {
        "atmosphere_gm": ">=0.0.1"
    }
}

Benchmarking

It is suggested the llvm D compiler be used for benchmarks. Project requires LDC version >= 0.15.0 or DMD >= 2.066, or corresponding GDC release. The DMD is easy way to start.

atmosphere_gm. Российский анонс библиотеки быстрых сеточных методов разделения смесей вероятностных распределений.

Всем добрый день и с Наступившим!

Библиотека

Основной алгоритм разделения смесей сейчас в значительной степени переработан и улучшен.
Обобщен шаблонами на случай смесей N-параметрических распределений. Пользователю нужно задать размерность и необходимую функцию плотности соответствующей размерности.
Обобщен для использования скользящего окна с неравномерным распределением событий по времени. Практически каждый этап итерационного шага управляем.

Библиотека разделена на пакеты:

1. kernel - весь BLAS и подобный API, который можно оптимизировать. Также утилиты работы с многомерной кластеризацией на гиперторах и многомерными индексами, поиском корней функции одной переменной.
2. fractal - описание фракталов.
3. cyclone - алгоритмы независимые от платформы
1. phiallocator - максимально оптимизированный и гибкий модуль управления памятью при поточном получении данных. Оптимально работает в том числе с данными, которые поступают неравномерными порциями. Данный пункт очень важен: количество перемещений не должно превышать количество новых данных, и при этом нужно иметь в памяти плоскую матрицу для оптимальной работы в том числе с ее транспонированной версией. Управление детерминировано по времени и может использоваться в программах реального времени.
2. generator - итератор по фракталам
3. eye - модуль структуры одной итерации. Работает с информацией с предыдущих шагов при принятии решений (то, что нам известно с предыдущего шага мы не пересчитываем, что во многом зависит от обстоятельств работы и аргументов итерации). В нем сосредоточена вся логика алгоритма. В зависимости от того, как им пользоваться можно получить совершенно разные алгоритмы, в том числе родительский алгоритм Алексея Назарова.
4. cyclone - примеры управления структурой одной итерации. В настоящий момент сделан двухпараметрический пример для смесей нормальных законов. Всего 150 строчек задают все параметры, весь алгоритм от чтения файла до полного управления структурой одной итерации и параметров сетки, и описания отображения собственно нормального распределения на параметрическую сеть.

BLAS API стандартный. Возможна работы без BLAS вовсе.

Технические возможности оптимизации и использования API

1. Библиотека написана на языке программирования D (компиляторы dmd (dmc-32 или win64), gdc (с версии 4.9 часть gcc), ldc (llvm компилятор, включен в анонс llvm 3.4)). Язык имеет С прямой API - можно интегрировать с openBLAS, cuBLAS, Intel MKL, MATLAB, Julia, R и т.д.
2. Возможно использовать библиотеку в системах реального времени. Память управляется самим алгоритмом и не участвует в сборке мусора, поэтому сборщик можно вообще отключить.
3. В плане оптимизации кода компиляторами библиотека полностью эквивалента аналогичной программе написанной на чистом С. Возможно использовать SIMD или встроеный ассемблер. Весь код, который может быть подвержен оптимизации на уровне железа вынесен в пакет kernel.
4. Библиотека платформнонезависима. Сейчас работает на всех основных ОС для ПК или серверов. Потом добавлю поддержку смартфонов :-)
5. Вы можете задать свою функцию правдоподобия в модуле kernel, для этого также нужно сделать согласованную версию кластеризации. Изменение алгоритма при этом не потребуется. По умолчанию алгоритм можно запускать с двумя версиями функции правдоподобия: логарифм произведения сумм (быстрая), или сумма логарифмов сумм (всегда точная).

Стабильность

Программа стабильна и протестирована до уровня RC (release candidate). Сама структура библиотеки располагает к быстрому обнаружению ошибок. Для тестирования использовалась специально созданная система визуализации. В прямом смысле алгоритм тестировался и визуализировался на каждом шаге внутри итерации. При этом API вероятно будет немного меняться согласно пожеланиям пользователей.

Сроки

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

Дополнительно

• Разработана библиотека на основе openGL для визуализации данных в реальном времени. Сейчас 2D, будет и 3D к концу учебного года (к примеру для обобщенных гамма, картинки обещают быть красивыми).
• Биндинг blas к D.

TODO

1. Поддержка для kernel различных драйверов:
1. Intel MKL (в том числе SIMD реализация функции логарифма для ФП и других)
2. cuBLAS
3. обобщенная реализация SIMD, в том числе AVX2 и AVX3.
2. 3D визуализация

Оптимизация EM-алгоритмов для вероятностного тематического моделирования на примере PLSA-EM-алгоритма

Как правило тестирование EM алгоритмов затратно по времени. Данная статья призвана обозначить основные способы оптимизации и готовые инструменты, применимые к широкому спектру алгоритмов, работающих с матрицами на примере PLSA-EM алгоритма. 

Матричная форма

Последующее изложение происходит в терминологии "Модификации EM-алгоритма для вероятностного тематического моделирования" К. В. Воронцов, А. А. Потапенко. В ней вы можете ознакомиться с постановкой задачи и с исходным алгоритмом, который далее будет представлен в матричной форме.

Матричное представление удобно по следующим причинам:

1. Простота формулировки, записи и программной реализации. 
2. Алгоритмы работы с матрицами хорошо изучены и оптимизированы.

Основные оптимизации

SIMD

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

Multithreading

Использование преимущества многоядерных процессоров само по себе довольно просто. В тестах оптимизированного тривиального алгоритма перемножения матриц, приведенных в данной статье, используется группировка перемножения матриц по внешнему циклу. Это позволяет создавать в 2 раза меньше новых потоков.

BLAS 

BLAS (Basic Linear Algebra Subprograms) является наиболее распространенным API для работы с матрицами, его реализации доступны как под свободными, так и под коммерческими лицензиями. BLAS разработано и оптимизировано для широкого спектра устройств, начиная со смартфонов и заканчивая суперкомпьютерами. Существуют реализации для графических ускорителей, их тесты можно найти в работе "BLAS Comparison on FPGA, CPU and GPU" Srinidhi Kestury, John D. Davis, Oliver Williams.

OpenBLAS

Наиболее удачным и производительным вариантом для настольных компьютеров (и ноутбуков) является реализация OpenBLAS, являющаяся потомком gotoBLAS2, разработанной Kazushige Goto . В основе реализации, кроме алгоритмов быстрого умножения матриц, использования SIMD инструкций и многопоточности, лежит оптимизация алгоритмов с учетом особенностей кеширования памяти для конкретной процесcорной архитектуры "Anatomy of high-performance matrix multiplication", "High-performance implementation of the level-3 BLAS" Kazushige Goto. Кроме того, ядро OpenBLAS имеет различные драйвера в зависимости от версии процессора. OpenBLAS опубликовано под лицензией BSD.

Оболочки для BLAS

MATLAB

MATLAB по умолчанию использует BLAS, представленный в Intel MKL. Существует реализация BLAS от AMD - ACML, ее также возможно подключить к MATLAB. При этом нужно обратить особое внимание на механизмы реализации многопоточности в конкретной версии MATLAB и ACML (зависит от компилятора).

R

Вероятно самой гибкой оболочкой для BLAS является язык R. В работе "Benchmarking Single- and Multi-Core BLAS Implementations and GPUs for use with R" Dirk Eddelbuettel можно найти сравнение реализаций GotoBLAS, Intel MKL, ATLAS и других. Тут стоит отметить, что, чем производительнее BLAS, тем больше оптимизации в ней сделано под конкретные процессоры.

Julia

Julia - новый язык с JIT компилятором на основе LLVM. В качестве ядра использует OpenBLAS. Пройдя по ссылке можно увидеть таблицу (приведена ниже) сравнения производительноcти, в том числе с MATLAB.

Тестирование

Для EM-алгоритмов, их модификаций иногда необходимо отойти от матричной формулуровки. Также интерес представляет возможность работы с разряженными матрицами. Преимущество тривиальной реализации заключается в возможности ее модифицирования под смежные алгоритмы. Поэтому в данной статье сравниваются тривиальная и BLAS реализации:

1. Тривиальная реализация на основе матричной формулировки.
2. Многопоточная тривиальная реализация на основе матричной формулировки, с оптимизированным скалярным произведением, суммой (используется при нормировании).
3. BLAS
4. BLAS c использованием оптимизированной суммы.

Платформа

Processor	Intel i5-4570, Haswell
Instruction set	AVX2
System	Ubuntu 13.10
Compiler	GDC-4.8.1
Compiler flags	-march=native -fno-bounds-check -frename-registers -frelease -O3

Реализация

Структура матрицы (язык D)

Для чистоты эксперимента здесь представлена простая плоская реализация двумерной матрицы. Производительность данной реализации равна аналогичной на C++.

Тривиальная реализация

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

Оптимизированная тривиальная реализация

Оптимизированная тривиальная реализация абсолютна аналогична тривиальной, за исключением:

1. Используется оптимизированное скалярное произведение
2. Используется распараллеливание по каждому внешнему циклу.

BLAS реализация

Используется многопоточный OpenBLAS 2.8, оптимизированный под архитектуру Intel Haswell.

BLAS + fast sum реализация

Аналогична BLAS, но используется оптимизированное суммирование.

Результаты

Было проведено 64 теста с различным количеством тем, документов и слов для каждой реализации. Результаты приведены в следующей таблице:
Отдельно вынесены 2 графика

Заключение

BLAS является достаточно мощным иструментом, с которым можно работать как нативно через компилируемые языки, так и через оболочки, такие как MATLAB, R и Julia. Однако, при оптимизации тривиального алгоритма можно достичь сопостовимой производительности, что позволяет оптимизировать смежные версии алгоритмов.

SIMD implementation of dot-product.Benchmarks

Dot product can be auto-vectorized only with fast-math or similar compiler option. fast-math allows compiler to use associative floating point transformations. Other math functions like exponent can be damaged consequently.
Possible solutions:

1. Compile fast_math code from other program separately and then link it. This is easy solution. However this is a step back to C.
2. To introduce a @fast_math attribute. This is hard to realize. But I hope this will be done for future compilers.
3. Do vectorization yourself. In that case you need to realize SIMD accessory functions like unaligned load.

So I wrote implementations of dot product for real and complex numbers.

Code

Dot product for real numbers:

Complex version is similar. The main loop:

There are gdc and ldc implementations, but ldc implementation was not tested. Source code is available at GitHub.

Benchmarks

Processor	Intel i5-4570, Haswell
Instruction set	AVX2
System	Ubuntu 13.10
Compiler	GDC-4.8.1
Compiler flags	-march=native -fno-bounds-check -frename-registers -frelease -O3