Причинный тест Грейнджера на Java

В своей прошлой заметке о причинном тесте Грейнджера я рассказал о его применении к данным о производстве яиц и количестве куриц-несушек в США. Я не рассказал, однако, как же этот тест провести самостоятельно и что значат все эти числа в таблицах. Рассказываю.

Вероятностный тест Грейнджера поддерживается многими статистическими пакетами, такими как Stata или R, но гораздо интересней написать реализацию самому. Заодно это поможет разобраться в тонкостях и деталях. Будем реализовывать его на Java. Для этого нам понадобиться Apache Commons Math в ней мы найдем реализации F-статистики и линейной регрессии.



Итак, у нас есть две переменные, Yи X изменяющиеся во времени и мы хотим оценить вероятность того, что X является причиной по-Грейнджеру для Y с временным лагом L. То есть, если, например L = 1, то вчерашние значения X должны предсказывать значения Y.

Выдвигаем две гипотезы, нулевая гипотеза H₀ что Y не зависит от прошлых значений X. И альтернативная гипотеза H₁, что не зависит. Для проверки этих гипотез, нам следует провести две линейные регрессии, и определить предполагаемые коэффициенты μ_i, α_k, β_k и ошибку ε_i в выражениях:

H₀: Y_i = μ_i + ∑α_kY_i-l + ε_i (где k = 1..L) (1)

H₁: Y_i = μ_i + ∑α_kY_i-l + ∑β_kX_i-l + ε_i (где k = 1..L) (2)


Линейная регрессия методом наименьших квадратов уже реализована в одном из классов Apache Commons Math org.apache.commons.math.stat.regression.OLSMultipleLinearRegression, чем мы и воспользуемся. Следует только помнить, что апачевская реализация считает регрессию для выражения:
Y=X*b+u
нам же, обязательно нужно учесть свободный член μ, для этого мы включим его в X, равно как и прошлые значения Y для регрессии (2).

Итак, код:

//Для нулевой гипотезы нам необходимы сдвинутые по времени значения y
        double[][] laggedY = createLagged(L, y);
//Для альтернативной гипотезы нам необходимы сдвинутые по времени значения y и x
        double[][] laggedXY = createLagged(L, x, y);

//Определяем две линейные регрессии
        OLSMultipleLinearRegression h0 = new OLSMultipleLinearRegression();
        OLSMultipleLinearRegression h1 = new OLSMultipleLinearRegression();
        h0.newSampleData(y, laggedY);
        h1.newSampleData(y, laggedXY);

Теперь когда все сделано, наступает самая сложная часть, понять результаты.
Для начала, давайте оценим ошибки в каждой из регрессий. Что такое ошибка? Ошибка это последний член выражений (1) и (2): ε_i
И оценить её можно весьма популярным в статистике и теории вероятности способом, суммой квадратов:
RSS = ∑ε_i²
В нашей Java версии это выглядит так:

double rs0[] = h0.estimateResiduals();
double RSS0 = sqrSum(rs0);

Сама по себе величина RSS мало о чем говорит, ведь мы же не знаем в каком диапазоне колеблется значение оцениваемой нами величины Y. Поэтому, оценим величину отклонения Y от среднего значения:
TSS = ∑(Y_i - Ymean)² , где Ymean - среднее значение Y

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

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



Получив это значение, статистики, обычно смотрят специальные таблицы, на основании которых делают выводы о применимости гипотез. Но настоящие программисты должны быть достаточно ленивы чтобы им было лень листать таблицы в поисках соответствующего значения. Тем более, что F-распределение реализовано в Apache Commons Math и мы можем просто вычислить статистическую значимость (P-value) для полученного нами значения:

FDistribution fDist = new FDistributionImpl(L, n-2*L-1);
        double pValue = 1.0 - fDist.cumulativeProbability(ftest);

Как понимать это p-value? Многие полагают, что p-value это вероятность того, что нулевая гипотеза H₀ - верна. На самом деле это не совсем так. Но для простоты, мы можем это принять. Действительно, значения p-value близкие к нулю говорят нам о том, что нулевую гипотезу следует отвергнуть и стало быть, принять альтернативную гипотезу, то есть в этом случае мы можем утверждать что значения X помогают предсказать Y.

Кроме того, нам следует проверить, насколько вообще регрессия (2) адекватно описывает изменения величины Y. Для этого вводится величина

Эта величина принимает значения 0 ≤ R² ≤ 1 и чем ближе это значение к 1 тем лучше наша регрессионаая модель описывает поведение величины Y.

И, на последок, замечание о финансах, у теста Грейнджера есть одно важное ограничение. Он работает только на стационарных данных, т.е. данных, распределение которых не меняется во времени. Котировки на бирже к таковым не относятся. А вот изменение котировок, например разница цен закрытия за два дня - относится. Поэтому при применении теста Грейнджера к рыночным ценам, необходимо сначала приобразовать данные.

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