Тест Бартельса на значимость цикла: секрет утерянной жемчужины
Статья Сергея Тарасова
Когда мы говорим о циклическом анализе или гармоническом анализе, обычно первое, что приходит нам на ум, - это преобразование Фурье. Применение этого метода путем вычисления спектрограммы (периодограммы) было предложено и разработано известным британским физиком сэром Артуром Шустером. Артур Шустер рассчитал первую спектрограмму в 1899 году. С тех пор этот подход получил большое развитие. Многие книги и монографии основаны на этой методике, и она является обязательной частью университетских курсов. Великие идеи Артура Шустера легли в основу современного гармонического анализа.
Артур Шустер и Юлиус Бартельс
Есть еще одно великое имя, которое заслуживает признания. Это немецкий геофизик Юлиус Бартельс. Юлиус Бартельс внес большой вклад в геомагнетизм; он изучал влияние на него Луны и Солнца. Он разработал альтернативный подход к циклическому анализу, который сейчас известен как критерий значимости Бартельса (Bartels Significance Test).
Следуя стандартам научных исследований своего времени, Юлиус Бартельс должен был применить некоторые методы, которые могли бы обеспечить некоторую статистическую оценку анализируемых циклов. Например, ученые, которые исследуют влияние экваториального вращения Солнца (24-29 дней для разных широт) на геомагнитную активность, - им нужно получить некоторую статистическую значимость, что-то помимо корреляции Пирсона. До Бартельса физики использовали подход Артура Шустера, сосредоточившись на поиске циклов. Бартельс наблюдал циклы, которые уже были установлены. В этом случае классический подход (периодограмма) не слишком помогает. Это была горячая тема и предмет обсуждения (я бы даже сказал, "напряженности") между физиками и статистиками в 1930-1940 годах. Джулиус Бартельс предложил свой подход при изучении геомагнитной активности. Тот же подход был использован для исследований в области экономики (Э. Дьюи). Первоначальная идея была действительно блестящей. Однако, по какой-то причине, сегодня, 90 лет спустя, мы все еще мало что знаем о тесте значимости Бартельса; это определенно не основное направление в современном гармоническом анализе. Что же случилось? Почему это произошло?
Похоже, что в оригинальной публикации теста значимости Бартельса была допущена ошибка; мы обсудим это позже. Из-за этой ошибки вероятность, рассчитанная этим методом, кажется нереалистично высокой: она равна Exp(N ^ 2); выше, чем правильно рассчитанная вероятность (N- размер выборки). В результате сейчас, 90 лет спустя, этот подход не используется так так широко, как следовало бы.
Моя личная история взаимоотношений с тестом значимости Бартельса - хорошая иллюстрация того, что мы упустили. Я слышал о тесте Бартельса практически с самого начала проекта Timing Solution, 17 лет назад. В течение всех этих лет пользователи программного обеспечения Timing Solution время от времени просили меня включить его в программу (в основном, они узнали об этом из деятельности Фонда изучения циклов). Этот вопрос появлялся практически каждый месяц в течение этих лет (это своего рода мой личный день сурка). Я не решался попробовать. Мои причины? Одна из них заключалась в том, что, будучи физиком в прошлом, я не слышал об этом, хотя много знал о периодограммах; добавьте к этому тот факт, что существует очень мало научных публикаций на эту тему. Другая причина заключалась в том, что тогда мое основное внимание было сосредоточено на поиске циклов, которые могли бы присутствовать в моих наборах данных. Я делал это путем рассчета классических периодограмм (метод Шустера). Периодограмма, которую я получил, пытаясь воспроизвести подход Бартельса, выглядела на тот момент очень неуверенно (почему - это будет объяснено далее). Итак, на первый взгляд, это выглядело как некое переосмысление гармонического анализа (который уже был одним из наиболее развитых и востребованных разделов математики).
Наконец, в этом году я решил положить конец этой истории. Я изучил доступные открытые источники и провел детальный анализ теста значимости Бартельса. Вот хроника моей битвы: https://www.timingsolution.com/TS/Articles/Bartels/index.htm .
И знаете что? Вдаваясь в подробности, я обнаружил нечто совершенно неожиданное для себя:
Я утверждаю, что в исходной формуле теста значимости Бартельса допущена ошибка, которая делает вероятность, рассчитанную предложенным способом, нереалистично высокой. Исправив эту ошибку, мы получаем значительно улучшенную периодограмму. С точки зрения качества (отношение сигнал/шум) периодограмма, рассчитанная по скорректированному методу Бартельса, является лучшей из тех, что я когда-либо наблюдал ранее.
Посмотрите сами, как это работает по сравнению с классическим подходом:
Взгляните на эти две периодограммы, рассчитанные для моего тестового набора данных при наложении синусоидальных волн в 35, 76 и 125 дней:
На увеличенном графике видно, что периодограмма Бартельса в 5-6 раз менее шумная, чем классическая периодограмма Шустера:
Что было потеряно
Вот еще примеры применения этих двух подходов к физическим и экономическим задачам (показаны модифицированные периодограммы Бартельса и классические периодограммы периода времени):
1. Это спектральный анализ индекса активности солнечных пятен с 1818 года по настоящее время:
Как видите, периодограмма Бартельса выглядит гораздо более чистой от шума, чем классическая периодограмма Шустера.
2. Давайте попытаемся выявить основные экономические циклы, которые присутствуют в промышленном индексе Доу-Джонса с 1885 года по настоящее время, принимая во внимание 135-летнюю историю цен:
Обе периодограммы выявили основные экономические циклы - 40-месячный цикл Китчина и 10-летний цикл Юглара, хотя, как вы видите, периодограмма Бартельса выглядит более определенной.
3. Преимущество периодограммы Бартельса особенно очевидно, когда она применяется к большим массивам данных. Вот две периодограммы, рассчитанные для глобального индекса геомагнитной активности Kp3H - 3 часа за 1932-1970 годы (в общей сложности более 100 тысяч пунктов).:
Разница очевидна.
Во всех приведенных выше примерах применялся модифицированный мною тест Бартельса, а не оригинальный. Оригинальный алгоритм содержит ошибку, которая затрудняет его применение к практическим задачам.
Чтобы увидеть, насколько велика эта разница, посмотрите на две периодограммы для индекса Kp3H, рассчитанные по оригинальному алгоритму Бартельса и по модифицированному:
Как видим, периодограмма, основанная на оригинальном алгоритме Бартельса, показывает, что практически любой цикл имеет статистическую значимость 50-60%. Это невозможно.
Модифицировав оригинальный алгоритм Бартельса, мы нашли кое-что действительно хорошее и полезное. Похоже, что эта печальная ошибка закрыла дверь в будущее для метода Бартельса. Подход Бартельса заслуживает гораздо большего. Это реальная альтернатива классическому подходу к гармоническому анализу. Я просто хочу обратить ваше внимание на этот подход, он стоит вашего времени!
Объяснение идей Юлиуса Бартельса
Классический циклический анализ, или гармонический анализ, начинается с преобразования Фурье. Его применение к задачам современной науки, особенно физики, было разработано Артуром Шустером. Бартельс был знаком с работами Шустера (см. http://ui.adsabs.harvard.edu/abs/1934TeMAE..39..345B/abstract ) и предложил свой собственный подход. Сравнивая классический подход и то, что сделал Бартельс, я был впечатлен математической красотой подхода Бартельса. Это действительно дает нам что-то новое.
Чтобы продемонстрировать, как тест значимости Бартельса работает на практике, мы применим этот подход к глобальному геомагнитному индексу KpH3 для небольшого фрагмента данных за 1939-1940 годы (вполне возможно, что Дж. Бартельс работал с этим набором данных). Там мы попытаемся обнаружить наличие цикла вращения Солнечного экватора. Цикл хорошо известен и описан, его усредненный период составляет 28 дней.
Обратите внимание на различный подход к одной и той же задаче:
- следуя Шустеру, мы будем искать просто наличие любого цикла,
- в то время как Бартельс предполагает, что цикл существует (28 дней), и мы будем искать доказательства его влияния на анализируемые данные.
Чтобы сделать это, давайте разделим весь набор данных на фрагменты данных за 28 дней и назовем каждый фрагмент интервалом с числом 1, 2,... Набор данных охватывает примерно год, так что в общей сложности у нас есть 12 интервалов.
Для каждого из этих интервалов мы проводим стандартную процедуру Фурье, определяя амплитуду и фазу 28-дневной волны, скорректированные для каждого интервала. Другими словами, для каждого интервала из этих 12 мы вычисляем коэффициенты "a" и "b" для этой волны: X*cos(2 Pi*t/28)+Y*sin(2 Pi*t/28).
Следующий шаг - поместить вектор (X,Y) на шкалу гармоник или на гармонический циферблат:
Например, посмотрите на первый интервал в первые 28 дней. Он охватывает период времени с 18 декабря 1939 года в 19:30 вечера по 15 января 1940 года в 19:30 вечера. Применяя стандартную процедуру Фурье, мы получаем следующие коэффициенты для формулы 28-дневной волны:
-1.98683*Cos(вт*т) -4.82192*Sin(вт*т)
где w=2*Pi*t/28,1 t- время в днях, начиная с полуночи 19:00 по Гринвичу
Мы можем представить эту первую волну в виде вектора на шкале гармоник с координатами X (коэффициент для Cos) и Y (коэффициент для Sin).
Мы повторяем ту же процедуру для других оставшихся интервалов, от 2 до 11, находя коэффициенты a и b для каждого из них. После этого мы помещаем все эти векторы на шкалу гармоник:
На рисунке выше показаны 12 векторов (размер выборки =12); они представляют, как работает 28-дневный цикл на каждом из этих 12 интервалов, как амплитуда и фаза изменялись от одного интервала к другому. Если 28-дневный цикл окажет какое-то влияние на этот набор данных, все эти векторы на шкале гармоник сформируют некоторый кластер, некоторое облако. В идеале, если этот цикл работает одинаково для всех интервалов, все векторы будут совпадать. Чем ближе векторы расположены друг к другу, тем плотнее облако - тем лучше этот цикл работает для данного набора данных.
Эта блестящая идея Юлиуса Бартельса - одна из самых красивых математических идей, которые я знаю. Мы заменяем стандартную процедуру гармонического анализа кластерным анализом, предоставляя тем самым пространство для статистического анализа, статистического анализа, который не совсем хорошо подходил для подхода, предложенного А. Шустером (как это обсуждалось в 1930-1940-х годах). Это действительно привносит что-то новое!
Итак, теперь нам нужно решить статистическую задачу: найти статистическую меру этого облака. Чем плотнее облако (что можно наблюдать, если все эти векторы расположены ближе друг к другу), тем выше вероятность того, что 28-дневный цикл играет определенную роль в нашем наборе данных (не является "случайным явлением", как выражается Бартельс) и тем больше статистическая значимость этого цикла.
Чтобы сделать это, давайте проведем несколько простых вычислений с помощью гармонического циферблата Бартельса:
Мы вычисляем центр масс этих 12 векторов, точку A (см. ниже). Это представляет собой ожидаемую амплитуду и фазу нашего 28-дневного цикла.
Теперь давайте посмотрим, как эти 12 векторов распределены относительно центра масс A.
Мы вычисляем ожидаемый радиус, который является среднеквадратичным значением всех 12 векторов с центром в A, среднеквадратичным значением векторов A1, A2, .., A12. Маленький желтый круг на шкале гармоник - это ожидаемый радиус с центром в точке A. Чем меньше это облако и дальше оно от центра O, тем статистически значимее 28-дневный цикл.
Мы вычисляем соотношение (радиус центра масс) / (ожидаемый радиус) и применяем его к показателю Гаусса (предполагая нормальное распределение для этих радиусов):
Статистическая значимость 28-дневного цикла для индекса геомагнитной активности Ko3H составляет 93%.
Команда Timing Solution провела множество тестов на многих наборах данных; мы обнаружили, что вероятность, рассчитанная таким образом, дает очень четкие периодограммы.
В чем была ошибка?
По какой-то неизвестной причине во всех доступных мне источниках - начиная с работы Дж. Бартелса, опубликованной в 1932 году, и заканчивая статьей Чарльза Э. Армстронга, опубликованной Фондом изучения циклов в 1973 году, - расчет ожидаемого радиуса показан как среднеквадратичное значение векторов, центрированных по облаку, деленное на квадрат корень из размера выборки. В данном конкретном случае (пример выше) применение опубликованной формулы приводит к тому, что ожидаемый радиус меньше квадратного корня (12)=в 3,46 раза.
Если я произведу расчет в соответствии с опубликованным описанием, то гармонический циферблат для нашего примера будет выглядеть следующим образом, при этом ожидаемый радиус здесь нереально мал:
Я проверил каждую цифру в публикации Журнала Американской статистической ассоциации (Journal of the American Statistical Association, июнь 1940 г.), где обсуждалась технология Бартельса:
В качестве модели была использована идея броуновского движения по причине, объясненной Джулиусом Бартельсом как "низкий уровень распространения случайной ошибки".
Это еще один недостаток. Броуновское моделирование предпочтительнее, если имеется некоторое количество последовательных шагов: сделав N последовательных случайных шагов, наше ожидаемое расстояние от начальной точки будет равно . Таким образом, вероятность того, что мы, сделав N случайных шагов, достигнем точки, в которой мы находимся на расстоянии K шагов от начальной точки, равна:
По-моему, это не наш случай. Вместо броуновского движения, где каждый шаг приводит к росту кумулятивной ошибки, у нас есть набор данных, который строится из меньших независимых наборов данных, - так мы получаем кластер отдельных векторов.
Приведенные ниже две периодограммы рассчитаны с использованием оригинального и скорректированного алгоритмов Бартельса для набора тестовых данных с 35-дневной волной. Почувствуйте разницу:
Оригинальный алгоритм Бартельса дает неопределенный результат и бесполезную информацию: практически любой цикл в диапазоне 34-36 дней имеет 100%-ную значимость:
Я провел много тестов, пытался применить различные формулы, чтобы воспроизвести график в оригинальной статье; я просто понятия не имею, что это может быть, если не ошибка. Если я сделал здесь что-то не так, я буду признателен, если кто-нибудь объяснит мне это.
Статистическая значимость в сравнении с информационной значимостью
Чарльз Э. Армстронг в своей статье "Применение критерия значимости Бартельса к анализу циклов временных рядов" ("Applying the Bartels Test of Significance to a Time Series Cycle Analysis) в 1944 году отмечал: "Тест Бартельса разработан не в первую очередь как средство определения периодов циклов, присутствующих в ряду. Его главная ценность заключается в его применении в качестве критерия значимости после того, как период был определен каким-либо образом" (The Bartels' test is not designed primarily as a mean of locating the periods of cycles present in a series. Its chief value lies in its application as a test of significance after the period has been located by some means).
Давайте посмотрим на эту проблему с точки зрения физика. Сейчас 1944 год. В мире науки методы математической статистики были окончательно усовершенствованы и стали обычной практикой в любых исследованиях и процессе принятия решений. Наш физик провел кое-какие исследования и хотел знать, можно ли доверять 35-дневному циклу, найденному в его наборе данных. На тот момент были доступны две периодограммы: классическая периодограмма Шустера и оригинальная периодограмма Бартельса (не исправленная). Вот они:
Бартельс утверждает, что практически любой цикл в диапазоне 34-36 дней имеет 100% статистическую значимость. Это означает, что у нас есть из чего выбирать. Вот почему опубликованный тест Бартельса (оригинальный) не может быть использован в качестве инструмента для выбора циклов. С другой стороны, классическая периодограмма не предоставляет статистической информации, в то время как пик на ней точно указывает на 35-дневный цикл. Первый случай (т.е. подход Бартельса) обеспечивает статистическую значимость по сравнению со вторым, который дает информационную значимость. Статистика и здравый смысл здесь вступают в противоречие, и угадайте, что выбрал наш физик? На мой взгляд, именно по этой причине методика Бартельса не представлена в университетских курсах.
Применяя исправленный алгоритм Бартельса, мы избегаем этого противоречия, используя лучшие части этих подходов. Теперь обе периодограммы указывают на один и тот же 35-дневный цикл:
Скорректированная периодограмма Бартельса анализирует цикл со статистической точки зрения, в то время как классическая периодограмма Шустера делает это с точки зрения соответствия этого цикла анализируемым данным. Как я вижу, в большинстве случаев оба метода работают вместе.
Идеальная периодограмма
Тест Бартельса имеет свои собственные побочные эффекты - ложные пики, которые появляются при небольших размерах выборки, то есть на больших периодах:
В приведенном выше примере классическая периодограмма выглядит намного лучше. Периодограмма Бартельса, начиная с определенного периода, показывает множество ложных пиков, которые представляют несуществующие циклы. Это артефакт, связанный с резонансом используемых интервалов с ложным циклом. Этот побочный эффект проявляется при размерах выборки 15 и менее, т.е. когда доступный набор данных охватывает менее 15 полных периодов.
Я провел некоторое исследование этого побочного эффекта и нашел решение. Поскольку мы делим исходные данные на меньшие интервалы, похоже, что здесь возникает резонанс.
Построение гармонического шкалы для одного из этих ложных пиков - вот что было найдено:
Все векторы выровнены, они синхронизированы только по фазе.
Для решения этой проблемы был применен другой алгоритм.
Посмотрите теперь на эти две периодограммы: Бартельс, скорректированный с помощью нового алгоритма, чтобы избежать ложных пиков, и классическая периодограмма Шустера. Я увеличиваю изображение, чтобы показать, насколько периодограмма Бартельса четкая по сравнению с классической:
Я в жизни не видел такой чистой периодограммы. Я работал со многими версиями периодограмм, основанных на классическом подходе Шустера, включая более продвинутые (преобразование Фурье для ковариации, различное сглаживание и т.д.). Не было ничего даже близкого по сравнению с модифицированной периодограммой Бартельса. Для меня это было большим сюрпризом!
Послесловие
Эта статья была написана только с одной целью: донести сообщение о том, что периодограмма, основанная на скорректированном тесте значимости Бартельса, работает очень хорошо, значительно лучше, чем другие алгоритмы. Я чувствую себя обязанным сделать это, особенно после того, как игнорировал этот подход в течение 17 лет.
Я действительно уверен, что, если бы не эта печальная ошибка, основное направление современного гармонического анализа выглядело бы по-другому. Алгоритм объяснен здесь - попробуйте. Я надеюсь, что мы сможем восполнить этот 90-летний пробел.
И последнее: мое глубочайшее уважение (в алфавитном порядке) Джулиусу Бартелсу и Артуру Шустеру. Они сделали лучшее из того, что было возможно в их некомпьютерную эпоху:
"Если я и видел дальше, то только стоя на плечах гигантов" (сэр Исаак Ньютон).
Сергей Тарасов
17 сентября 2021 года
Оригинал: https://www.timingsolution.com/TS/Articles/Bartels!/