Поддержание точности с плавающей запятой с помощью скользящего среднего - PullRequest
Новая
       2

Поддержание точности с плавающей запятой с помощью скользящего среднего

13 голосов
/ 27 января 2011

Мне нужно рассчитать среднеквадратичную ошибку 16-битной операции для произвольного числа точек данных (свыше 100 миллионов).Я решил использовать скользящее среднее, чтобы мне не пришлось беспокоиться о переполнении при добавлении большого количества квадратов ошибок.При 100 миллионах выборок у меня были проблемы с точностью с плавающей запятой (неточные результаты), поэтому я перешел на удвоение.

Вот мой код 

int iDifference = getIdeal() - getValue();

m_iCycles++;


// calculate the running MSE as

// http://en.wikipedia.org/wiki/Moving_average

// MSE(i + 1) = MSE(i) + (E^2 - MSE(i))/(i + 1)

m_dMSE = m_dMSE + ((pow((double)iDifference,2) - m_dMSE) / (double)m_iCycles);

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

Ответы [ 3 ]

13 голосов
/ 27 января 2011

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

6 голосов
/ 27 января 2011

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

Чтобы сохранить точность в промежуточном итоге, сохраните промежуточные итоги вместо единого итога.Просто продолжайте добавлять к промежуточному итогу, пока добавление еще одного не вызовет переполнение.Затем перейдите к следующему промежуточному итогу.Поскольку все они имеют одинаковый порядок величины (в базе 2), оптимальная точность может быть достигнута путем преобразования в плавающую точку и использования попарного накопления в одну итоговую сумму. 

// first = errors, second = counter
typedef pair< vector< uint32_t >, uint32_t > running_subtotals;

void accumulate_error( uint32_t error, running_subtotals &acc ) {
    ( numeric_limits< uint32_t >::max() - error < acc.first.back()?
        * acc.first.insert( acc.first.end(), 0 ) : acc.first.back() )
        += error; // add error to current subtotal, or new one if needed
    ++ acc.second; // increment counter
}

double get_average_error( running_subtotals const &total ) {
    vector< double > acc( total.first.begin(), total.first.end() );
    while ( acc.size() != 1 ) {
        if ( acc.size() % 2 ) acc.push_back( 0 );
        for ( size_t index = 0; index < acc.size() / 2; ++ index ) {
            acc[ index ] = acc[ index * 2 ] + acc[ index * 2 + 1 ];
        }
        acc.resize( acc.size() / 2 );
    }
    return acc.front() / total.second;
}
2 голосов
/ 27 января 2011

Если другие ваши решения не работают, вы можете исследовать библиотеку Bignum

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

...