Диссертация (1150844), страница 19
Текст из файла (страница 19)
В частном случае shallow = 1 модель принимает вид«экспонента плюс константа». anterior и shallow — коэффициенты перед соответствующим экспонентам.В качестве модели исходных данных использована модель () + (), где() представляет из себя шум.5.3.2. Оценивание параметров сигналаВ [62] для оценивания параметров сигнала использован метод ESPRIT(точнее говоря, LS-ESPRIT [4, 5]). ESPRIT даёт оценки для anterior и shallow ,после чего anterior и shallow оцениваются по стандартному методу наименьших133квадратов, так как они линейно входят в модель (5.2).
Так как считается, что головная компонента убывает быстро (anterior < 1), а пологая компонента близкак константе (shallow находится возле единицы), то оценки параметров упорядочены соответствующим образом.В дополнение к оценкам ESPRIT для тех же самых профилей построенаоценка anterior и shallow с помощью Модифицированного метода Гаусса-Ньютона (MGN) с использованием алгоритма Cadzow в качестве метода поиска начального значения. Это возможно сделать, так как значения модели (5.2) наравномерной решётке представляют из себя ряд ранга 2.
Полученная в результате ОЛРФ() имела ненулевую последнюю координату, соответственно, anteriorи shallow были вычислены как корни соответствующего характеристическогоЛРФ полинома, см. определение 1.1.4. Шум () имеет сложную структуру,так как состоит из шума от препроцессинга и вариабельной активности генов.Поэтому шум не удовлетворяет модели стационарного процесса, хотя не слишком сильно отличается от белого шума. В качестве более простого вариантабыла взята матрица W = I , соответствующая случаю белого шума.Пример сравнения оценки сигнала для одного и того же наружного профиля методами ESPRIT и MGN представлен на рисунке 5.13, соответствующиеоценкам головной и пологой части для каждого из методов представлены нарисунке 5.14.
Заметим, что полученная с помощью метода MGN головная компонента более резко убывает, чем с помощью оценки по ESPRIT. Это связано стем, что полученная MGN оценка лучше аппроксимировала резкий спад на левом краю ряда (среднеквадратичная ошибка аппроксимации по ESPRIT равна13.62, для оценки по MGN равна 12.95).13420304050607080SeriesESPRITMGN20406080APРис. 5.13. Различные способы оценки сигнала в apical профиле.01020304050ESPRIT anteriorESPRIT shallowMGN anteriorMGN shallow20406080APРис. 5.14. Различные способы оценки частей сигнала в apical профиле.1355.3.3. Сравнение характеристик полученных оценокПосле очистки данных, описанной в [62], были построены доверительныеинтервалы по каждой из трёх групп эмбрионов: Cleavage, nc10-13 и nc14. В ка(apical)честве основных рассматриваемых характеристик в [62] используется anterior —(apical)(basal)оценка головной экспоненты наружного профиля и = log(anterior /anterior )— независимое от линейных преобразований отношение коэффициентов передэкспонентой наружного и внутреннего профиля.
Для каждой из характеристик были построены 90-процентные доверительные интервалы для оценки поESPRIT и для оценки по MGN. Результаты представлены в таблицах 5.1 и 5.2(apical)(для краткости anterior записано как ).Таблица 5.1. Средние значения и границы 90-процентных доверительных интервалов харак(apical)теристик anterior и , оценённых по методу ESPRIT.Среднее 5% 95% Среднее 5% 95%10-13c0.8790.8620.895-0.170-0.300-0.04114c0.8710.8620.8800.6520.5080.797Cleavage0.8150.7850.8460.2830.1240.443Таблица 5.2. Средние значения и границы 90-процентных доверительных интервалов харак(apical)теристик anterior и , оценённых по методу MGN.Среднее 5% 95% Среднее 5% 95%10-13c0.8450.8270.863-0.151-0.249-0.05314c0.8490.8360.8620.5890.4670.710Cleavage0.7640.7290.7990.1670.0530.282Те же самые результаты в графическом виде с построенными для каждой из групп 80-процентными доверительными эллипсоидами представлены нарисунке 5.15.136221CAB_esprit1CAB_mgn●fgroup●●10−13c14cCleavage●●●0●●●●●●●●●●●0●●●●●●fgroup●●●●●●●●●●●●●●● ●●●●●10−13c14cCleavage●●●●●−1−10.60.7А0.80.9apic_esprit_ant0.61.0Б0.70.80.91.0apic_mgn_ant(apical)Рис.
5.15. anterior и для оценок по методу (А) ESPRIT (Б) MGN.Более того, были сравнены ошибки аппроксимации сигнала, построенныхразличными методами. Результаты представлены в таблице 5.3. Видно, что оценённый с помощью метода MGN сигнал более точно аппроксимирует данные(среднеквадратичная ошибка аппроксимации по каждой паре группа/профильуменьшилось). Сравнение средних оценок параметров в таблицах 5.1 и 5.2 по(apical)казывает, что среднее значение характеристики anterior уменьшилось, и, следовательно, привело к более резкому поведению оценённой головной компоненты.Таблица 5.3. Средние значения ошибок аппроксимации для оценок сигнала, полученных методами ESPRIT и MGN.apical, ESPRIT apical, MGN basal, ESPRIT basal, MGN10-13c18.1617.678.978.6614c17.0216.794.103.95Cleavage22.6421.6424.2622.701375.4.
Аппроксимация временными рядами конечного рангапри неизвестных параметрах авторегрессионногошумаРассмотрим задачу оценки сигнала по наблюдаемому ряду X = S + длины , где S — сигнал конечного ранга, — шум, например, представляющий из∑︀себя процесс авторегрессии порядка : = =1 − + , где — гауссовскийбелый шум с нулевым средним и дисперсией 2 . На практике, коэффициентыпроцесса авторегрессии неизвестны, то есть неизвестен точный вид ковариационной матрицы Σ, что не позволяет сразу решать задачу (3). Поэтому, задавпорядок процесса авторегрессии и ранг сигнала , оценку сигнала можно построить с помощью максимизации функции правдоподобия:maxY∈ , =(1 ,..., )T , >0log X,DΣ(,)D (Y),(5.3)где X,DΣ(,)D — функция правдоподобия многомерного нормального распределения со средним X и ковариационной матрицей DΣ(, )D, где Σ(, ) ∈R × — автоковариационная матрица процесса AR() с параметрами =(1 , .
. . , )T и (в явном виде обратная автоковариационная матрица указанав [49]), D — заранее известная диагональная положительно определённая матрица, задающая огибающую шума (вид обратной матрицы приведён в (4.26)).Для решения этой задачи воспользуемся покоординатным спуском, решая поочерёдно задачу аппроксимацию рядом конечного ранга и поиска коэффициентов процесса авторегрессии.
Ниже данный подход записан в виде алгоритма.Алгоритм 5.4.1 (Покоординатный спуск в случае неизвестных коэффициентов авторегрессии). Вход: 0 — начальный вектор параметров авторегрессии, — число итераций.1:Положить = 0.1382:while < или не выполнен критерий остновки STOP do3:Положить = + 1.4:Вычислить Y = arg maxY∈ log X,DΣ(−1 ,−1 )D (Y), что сводится к решению задачи (3) для ряда X и W = D−1 Σ−1 (−1 , −1 )D−1 .5:Вычислить ( , ) = arg max=(1 ,..., )T , log X,DΣ(,)D (Y ), что сводится к задаче поиска коэффициентов процесса AR() для ряда D−1 (X − Y )(например, с помощью алгоритма [63]).6:end while7:return Y .Пример использования подходаВ качестве примера был взят ряд US Unemployment, male — ежемесячныеданные по безработице в тысячах человек с 1948 по 1991 год, длина ряда равна = 408, данные взяты из [64]. Тот же самый ряд рассматривался в качествепримера в статье [65], где для его исследования применялись техники улучшения разделимости метода SSA [8].К этому ряду был применён алгоритм 5.4.1 в предположении, что ряд из себя представляет сумму сигнала конечного ряда и гетероскедастичного красногошума.
Сделано предположение о том, что шум является «мультипликативным»,то есть огибающая шума (диагональ матрицы D) примерно пропорциональнатренду. Поэтому была сделана «грубая» оценка тренда с помощью примененияметода SSA [8, 21] с небольшой длиной окна = 84 и взятия первой собственнойтройки, таким образом, метод SSA использовался для сглаживания исходногоряда [8]. Ряд вместе с оценкой тренда представлен на рисунке 5.16.Так как задача (5.3) может иметь много локальных минимумов, алгоритм5.4.1 может сходиться к различным решениям при различных начальных коэффициентах авторегрессии. Было взято 0 = (0.9, 0, .
. . , 0)T и сделано 10 итера139300020001000USUnemployment, MALE4000SeriesRaw trend estimate0100200300400IndexРис. 5.16. Ряд US Unemployment, male и оценка тренда по методу SSA.ций алгоритма 5.4.1, то есть изначально предполагалось, что шум красный, чтобы получить как можно более простую модель при подборе . Выбор параметровбыл произведён по байесовскому информационному критерию (BIC) [66], который в данном случае выглядит как −2 log(X,DΣ( * ,* )D (Y* )) + (2 + ) log ,где * , Y* , * — решения задачи (5.3), 2 + — общее число степеней свободы. Сравнением BIC с соседними значениями параметров был получен рангряда = 16 и порядок процесса авторегрессии = 3.