Avoiding Possible Precision Loss with a Simple Moving Average

Question

Suppose we had a basic moving average function that was keeping track of the sum. For example:

Queue values;
double sum;

double CalcSMA(double next) {
    values.push(next);
    sum -= values.pop();
    sum += next;
    return sum / SMA_LENGTH;
}

One example of how this could break down would be if our window was 5 wide that we fed it something like: 2, 2, 2, 2, 2, 2, 2, 1E100, 2, 2, 2, 2, 2, 2, 2, 2. The output would then be 2, 2, 2, 2E99, 2E99, 2E99, 2E99, 2E99, 0, 0, 0.

Even if the sum isn't quite that dramatically off, there still could be quite reasonable instances where a small loss in precision could make the sum artificially increase by a tiny amount. Over a long period of time, this could add up and become an issue.

Does anyone have any ideas of how to work around the loss in precision?

EDIT: note that this function is designed to work O(1). That is necessary. So, recalculating each time won't work: the window is too large.

Answer 1

You could recalculate a fresh sum over every SMA_LENGTH values to stop the errors accumulating:

Queue values;
double sum = 0.0;
double freshsum = 0.0;
int count = 0;

double CalcSMA(double next) {
    values.push(next);
    sum -= values.pop();
    sum += next;
    freshsum += next;
    if(++count == SMA_LENGTH)
    {
        sum = freshsum;
        freshsum = 0.0;
        count = 0;
    } 
    return sum / SMA_LENGTH;
}

Answer 2

What samgak proposed doesn't actually guarantee good averages if you constantly provide it with evil values.

You can use Neumaier's algorithm to generate accurate results in O(1) time. Something like this:

const double SMA_LENGTH = 5;
Queue values;
double sum = 0.0;
double correction = 0.0;

static void Neumaier(double value, ref double sum, ref double correction)
{
    var t = sum + value;
    if (Math.Abs(sum) >= Math.Abs(value))
        correction += (sum - t) + value;
    else
        correction += (value - t) + sum;
    sum = t;
}

double CalcSMA(double next)
{
    Neumaier(-values.pop(), ref sum, ref correction);
    Neumaier(next, ref sum, ref correction);
    values.push(next);
    return (sum + correction) / SMA_LENGTH;
}

If you have huge sequences, you can reset the window every 10^15 or so additions. This is because, with double precision, the algorithm starts losing accuracy after about 10^16 additions.

On the other hand, Neumaier is more complex so it is probably slower.