Математическая оптимизация в C#

Если это для функции активации, имеет ли это большое значение, если расчет e^x полностью точен?

double eapprox = (1d + x / 256d);
eapprox *= eapprox;
eapprox *= eapprox;
eapprox *= eapprox;
eapprox *= eapprox;
eapprox *= eapprox;
eapprox *= eapprox;
eapprox *= eapprox;
eapprox *= eapprox;

продолжайте удвоение или деление пополам 256 (и добавление/удаление умножения) в зависимости от того, насколько точным вы хотите, чтобы приближение было. Даже при n=4 он по-прежнему дает около 1,5 десятичных знаков точности для значений x между -0,5 и 0,5 (и появляется в 15 раз быстрее, чем математика.опыт.))(

P. S. Я забыл упомянуть, -- вы, очевидно, не следует действительно разделить на 256: умножить на константу 1/256. JIT-компилятор Java делает эту оптимизацию автоматически (по крайней мере, Hotspot делает), и я предполагал, что C# тоже должен это делать.

посмотреть этот пост. он имеет приближение для e^x, написанного на Java, это должен быть код C# для него (непроверенный):

public static double Exp(double val) {  
    long tmp = (long) (1512775 * val + 1072632447);  
    return BitConverter.Int64BitsToDouble(tmp << 32);  
}

в моих бенчмарках это больше чем в 5 раз быстрее, чем математика.exp () (в Java). Аппроксимация основана на документе"быстрое, компактное приближение экспоненциальной функции", который был разработан именно для использования в нейронных сетях. Это в основном то же самое, что и таблица поиска 2048 записей и линейное приближение между записями, но все это с помощью трюков с плавающей точкой IEEE.

EDIT: по данным Специальный Соус это ~3.25 x быстрее, чем реализация среды CLR. Спасибо!

1. помните, что любые изменения в этой функции активации происходят за счет различного поведения. Это даже включает переключение на поплавок (и, таким образом, снижение точности) или использование заменителей активации. Только экспериментируя с вашим вариантом использования покажет правильный путь.
2. В дополнение к простой оптимизации кода, я бы также рекомендовал рассмотреть распараллеливания вычислений (т. е.: использовать несколько ядер вашего компьютера или даже машины в облаках Windows Azure) и улучшение алгоритмов обучения.

обновление:сообщение в таблицах поиска для функций активации ANN

обновление 2: Я удалил точку на LUTs, так как я перепутал их с полным хэшированием. Спасибо перейти к Хенрик Густафссон за то, что вернул меня на трассу. Таким образом, память не является проблемой, хотя пространство поиска все еще перепутано с локальным крайности немного.

при 100 миллионах вызовов я бы начал задаваться вопросом, не искажает ли накладные расходы профилировщика ваши результаты. Замените расчет на no-op и посмотрите, является ли он еще сообщается, что потребляет 60% времени выполнения...

или еще лучше, создайте некоторые тестовые данные и используйте таймер секундомера для профилирования миллиона или около того вызовов.

Если вы можете взаимодействовать с C++, вы можете рассмотреть возможность хранения всех значений в массиве и перебирать их с помощью SSE следующим образом:

void sigmoid_sse(float *a_Values, float *a_Output, size_t a_Size){
    __m128* l_Output = (__m128*)a_Output;
    __m128* l_Start  = (__m128*)a_Values;
    __m128* l_End    = (__m128*)(a_Values + a_Size);

    const __m128 l_One        = _mm_set_ps1(1.f);
    const __m128 l_Half       = _mm_set_ps1(1.f / 2.f);
    const __m128 l_OneOver6   = _mm_set_ps1(1.f / 6.f);
    const __m128 l_OneOver24  = _mm_set_ps1(1.f / 24.f);
    const __m128 l_OneOver120 = _mm_set_ps1(1.f / 120.f);
    const __m128 l_OneOver720 = _mm_set_ps1(1.f / 720.f);
    const __m128 l_MinOne     = _mm_set_ps1(-1.f);

    for(__m128 *i = l_Start; i < l_End; i++){
        // 1.0 / (1.0 + Math.Pow(Math.E, -value))
        // 1.0 / (1.0 + Math.Exp(-value))

        // value = *i so we need -value
        __m128 value = _mm_mul_ps(l_MinOne, *i);

        // exp expressed as inifite series 1 + x + (x ^ 2 / 2!) + (x ^ 3 / 3!) ...
        __m128 x = value;

        // result in l_Exp
        __m128 l_Exp = l_One; // = 1

        l_Exp = _mm_add_ps(l_Exp, x); // += x

        x = _mm_mul_ps(x, x); // = x ^ 2
        l_Exp = _mm_add_ps(l_Exp, _mm_mul_ps(l_Half, x)); // += (x ^ 2 * (1 / 2))

        x = _mm_mul_ps(value, x); // = x ^ 3
        l_Exp = _mm_add_ps(l_Exp, _mm_mul_ps(l_OneOver6, x)); // += (x ^ 3 * (1 / 6))

        x = _mm_mul_ps(value, x); // = x ^ 4
        l_Exp = _mm_add_ps(l_Exp, _mm_mul_ps(l_OneOver24, x)); // += (x ^ 4 * (1 / 24))

#ifdef MORE_ACCURATE

        x = _mm_mul_ps(value, x); // = x ^ 5
        l_Exp = _mm_add_ps(l_Exp, _mm_mul_ps(l_OneOver120, x)); // += (x ^ 5 * (1 / 120))

        x = _mm_mul_ps(value, x); // = x ^ 6
        l_Exp = _mm_add_ps(l_Exp, _mm_mul_ps(l_OneOver720, x)); // += (x ^ 6 * (1 / 720))

#endif

        // we've calculated exp of -i
        // now we only need to do the '1.0 / (1.0 + ...' part
        *l_Output++ = _mm_rcp_ps(_mm_add_ps(l_One,  l_Exp));
    }
}

однако помните, что массивы, которые вы будете использовать, должны быть выделены с помощью _aligned_malloc(some_size * sizeof(float), 16), потому что SSE требует памяти, выровненной по границе.

используя SSE, я могу вычислить результат для всех 100 миллионов элементов примерно за полсекунды. Однако, выделяя столько памяти за один раз это будет стоить вам почти две трети гигабайта, поэтому я бы предложил обрабатывать больше, но меньшие массивы за раз. Возможно, вы даже захотите использовать подход двойной буферизации с элементами 100K или более.

кроме того, если количество элементов начинает значительно расти, вы можете захотеть обработать эти вещи на GPU (просто создайте текстуру 1D float4 и запустите очень тривиальный шейдер фрагментов).

FWIW, вот мои тесты C# для уже опубликованных ответов. (Пустая функция, которая просто возвращает 0, чтобы измерить накладные расходы вызова функции)

Empty Function:       79ms   0
Original:             1576ms 0.7202294
Simplified: (soprano) 681ms  0.7202294
Approximate: (Neil)   441ms  0.7198783
Bit Manip: (martinus) 836ms  0.72318
Taylor: (Rex Logan)   261ms  0.7202305
Lookup: (Henrik)      182ms  0.7204863

public static object[] Time(Func<double, float> f) {
    var testvalue = 0.9456;
    var sw = new Stopwatch();
    sw.Start();
    for (int i = 0; i < 1e7; i++)
        f(testvalue);
    return new object[] { sw.ElapsedMilliseconds, f(testvalue) };
}
public static void Main(string[] args) {
    Console.WriteLine("Empty:       {0,10}ms {1}", Time(Empty));
    Console.WriteLine("Original:    {0,10}ms {1}", Time(Original));
    Console.WriteLine("Simplified:  {0,10}ms {1}", Time(Simplified));
    Console.WriteLine("Approximate: {0,10}ms {1}", Time(ExpApproximation));
    Console.WriteLine("Bit Manip:   {0,10}ms {1}", Time(BitBashing));
    Console.WriteLine("Taylor:      {0,10}ms {1}", Time(TaylorExpansion));
    Console.WriteLine("Lookup:      {0,10}ms {1}", Time(LUT));
}

с моей головы,в этой статье объясняется способ аппроксимации экспоненты путем злоупотребления плавающей точкой, (нажмите на ссылку в правом верхнем углу для PDF), но я не знаю, будет ли это очень полезно для вас в. NET.

кроме того, еще один момент: для быстрого обучения больших сетей логистическая сигмоида, которую вы используете, довольно ужасна. См. раздел 4.4 из эффективная сеть с обратным распространением ошибки по Лекун и соавт и использовать что-то с нулевым центром (на самом деле, читать вся эта статья, она очень полезна).

F# имеет лучшую производительность, чем C# в математических алгоритмах .NET. таким образом, переписывание нейронной сети В F# может улучшить общую производительность.

Если мы повторно реализовать LUT benchmarking snippet (Я использую слегка измененную версию) в F#, затем полученный код:

• выполняет тест sigmoid1 в 588,8 МС вместо 3899,2 МС
• выполняет тест sigmoid2 (LUT) в 156,6 МС вместо 411.4 ms

более подробную информацию можно найти в блоге. Вот фрагмент F# JIC:

#light

let Scale = 320.0f;
let Resolution = 2047;

let Min = -single(Resolution)/Scale;
let Max = single(Resolution)/Scale;

let range step a b =
  let count = int((b-a)/step);
  seq { for i in 0 .. count -> single(i)*step + a };

let lut = [| 
  for x in 0 .. Resolution ->
    single(1.0/(1.0 +  exp(-double(x)/double(Scale))))
  |]

let sigmoid1 value = 1.0f/(1.0f + exp(-value));

let sigmoid2 v = 
  if (v <= Min) then 0.0f;
  elif (v>= Max) then 1.0f;
  else
    let f = v * Scale;
    if (v>0.0f) then lut.[int (f + 0.5f)]
    else 1.0f - lut.[int(0.5f - f)];

let getError f = 
  let test = range 0.00001f -10.0f 10.0f;
  let errors = seq { 
    for v in test -> 
      abs(sigmoid1(single(v)) - f(single(v)))
  }
  Seq.max errors;

open System.Diagnostics;

let test f = 
  let sw = Stopwatch.StartNew(); 
  let mutable m = 0.0f;
  let result = 
    for t in 1 .. 10 do
      for x in 1 .. 1000000 do
        m <- f(single(x)/100000.0f-5.0f);
  sw.Elapsed.TotalMilliseconds;

printf "Max deviation is %f\n" (getError sigmoid2)
printf "10^7 iterations using sigmoid1: %f ms\n" (test sigmoid1)
printf "10^7 iterations using sigmoid2: %f ms\n" (test sigmoid2)

let c = System.Console.ReadKey(true);

и вывод (компиляция выпуска против F# 1.9.6.2 CTP без отладчика):

Max deviation is 0.001664
10^7 iterations using sigmoid1: 588.843700 ms
10^7 iterations using sigmoid2: 156.626700 ms

обновление: обновлен бенчмаркинг для использования 10^7 итераций, чтобы сделать результаты сопоставимыми с C

обновление 2: вот результаты работы Си из той же машина для сравнения с:

Max deviation is 0.001664
10^7 iterations using sigmoid1: 628 ms
10^7 iterations using sigmoid2: 157 ms

Примечание: Это продолжение этой пост.

Edit: обновление для расчета то же самое, что и этой и этой, принимая вдохновение от этой.

первая мысль: Как насчет некоторой статистики по переменной values?

• значения "value" обычно малы -10

Если нет, то вы, вероятно, можете получить импульс путем тестирования для значений за пределами границ

if(value < -10)  return 0;
if(value > 10)  return 1;

• часто ли повторяются значения?

Если это так, вы, вероятно, можете получить некоторую выгоду от Memoization (вероятно, нет, но это не помешает проверять....)

if(sigmoidCache.containsKey(value)) return sigmoidCache.get(value);

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

сопрано было несколько хороших оптимизаций ваш звонок:

public static float Sigmoid(double value) 
{
    float k = Math.Exp(value);
    return k / (1.0f + k);
}

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

например, попробуйте кэшировать последнее значение и результат. Если следующий вызов имеет то же значение, что и предыдущий, вам не нужно вычислять его, поскольку вы бы кэшировали последний результат. Если текущий вызов был таким же, как предыдущий вызов даже 1 из 100 раз, вы могли бы потенциально сэкономить 1 миллион вычислений.

или вы можете обнаружить, что в течение 10 последовательных вызовов параметр value в среднем одинаковый 2 раза, поэтому вы можете попробовать кэшировать последние 10 значений/ответов.

идея: возможно, вы можете сделать (большую) таблицу поиска с предварительно рассчитанными значениями?

Это немного не по теме, но просто из любопытства, я сделал ту же реализацию, что и в C,C# и F# в Java. Я просто оставлю это здесь на случай, если кто-то еще заинтересуется.

результат:

$ javac LUTTest.java && java LUTTest
Max deviation is 0.001664
10^7 iterations using sigmoid1() took 1398 ms
10^7 iterations using sigmoid2() took 177 ms

Я полагаю, что улучшение по сравнению с C# в моем случае связано с тем, что Java лучше оптимизирована, чем Mono для OS X. На аналогичной MS .NET-реализации (против Java 6, если кто-то хочет опубликовать сравнительные числа) я полагаю результаты будут другими.

код:

public class LUTTest {
    private static final float SCALE = 320.0f;
    private static final  int RESOLUTION = 2047;
    private static final  float MIN = -RESOLUTION / SCALE;
    private static final  float MAX = RESOLUTION / SCALE;

    private static final float[] lut = initLUT();

    private static float[] initLUT() {
        float[] lut = new float[RESOLUTION + 1];

        for (int i = 0; i < RESOLUTION + 1; i++) {
            lut[i] = (float)(1.0 / (1.0 + Math.exp(-i / SCALE)));
        }
        return lut;
    }

    public static float sigmoid1(double value) {
        return (float) (1.0 / (1.0 + Math.exp(-value)));
    }

    public static float sigmoid2(float value) {
        if (value <= MIN) return 0.0f;
        if (value >= MAX) return 1.0f;
        if (value >= 0) return lut[(int)(value * SCALE + 0.5f)];
        return 1.0f - lut[(int)(-value * SCALE + 0.5f)];
    }

    public static float error(float v0, float v1) {
        return Math.abs(v1 - v0);
    }

    public static float testError() {
        float emax = 0.0f;
        for (float x = -10.0f; x < 10.0f; x+= 0.00001f) {
            float v0 = sigmoid1(x);
            float v1 = sigmoid2(x);
            float e = error(v0, v1);
            if (e > emax) emax = e;
        }
        return emax;
    }

    public static long sigmoid1Perf() {
        float y = 0.0f;
        long t0 = System.currentTimeMillis();
        for (int i = 0; i < 10; i++) {
            for (float x = -5.0f; x < 5.0f; x+= 0.00001f) {
                y = sigmoid1(x);
            }
        }
        long t1 = System.currentTimeMillis();
        System.out.printf("",y);
        return t1 - t0;
    }    

    public static long sigmoid2Perf() {
        float y = 0.0f;
        long t0 = System.currentTimeMillis();
        for (int i = 0; i < 10; i++) {
            for (float x = -5.0f; x < 5.0f; x+= 0.00001f) {
                y = sigmoid2(x);
            }
        }
        long t1 = System.currentTimeMillis();
        System.out.printf("",y);
        return t1 - t0;
    }    

    public static void main(String[] args) {

        System.out.printf("Max deviation is %f\n", testError());
        System.out.printf("10^7 iterations using sigmoid1() took %d ms\n", sigmoid1Perf());
        System.out.printf("10^7 iterations using sigmoid2() took %d ms\n", sigmoid2Perf());
    }
}

Я понимаю, что прошел год с тех пор, как этот вопрос появился, но я столкнулся с ним из-за обсуждения производительности F# и C относительно C#. Я играл с некоторыми образцами из других респондентов и обнаружил, что делегаты, похоже, выполняются быстрее, чем обычный вызов метода, но нет явного преимущества peformance для F# над C#.

• C: 166ms
• C# (делегат): 275ms
• на C# (метод): 431ms
• C# (метод, счетчик с плавающей точкой): 2,656 МС
• F#: 404ms

C# с поплавковым счетчиком был прямым портом кода C. Гораздо быстрее использовать int в цикле for.

вы также можете рассмотреть возможность экспериментировать с альтернативными функциями активации, которые дешевле оценить. Например:

f(x) = (3x - x**3)/2

(который может быть учтена как

f(x) = x*(3 - x*x)/2

для одного умножения меньше). Эта функция имеет нечетную симметрию, а ее производная тривиальна. Использование его для нейронной сети требует нормализации суммы входов путем деления на общее количество входов (ограничение области до [-1..1], который также является диапазоном).

мягкая вариация на тему Сопрано:

public static float Sigmoid(double value) {
    float v = value;
    float k = Math.Exp(v);
    return k / (1.0f + k);
}

Так как вы только после одного результата точности, зачем делать математику.Exp функция вычислить двойной? Любой калькулятор экспонент, который использует итерационное суммирование (см. расширение e^x) займет больше времени для большей точности, каждый раз. А двойной-это в два раза больше работы одиночного! Таким образом, вы сначала конвертируете в сингл,затем сделать свой экспоненциальный.

но функция expf должна быть еще быстрее. Я не вижу необходимости в передаче soprano (float) в expf, хотя, если C# не выполняет неявное преобразование float-double.

в противном случае, просто использовать реальные язык, как FORTRAN...

есть много хороших ответов здесь. Я бы предложил запустить его через эта техника, просто чтобы убедиться, что

• вы не называете его больше раз, чем вам нужно.
  (Иногда функции вызываются больше, чем необходимо, просто потому, что их так легко вызвать.)
• вы не вызываете его повторно с теми же аргументами
  (где вы могли бы использовать мемоизации)

кстати функция у вас обратная логит-функция,
или обратное логарифмической функции отношения шансов log(f/(1-f)).

(обновлено с измерениями производительности) (Обновлено снова с реальными результатами :)

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

следующий фрагмент является примером реализации в C (я не говорю на c# достаточно свободно, чтобы высушить его). Он работает и работает достаточно хорошо, но я уверен, что в нем есть ошибка:)

#include <math.h>
#include <stdio.h>
#include <time.h>

#define SCALE 320.0f
#define RESOLUTION 2047
#define MIN -RESOLUTION / SCALE
#define MAX RESOLUTION / SCALE

static float sigmoid_lut[RESOLUTION + 1];

void init_sigmoid_lut(void) {
    int i;    
    for (i = 0; i < RESOLUTION + 1; i++) {
        sigmoid_lut[i] =  (1.0 / (1.0 + exp(-i / SCALE)));
    }
}

static float sigmoid1(const float value) {
    return (1.0f / (1.0f + expf(-value)));
}

static float sigmoid2(const float value) {
    if (value <= MIN) return 0.0f;
    if (value >= MAX) return 1.0f;
    if (value >= 0) return sigmoid_lut[(int)(value * SCALE + 0.5f)];
    return 1.0f-sigmoid_lut[(int)(-value * SCALE + 0.5f)];
}

float test_error() {
    float x;
    float emax = 0.0;

    for (x = -10.0f; x < 10.0f; x+=0.00001f) {
        float v0 = sigmoid1(x);
        float v1 = sigmoid2(x);
        float error = fabsf(v1 - v0);
        if (error > emax) { emax = error; }
    } 
    return emax;
}

int sigmoid1_perf() {
    clock_t t0, t1;
    int i;
    float x, y = 0.0f;

    t0 = clock();
    for (i = 0; i < 10; i++) {
        for (x = -5.0f; x <= 5.0f; x+=0.00001f) {
            y = sigmoid1(x);
        }
    }
    t1 = clock();
    printf("", y); /* To avoid sigmoidX() calls being optimized away */
    return (t1 - t0) / (CLOCKS_PER_SEC / 1000);
}

int sigmoid2_perf() {
    clock_t t0, t1;
    int i;
    float x, y = 0.0f;
    t0 = clock();
    for (i = 0; i < 10; i++) {
        for (x = -5.0f; x <= 5.0f; x+=0.00001f) {
            y = sigmoid2(x);
        }
    }
    t1 = clock();
    printf("", y); /* To avoid sigmoidX() calls being optimized away */
    return (t1 - t0) / (CLOCKS_PER_SEC / 1000);
}

int main(void) {
    init_sigmoid_lut();
    printf("Max deviation is %0.6f\n", test_error());
    printf("10^7 iterations using sigmoid1: %d ms\n", sigmoid1_perf());
    printf("10^7 iterations using sigmoid2: %d ms\n", sigmoid2_perf());

    return 0;
}

предыдущие результаты были связаны с оптимизатором делает свою работу и оптимизировать расчеты. Выполнение кода на самом деле дает немного другие и гораздо более интересные результаты (на моем пути медленный MB Air):

$ gcc -O2 test.c -o test && ./test
Max deviation is 0.001664
10^7 iterations using sigmoid1: 571 ms
10^7 iterations using sigmoid2: 113 ms

---

TODO:

есть вещи, чтобы улучшить и способы устранения недостатков; как это сделать, остается в качестве упражнения для читателя:)

• настройка диапазона функции, чтобы избежать скачка, где начинается таблица и концы.
• добавьте функцию небольшого шума, чтобы скрыть артефакты сглаживания.
• как сказал Рекс, интерполяция может дать вам совсем немного больше точности, будучи довольно дешевой с точки зрения производительности.

есть гораздо более быстрые функции, которые делают очень похожие вещи:

x / (1 + abs(x)) - быстрая замена для TAHN

а также:

x / (2 + 2 * abs(x)) + 0.5 быстрая замена сигмовидной

сравнить графики с фактическим сигмоидом

делая поиск в Google, я нашел альтернативную реализацию сигмоидной функции.

public double Sigmoid(double x)
{
   return 2 / (1 + Math.Exp(-2 * x)) - 1;
}

это правильно для ваших нужд? Это быстрее?

http://dynamicnotions.blogspot.com/2008/09/sigmoid-function-in-c.html

1) вы называете это только из одного места? Если это так, вы можете получить небольшую производительность, переместив код из этой функции и просто поместив его туда, где вы обычно вызывали Сигмоидную функцию. Мне не нравится эта идея с точки зрения удобочитаемости кода и организации, но когда вам нужно получить каждый последний прирост производительности, это может помочь, потому что я думаю, что вызовы функций требуют push/pop регистров в стеке, чего можно было бы избежать, если бы ваш код был весь линейный.

2) я понятия не имею, может ли это помочь, но попробуйте сделать ваш параметр функции параметром ref. Смотрите, если это быстрее. Я бы предложил сделать его const (что было бы оптимизацией, если бы это было в c++), но c# не поддерживает параметры const.

Если вам нужен гигантский прирост скорости, вы, вероятно, можете посмотреть на распараллеливание функции с помощью силы (ge). IOW, используйте DirectX для управления видеокартой, чтобы сделать это за вас. Я понятия не имею, как это сделать, но я видел, как люди используют графические карты для всех видов вычислений.

Я видел, что многие люди здесь пытаются использовать приближение, чтобы сделать сигмоид быстрее. Однако важно знать, что сигмоид также может быть выражен с помощью tanh, а не только exp. Вычисление сигмоида таким образом примерно в 5 раз быстрее, чем с экспоненциальным, и с помощью этого метода вы ничего не приближаете, поэтому исходное поведение сигмоида сохраняется как есть.

    public static double Sigmoid(double value)
    {
        return 0.5d + 0.5d * Math.Tanh(value/2);
    }

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