АСТРОНОМИЯ
УДК 519.23, 521.17
Р. В. Балуев
ИССЛЕДОВАНИЕ СТАТИСТИЧЕСКИХ СВОЙСТВ ЭКЗОПЛАНЕТ МЕТОДАМИ НЕПРЕРЫВНОГО ВЕЙВЛЕТ-АНАЛИЗА1
Введение. К настоящему времени открыто уже более 200 планет у звезд солнечного типа [8]. Статистика выборки такого объема позволяет достигнуть точности порядка 1 /л/200 « 0.07 (это число грубо характеризует, скажем, относительную ошибку некоторого параметра, но его не следует связывать с вероятностями, возникающими в статистических критериях). Такой точности достаточно для некоторых заключений относительно распределений различных параметров планетных систем, хотя она и соответствует «пределу чувствительности», предъявляя высокие требования к применяемым статистическим методам и критериям.
В этой статье мы будем использовать для анализа различных распределений вне- солнечных планет метод вейвлет-преобразований. Этот математический инструмент является весьма эффективным в самых разных приложениях, связанных с обработкой данных [3]. Вейвлет-преобразование позволяет разложить исходные «сырые» данные на совокупность структур определенной формы (которая задается видом применяемого вейвлета), но имеющих различные положения в пространстве исследуемых параметров и разные масштабы. При этом чувствительность алгоритма обработки к структурам именно этого вида значительно возрастает по сравнению с классическими методами анализа. Первые результаты этой работы были мною представлены на международной конференции [2]. В этой статье я более подробно опишу использованную методику обработки данных и приведу обновленные и более полные результаты такой обработки выборки открытых экзопланет.
Основные теоретические моменты. В ряде работ уже предпринимались попытки адаптации вейвлет-анализа к задаче исследования распределений. Предлагалось использовать дискретные вейвлет-преобразования, проводя разложение плотности распределения по счетной системе ортогональных вейвлетов (см. обзор [1]). Ортогональность дает преимущества в распространенной ситуации «сигнал + шум», когда шум
1Работа выполнена при финансовой поддержке РФФИ (гранты №05-02-17408 и 06-02-16795) и Совета по грантам Президента РФ для государственной поддержки молодых российских ученых и ведущих научных школ (грант №НШ -4929.2006.2).
© Р. В. Балуев, 2008
является аддитивным (и даже гауссовым). Но в указанной задаче шум не аддитивен и имеет принципиально иную природу. В этой ситуации ортогональность является излишним требованием, неоправдано ограничивающим форму вейвлетов. Дискретное вейвлет-преобразование зависит от параметров (начальных масштаба и смещения), не имеющих физического смысла для непрерывного распределения. По этим и другим причинам здесь представляются более удобными непрерывные вейвлет-преобразования.
Непрерывный вейвлет-анализ одномерных распределений. Сформулируем основные положения общей теории вейвлет-преобразований [3], учитывая специфику задачи анализа функций распределений случайных величин. Мы будем использовать дифференциальные вейвлеты
■ iP 'I'
j^\.
(1)
где функция Ф(в) пригодна для локального сглаживания, т.е. хорошо локализована
+ Л
около нуля вместе со своим Фурье-образом Ф ( ш ) = / л^е^^Б и удовлетворяет нор-
_ Л
мировке Ф(0) = 1. Для некоторой функции í (х) эти вейвлеты продуцируют последовательность дифференциальных вейвлет-преобразований
Замечу, что ради дальнейшего удобства мы использовали в этом определении вейвлет- преобразования нестандартную нормировку. Если функция f (х) есть плотность вероятности некоторой случайной величины £, то мы будем называть функции (2) вейвлет- преобразованиями этой случайной величины и обозначать также, как (а,Ь). Для краткости индексы Ф,А ниже могут опускаться. Функция í (х) восстанавливается с помощью формулы обращения
Требование СЦ < +то выражает собой условие допустимости вейвлетов фА (б). Для дифференциальных вейвлет-преобразований выполнено равенство
Щ (а, Ъ) = Щ1,о[£ «](а, Ъ), (4)
весьма важное при анализе статистических данных. Действительно, для подавления шума необходимо что-либо сглаживать (усреднять). Сглаживать саму исследуемую функцию (как это делает гистограмма) нерационально, при этом подавляются и реальные детали ее структуры. Информация о таких деталях сохраняется в производных исследуемой функции. Например, при поиске структур типа обрыва выгодно исследовать первую производную í (х) и использовать вейвлет л^б), а при поиске экстремумов (где велика кривизна графика { (х)) —вторую производную с помощью вейвлета л(б). Выбор сглаживающей функции Ф(в) может быть в значительной степени произволен и диктуется скорее вопросами математической простоты и удобства в обращении. Мы будем использовать гауссиану ф(в) = е~Б /2/У27т, которая порождает семейство гауссовых вейвлетов ф (б) с константами в формуле (3) вида С] = (] — 1)!
Оценка вейвлет-преобразования и критерии согласия. Пусть дана случайная выборка X = {х] |л=1. Значения XI считаются не отягощенными ошибками наблюдений (измерений). Очевидно, что оценкой вейвлет-преобразования (скалограммой) как математического ожидания (2) является выборочное среднее
W> = 7vh*tb(A)- <«
i=1 v 7
Необходимо описать погрешность этой оценки. В информативном случае (при наличии точной гипотезы о WjAi (a, b)) нужно построить некоторый критерий согласия. Для этого нам нужно построить такую статистику Z, что гипотеза WjAl(a,b) = W(a,b) отвергается при Z > z, а также доверительную вероятность (или ее нижнюю границу) как функцию a(z). Вероятно, наиболее удобно будет принять Z = maxa,b |Y(a,b)|, где Y(a,b) = aj+1 (Ix(a,b) — W (a,b)). В этом случае значение Z сохраняется при линейных преобразованиях £. К сожалению, добиться точной инвариантности относительно существенно более широкого класса преобразований вряд ли возможно. Однако в приложении статистика Z оценивается сверху случайными величинами с известными распределениями, не зависящими от f (x). Подобные оценки дают доказательство существования критерия, основанного на статистике Z. К сожалению, эти оценки оказались слишком грубыми для использования при нашем N ~ 100. Найти существенно лучшую оценку элементарными методами уже не представляется возможным. Было проведено Монте-Карло-моделирование статистики Z для гауссовых вейвлетов фА(А), различных f (x). Выяснилось, что реальные распределения Z почти не зависят от f (x) и хорошо аппроксимируются распределениями Гумбеля вида exp A- exp -Aj AzyHV — BjAj a . Эти оценки Aj, Bj приведены в таблице.
Результаты моделирования распределений статистики ^. _________________________________________________________________________Число модельных измерений N = 100.
Распределение f(x) к-во испыт. A1 Bi K-BO испыт. A2 B2
Нормальное 8697 14.75 0.230 18130 11.16 0.318
Равномерное 7455 14.23 0.231 24171 11.25 0.316
Коши 8897 16.56 0.224 16000 12.32 0.312
Смесь нормальныха 7814 14.40 0.231 17470 11.30 0.318
Итоговая оценка 14.2 0.231 11.1 0.318
а Отношение стандартов гауссиан а : 2а, весов 3 : 7, взаимное смещение на 2а.
Оценка плотности распределения. Для применения критериев согласия необходима некоторая гипотеза, полностью определяющая распределение. Часто такая априорная модель отсутствует. В этом, неинформативном, случае требуется найти оценку f(x), имеющую возможно меньше статистически необоснованных структур (наболее простую, еще согласующуюся с выборкой). Наш критерий согласия определяет точки (а,Ъ), в которых эмпиричиская скалограмма значима (значимо отлична от нуля). Дополнив ее в остальных точках нулями, получим функцию переменных (а, Ъ), имеющую наименьшее количество значимых деталей, еще согласующуюся с выборкой. Применив к этой очищенной скалограмме формулу обращения (3), получим функцию, вейвлет- преобразование которой ближе (в смысле применяемого критерия согласия) к скало- грамме. Такие итерации повторяются до получения функции, вейвлет-преобразование которой значимо не отличается от скалограммы. Эта функция и представляет собой требуемую оценку. Можно заметить, что в ходе работы этого алгоритма данным не навязывается какой-либо функциональной модели. Здесь, как и в дискретном вейвлет- анализе, может применяться «жесткая» или «мягкая» чистки [1]. По большому счету, эти виды чистки скалограммы представляют собой лишь разные способы ее интерпретации и не дают качественно различных результатов. Данный алгоритм родственен известному методу CLEAN, применяемому при анализе временных рядов и в радиоастрономии при обработке изображений.
Одномерные распределения экзопланет. Метод вейвлет-преобразований был использован для анализа распределений некоторых важных параметров, характеризующих экзопланетные системы. В вычислениях участвовали 187 планет, открытых на март 2007 г. методом лучевых скоростей [8]. Использовался вейвлет A2(s), позволяющий выявлять концентрации и разрежения плотности вероятности.
На гистограмме логарифмов орбитальных периодов P экзопланет (рис. 2) помимо выраженного семейства долгопериодических планет наблюдается концентрация около значения 4 сут. Аналогичная концентрация присутствует и на гистограмме логарифмов больших полуосей A орбит экзопланет. Этот максимум «горячих юпитеров» был отмечен довольно давно [5]. Однако стандартные оценки ошибок в бинах гистограммы не позволяют признать его достоверным. Кроме того, вид гистограммы зависит от размера бинов, что влечет неустойчивость результата. Интегральное распределение также дает низкий уровень достоверности структуры «горячих юпитеров». Кроме того, оно значительно труднее поддается визуальному анализу.
Скалограммы распределений величин ln P и ln A похожи (рис. 1). Этого следовало ожидать: массы большинства центральных звезд близки к солнечной, что влечет линейное соотношение ln P « 1.5ln A + const. На указанных графиках на фоне глобальной структуры распределения выделяется пик, образованный семейством короткопери- одических планет — в основном «горячих юпитеров» (относительно недавно этот класс пополнился некоторым количеством менее массивных планет—«горячих нептунов»), занимающий интервал 2 — 7 сут по периодам и 0.03 — 0.07 а.е. по большим полуосям. Относительный масштаб этого образования порядка единицы. Кроме того, имеется разрежение («долина периодов», «period valley») сравнимого масштаба, занимающее промежуток 20 — 100 сут и 0.2 — 0.4 а.е. для периодов и больших полуосей соответственно. Значимость этих структур составляет около 0.96. Но для каких-либо уверенных заключений относительно структуры этого семейства имеющейся выборки экзопланет уже недостаточно. На рис. 1 приводятся результаты работы описанного выше алгоритма чистки. Степень концентрации «горячих юпитеров» сильно зависит от уровня чистки. В класс короткопериодических планет входит около 20% открытых экзопланет.
На скалограмме логарифмов минимальных масс внесолнечных планет (рис. 3) не обнаруживается никаких деталей, кроме крупномасштабной структуры. Наблюдается максимум около значения m « 1.5 Mjup и плавный спад в сторону как меньших,
0 50 ----- 0.80 —- 0.90 ----- 0.95 . 0.98 --- 0.99 -
Рис. 1. Графики слева относятся к распределению логарифмов периодов (сут), справа — больших полуосей (а.е.). Верхняя пара изображает их скалограммы (контуры равных доверительных вероятностей, отрицательные вероятности условно соответствуют отрицательным значениям скалограммы). Смещение Ь (имеющее размерность исследуемой величины) отложено по горизонтальной оси, масштаб а (безразмерен)—по вертикальной. Ниже —
результаты «жесткой» и «мягкой» чистки распределений для ряда уровней значимости, на фоне соответствующих гистограмм.
так и больших масс. Падение плотности с увеличением массы отражает реальную особенность распределения, в то время как завал на меньших массах по крайней мере частично объясняется наблюдательной селекцией. Важно отметить, что здесь участвовали только наблюдаемые массы m = M sin i (i —угол наклона плоскости орбиты планеты к картинной плоскости). Чтобы получить распределение истинных масс M А m, необходимо провести коррекцию на эффект неизвестного наклона.
На рис. 3 хорошо виден резкий узкий максимум распределения эксцентриситетов вблизи e = 0. Эти почти круговые орбиты относятся, как правило, к «горячим юпитерам» (хотя частично эта деталь объясняется тем, что при недостаточности временного ряда наблюдений эксцентриситет полагается нулевым априорно и не оценивается). Далее плотность распределения остается приблизительно постоянной, спадая только после значения е « 0.4. Эксцентриситет — неотрицательная величина, но при анализе эта до-
1 10 100 1000 10000 1 10 100 1000 10000
Рис. 2. Гистограмма (слева) и эмпирическая функция распределения (справа) логарифмов орбитальных периодов экзопланет (единица периода — сутки). Гистограмма нормирована на полное число планет. Приведены доверительные интервалы уровня значимости «1а» (доверительная вероятность 0.68) для соответствующих бинам вероятностей и доверительные границы интегрального распределения того же уровня значимости. Заметим, что
доверительные границы для всей совокупности бинов гистограммы существенно шире, чем показанные границы для одиночных бинов.
полнительная информация не учитывалась, в связи с чем график на рис. 3 продолжен в область e < 0. Также виден некоторый дефект алгоритма восстановления плотности распределения: рядом с резким максимумом e = 0 наблюдаются слабый боковой лепесток, признаков которого нет на скалограмме.
Выводы и перспективы. Можно утверждать, что вейвлет-анализ оказался значительно эффективнее классических методов оценки функций распределений, по крайней мере при обработке выборки из N ~ 200 экзопланет. Например, удалось установить, что статистическая значимость выделения «горячих юпитеров» в отдельное семейство довольно высока (не менее 95%). Однако при этом мы не учитывали факторы, искажающие исследуемые распределения. Приведем два наиболее важных из них.
1. Эффект неизвестного наклона. Подавляющее большинство экзопланет открыты методом лучевых скоростей. Этот метод не позволяет определить массу планеты M, наблюдаемой величиной является m = M sin i (i —неизвестный угол наклона орбиты планеты к картинной плоскости).
2. Наблюдательная селекция. Планеты, дающие малую амплитуду кривой лучевых скоростей K, открываются труднее, чем планеты с большей K. Для обнаружения планеты необходимо проследить кривую лучевых скоростей на достаточно длинном интервале времени, порядка ее орбитального периода.
Математически наблюдательная селекция сводится к умножению плотности распределения масс, полуосей, эксцентриситетов на функцию селекции, определяющую вероятность открытия планеты с данными характеристиками. Эффект неизвестного наклона, в свою очередь, описывается интегральным уравнением, родственным уравнению Абеля ([7] и др.).
Вейвлет-анализ, основанный на применении к выборке некоторых интегральных преобразований, позволяет учесть такие искажающие эффекты. Интегральное преобразование искаженной функции распределения можно представить как аналогичное интегральное преобразование исходной функции распределения, но с другим ядром.
Рис. 3. Распределение логарифмов минимальных масс M sin i экзопланет (слева) и эксцентриситетов их орбит e (справа). Пояснения аналогичны пояснениям к рис. 1.
Можно сказать, что селекция переносится с интересующей нас функции распределения на ядро этого преобразования, которое, в отличие от распределения, известно точно. Решение этих задач предполагается осуществить в будущем.
Приложение: Оценки статистики Z. В рамках непараметрического критерия Колмогорова гипотеза о том, что истинное интегральное распределение описывается функцией К(х), отвергается при к = шаххЕЙ |К(х) — ^[х](х) > z (здесь F[X](x) — эмпирическая функция распределения). Соответствующая доверительная вероятность a(z) асимтотически (М А оо) выражается функцией K(zy/N) ~ 1 — 2в~2Мг . Заметим, что
Интеграл в (7) —константа, выражающая полную вариацию Уф вейвлета ф. Если базовая гипотеза верна (F(x) есть истинное распределение £), то с вероятностью, не меньшей K{zAfN/У.ф), для всех а, b выполнено IY(a, 6)| Л z. Таким образом, базовая гипотеза отвергается на доверительном уровне а, если в какой-либо точке оказалось IY(a, b) | > z, где z определяется из K{zЛfN/Vл) = а. При этом а не зависит от F(x).
Записав производную в (6) как разность положительной [*]+ = max{0, *} и отри-
+ж +ж
цательной [*]_ = — min{0, *} частей, с учетом равенств / A'(s)] + ds =J л ' (s)]_ds =
\Г(к,с).ф\^+^-уф. (8)
Здесь к+ = тахжек{ FX](x) — F (x)} и к_ = таххек{ F(x) — F[X](x)} —односторонние статистики Смирнова. Совместное предельное распределение статистик Смирнова известно [4]. Поэтому несложно получить функцию распределения их полусуммы и сравнить его с распределением
Колмогорова:
1 - 2^(1GA;V -/, 1
i ь; ^
i: = 1
Уф / 2, имеем более точную оценку
Очевидно, к±Ак, поэтому неравенство (8) является более сильным, чем (7), в связи с чем предпочтительнее использовать функцию P(z) вместо K(z).
Summary
R. V. Baluyev. Exploration of the Exoplanetary Statistics by Means of the Continuous Wavelet Analysis.
An effective algorithm of the analysis of statistical distributions, based on usage of continuous wavelet transforms is constructed. Some features of the wavelet analysis of statistical data are described, the practically important statistical tests are built. This technic was applied to exoplanetary statistics. A number of univariate distributions of important parameters (exoplanetary minimum masses, orbital periods, semimajor axes and eccentricities) are examined. The possible prospects for this technic development (like accounting for the observational selection and for the effect of unknown inclination) are discussed.
Литература
1. Abramovich F., Bailey T. C., Sapatinas T. Wavelet analysis and its statistical applications // The Statistitian 2000. Vol.49, N 1. P. 1-29.
2. Baluyev R. V. Statistics of masses and orbital parameters of extrasolar planets using continuous wavelet transforms // Arnold L., Bouchy F., Moutou C. (eds.), Tenth Anniversary of 51 Peg b: Status of and Prospects for Hot Jupiter Studies. Frontier Group. Paris, 2006. P. 103-110.
3. Daubechies I. Ten lectures on wavelets // S.I.A.M., Philadelphia, 1992.
4. Королюк В. С. и др. Справочник по теории вероятностей и математической статистике. Киев: Наукова думка, 1978.
5. Ксанфомалити Л. В. Закономерности внесолнечных планетных систем и роль метал- личности звезд в образовании планет // Астрон. вестн. 2004. Т.38, №5. С.428-439.
6. Холшевников К.В., Кузнецов Э.Д. Эффект селекции в больших полуосях орбит вне- солнечных планет // Астрон. вестн. 2002. Т. 36, №6. С. 504-515.
7. Zucker S., Mazeh T. Derivation of the mass distribution of extrasolar planets with MAXLIMA, a maximum likelihood algorithm // ApJ 2001. Vol. 562. P. 1038-1044.
8. Schneider J. The extrasolar planets encyclopaedia // http://www.exoplanet.eu, 2007.
Статья поступила в редакцию 11 сентября 2007 г.