Как эффективно рассчитать движущееся стандартное отклонение

Ответ 1

Ответ: да, вы можете. В середине 80-х я разработал именно такой алгоритм (возможно, не оригинальный) в FORTRAN для приложения мониторинга и управления процессом. К сожалению, это было более 25 лет назад, и я не помню точных формул, но эта техника была продолжением той, что была для скользящих средних, с расчетами второго порядка, а не только линейными.

Посмотрев на свой код, я думаю, что могу судить, как я это делал тогда. Обратите внимание, как ваш внутренний цикл создает сумму квадратов?:

            for (int x = i; x > (i - period); x--)
            {
                total_bollinger += Math.Pow(data.Values[x]["close"] - average, 2);
            }

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

SUM(i=1; n){ (v[i] - k)^2 }

является

SUM(i=1..n){v[i]^2 -2*v[i]*k + k^2}

который становится

SUM(i=1..n){v[i]^2 -2*v[i]*k} + k^2*n

который

SUM(i=1..n){v[i]^2} + SUM(i=1..n){-2*v[i]*k} + k^2*n

который также является

SUM(i=1..n){v[i]^2} + SUM(i=1..n){-2*v[i]}*k + k^2*n

Теперь первый термин - это всего лишь сумма квадратов, вы обрабатываете это так же, как и сумму значений для среднего. Последний член (k^2*n) - это только среднее квадратическое число period. Поскольку вы все равно разделите результат на период, вы можете просто добавить новый средний квадрат без дополнительного цикла.

Наконец, во втором члене (SUM(-2*v[i]) * k), так как SUM(v[i]) = total = k*n вы можете изменить его на это:

-2 * k * k * n

или просто -2*k^2*n, что в 2 раза превышает среднюю квадрат, как только период (n) будет снова разделен. Таким образом, окончательная комбинированная формула:

SUM(i=1..n){v[i]^2} - n*k^2

или

SUM(i=1..n){values[i]^2} - period*(average^2)

(обязательно проверьте правильность этого, так как я вывожу его с головы)

И включение в ваш код должно выглядеть примерно так:

public static void AddBollingerBands(ref SortedList<DateTime, Dictionary<string, double>> data, int period, int factor)
{
    double total_average = 0;
    double total_squares = 0;

    for (int i = 0; i < data.Count(); i++)
    {
        total_average += data.Values[i]["close"];
        total_squares += Math.Pow(data.Values[i]["close"], 2);

        if (i >= period - 1)
        {
            double total_bollinger = 0;
            double average = total_average / period;

            double stdev = Math.Sqrt((total_squares - Math.Pow(total_average,2)/period) / period);
            data.Values[i]["bollinger_average"] = average;
            data.Values[i]["bollinger_top"] = average + factor * stdev;
            data.Values[i]["bollinger_bottom"] = average - factor * stdev;

            total_average -= data.Values[i - period + 1]["close"];
            total_squares -= Math.Pow(data.Values[i - period + 1]["close"], 2);
        }
    }
}

Ответ 2

Проблема с подходами, которые вычисляют сумму квадратов, состоит в том, что она и квадрат сумм могут стать довольно большими, а вычисление их разности может ввести very большая ошибка, так что подумайте о чем-то лучше. Для чего это необходимо, см. Статью Википедии о Алгоритмы вычисления дисперсии и Джон Кук на Теоретическое объяснение численных результатов)

Во-первых, вместо вычисления stddev пусть сосредоточится на дисперсии. Как только у нас есть дисперсия, stddev - это просто квадратный корень дисперсии.

Предположим, что данные находятся в массиве с именем x; сканирование окна размером n можно рассматривать как удаление значения x[0] и добавление значения x[n]. Обозначим средние значения x[0]..x[n-1] и x[1]..x[n] через μ и μ соответственно. Разница между дисперсиями x[0]..x[n-1] и x[1]..x[n] заключается в том, что после отмены некоторых условий и применения (a²-b²) = (a+b)(a-b):

Var[x[1],..,x[n]] - Var[x[0],..,x[n-1]] 
= (\sum_1^n x[i]² - n µ’²)/(n-1) - (\sum_0^{n-1} x[i]² - n µ²)/(n-1)
= (x[n]² - x[0]² - n(µ’² - µ²))/(n-1) 
= (x[n]-µ’ + x[0]-µ)(x[n]-x[0])/(n-1)

Поэтому дисперсия возмущена тем, что не требует от вас поддержки суммы квадратов, что лучше для численной точности.

Вы можете рассчитать среднее значение и дисперсию один раз в начале с помощью правильного алгоритма (метод Welford). После этого каждый раз, когда вам нужно заменить значение в окне x[0] другим x[n], вы обновляете среднее значение и дисперсию следующим образом:

new_Avg = Avg + (x[n]-x[0])/n
new_Var = Var + (x[n]-new_Avg + x[0]-Avg)(x[n] - x[0])/(n-1)
new_StdDev = sqrt(new_Var)

Ответ 3

Я использовал commons-math (и внес вклад в эту библиотеку!) для чего-то очень похожего на это. Это с открытым исходным кодом, портирование на С# должно быть легким, как купил на складе пирог (вы пытались сделать пирог с нуля!?). Проверьте это: http://commons.apache.org/math/api-3.1.1/index.html. У них есть класс StandardDeviation. Пойдите в город!

Ответ 4

Самая важная информация уже была приведена выше, но, возможно, это все еще представляет общий интерес.

Небольшая библиотека Java для расчета скользящей средней и стандартного отклонения доступна здесь: https://github.com/tools4j/meanvar

Реализация основана на варианте метода Welford, упомянутом выше. Выведены методы удаления и замены значений, которые можно использовать для перемещения значений окон.