Диссертация (1150844), страница 17
Текст из файла (страница 17)
Более того, былиизмерены число обусловленности (A ) матрицы A , отношение первого к -му̂︀ 2 , обозначенные (L ) и (L̂︀ 2 ) соответственно.сингулярному числу для L и LРезультаты эксперимента представлены на рисунке 5.1. Для наглядности числа обусловленности домножены на 10−18 , чтобы привести к масштабу другихграфиков.В легенде на рисунке 5.1 вдобавок к измеренным величинам указана оценка порядка полиномиального роста, построенная с помощью линейной регрессии. Например, порядок 3.00 означает кубический рост.
Видим, что порядокΔ (3.01), Δ2 (2.88) и (A ) (3.04) примерно равен = 3, что согласуется стеоремой 3.2.1. Из рисунка 5.1 также видно, что, предположительно, порядок̂︀ 2 ) равен ( −1 ). Таким образом, порядок ошибкиобусловленности (L ) и (Lвычисления проекции для алгоритмов 3.2.1 и 3.2.2 определяется обусловленностью A .Аналогично, ошибки вычисления проекции были сделаны и для алгоритмов 3.2.3 и 3.2.4, результаты представлены на рисунке 5.2. Видно, что при рассмотренных компенсированная схема Горнера даёт настолько большое улучшение по точности, что порядок полиномиального роста (A ) при указанных11810−18κ(Ag) (3.04)1e−0910−18κ(Lr) (2.01)^10−18κ(L2r) (2.01)1e−151e−12value1e−06∆A (3.01)∆A2 (2.88)501002005001000200050001000050000NРис. 5.1. Сравнение ошибок вычисления проекции для различных и чисел обусловленности. был оценён как линейный.5.1.2.
Алгоритмы локального поискаДля того чтобы сравнить локальные методы решения задачи (3) по точности оценки решения, необходимо построить такую пару из исходного временногоряда X и ряда Y* = S0 , для которой было бы известно, что S0 является точкой локального минимума задачи (3). Потом необходимо запустить алгоритмлокального поиска из окрестности точки локального минимума и сравнить полученный результат с известным минимумом.⋆В качестве Y= (1 21 , . . .
, 1 2 )T взята квадратичная функция, управляемая ОЛРФ(0 ), 0 = (1, −3, 3, −1)T , где — равномерная решётка на отрезке⋆̂︀ =[−1; 1], положительное 1 взято так, чтобы ‖Y‖ = 1. Далее, был взят ряд N̂︀ ‖ = 1, вычислен N = N̂︀ − Π(2 ),W N̂︀ , и(3 |1 |, . . . , 3 | |)T такой, что ‖N0⋆*= S0 и X = X удовлетворяетвзят X = Y+ N . Таким образом, пара Y1e−09∆A (3.01)∆A, compensated (0.99)∆A2 (2.88)∆A2, compensated (0.85)1e−151e−12value1e−06119501002005001000200050001000050000NРис. 5.2.
Сравнение ошибки вычисления проекции при различных с компенсированнойсхемой и без неё.условиям леммы 2.4.2. Однако, лемма 2.4.2 даёт только необходимые условия локального минимума. Достаточное условие локального минимума (положительная определённость матрицы Гёссе целевой функции ‖X − S(( ) , ( ) )‖2W [35,Theorem 2.3]) было проверено численно при < 1000.Для численного построения примера необходимо строить проекцию на(20 ).
Для этого нужен алгоритм, который очень точно строит базис (20 ).Подпространство (0 ) состоит из полиномиальных последовательностей степени не более 2, (20 ) из степени не более 5, что объяснено в подразделе 5.1.1.Для вычисления Z(20 ) с наибольшей точностью были вычислены полиномыЛежандра [60] от 0 до 5 порядка в точках . Затем полученный базис был подвергнут дополнительной ортогонализации. Таким образом, пример построен.Теперь перейдём к сравнению методов.Для простоты ограничимся случаем, когда W — единичная матрица. Мыпровели сравнение методов VPGN, S-VPGN и MGN, описанных в разделе 3.3.5,120при различных от 20 до 50000 с использованием компенсированной схемыГорнера для вычисления базисов подпространств () и (2 ).
Детали применения алгоритмов для вычисления базисов описаны в разделе 5.1.1. В качестве меры ошибки было взято евклидово расстояние между результатом ра̃︀ ⋆ и точкой локального минимума Y⋆ . Дополнительно, былаботы алгоритма Ỹ︀ ⋆ управиспользована мера для расчёта того, насколько полученное решение Yляется последней вычисленной ОЛРФ(⋆ ).
Для этого была вычислена невязка̃︀ ⋆ ‖‖QT (⋆ )Y.⋆‖ ‖В качестве алгоритма локального поиска по направлению Δ на шаге 8алгоритмов 3.3.4 и 3.3.5 возьмём простой перебор = 1, 1/2, 1/4, . . . , 2−16 , 0 дотех пор, пока не выполнится улучшение решения относительно предыдущей()()итерации: ‖X − S* (( ) + Δ )‖W ≤ ‖X − S* (( ) )‖W . В этом случае, в качествекритерия остановки можно использовать условие = 0, т.е. невозможностьполучить более близкое решение, чем на предыдущей итерации. Проблема состоит в том, что если мы вначале возьмём = 1 (соответствующее полномушагу метода Гаусса-Ньютона), а относительное изменение будет мало, скажем⃦⃦⃦ * ()() ⃦*⃦ S (( ) + Δ ) − S (( ) ) ⃦⃦⃦ < ,(5.1)⃦⃦()*S (( ) )⃦⃦где имеет порядок корня из машинного эпсилон (мы используем = 5 · 10−8 ),то указанный локальный поиск нельзя применять вследствие недостаточной()точности вычисления целевой функции ‖X − S* (( ) + Δ )‖2W , сводящейся квычислению плохо обусловленных скалярных произведений.Давайте модифицируем сам метод локального поиска по направлению дляслучая, когда выполнено условие (5.1).
Как метод MGN, так и метод VPGN()можно рассматривать как последовательность параметров ( ) , так и после()довательность рядов S* (( ) ). Если критерий по последовательности рядов неработает, можно использовать критерий по набору параметров. Положим = 1,121если ≥ 0 и ‖Δ−1 ‖ > ‖Δ ‖, и = 0 в противном случае, т.е. завершим алго(+1)ритм тогда, когда очередной шаг по параметрам ‖Δ ‖ = ‖( )()− ( ) ‖ сталбольше предыдущего из-за достижения локального минимума с максимальновозможной точностью. Таким образом, мы получим комбинацию из метода локального поиска по направлению и критерия остановки, которая работает даже тогда, когда численное вычисление целевой функции может не показыватьулучшения из-за недостаточной точности его вычисления.В качестве начальных коэффициентов ОЛРФ был использован вектор 0 +⋆представлены на10−6 (1, 1, 1, 1)T .
Результаты для метрики-расстояния до Yрисунке 5.3, для невязки — на рисунке 5.4.⋆. Заметим, что алгоВсе три метода сошлись в одну и ту же окрестность Yритм S-VPGN не даёт существенного улучшения по расстоянию до решения посравнению с методом VPGN, в то время как решения метода MGN оказываютсясущественно ближе решений методов S-VPGN и VPGN из-за более устойчивоговычисления шага. К тому же заметим, что используемая в локальном поискеи изображённая на рисунке 5.3 граница из условия (5.1) повлияла на видошибки для алгоритма MGN — благодаря перебору длин шагов , который,условно говоря, работает над этой границей, удалось сделать решение близким,несмотря на падение точности при вычислении шага, но эта граница не повлияла на результат методов VPGN и S-VPGN из-за худшей точности вычисленияшага. С точки зрения невязки решения алгоритмы S-VPGN и MGN идентичныи имеют значительное преимущество по устойчивости над методом VPGN.1e−081e−04MGNVPGNS−VPGNCriterion threshold1e−161e−12Distance to solution1e+0012220501002005001000200050002000050000N1e−06Рис.
5.3. Сравнение алгоритмов по расстоянию до ответа при различных .1e−121e−151e−18Disperancy1e−09MGNVPGNS−VPGN20501002005001000200050002000050000NРис. 5.4. Сравнение алгоритмов по невязке при различных .1235.2. Исследование свойств оценок сигнала с помощьюстатистического моделирования5.2.1. Исследование оценок, полученных методом КэдзоуДля исследования работы предложенных в главе 4 алгоритмов были проведено численное моделирование. Сравнение происходило на двух рядах длины = 40 и ранга = 2:(1)(1)(1)1. Синусоидальный сигнал S(1) = (1 , .
. . , )T , где = 5 sin 26 ,(2)(2)(2)2. Линейный сигнал S(2) = (1 , . . . , )T , где = 100 .()Были промоделированы ряды X = S() + (1 ), где — независимые реализации длины процесса авторегрессии порядка 1 с параметром 1 : = 1 −1 + ,где — гауссовский белый шум с нулевым средним и дисперсией 2 , — номерреализации, = 1, . . . , . Мы рассматривали 1 = 0 — случай белого шума и1 = 0.9 — случай красного шума. Затем к полученным промоделированнымрядам был применен алгоритм Cadzow(quadratic) из раздела 4.6.1, и для полу()ченных оценок ̂︀S была оценена их точность с помощью MSE (mean squarederror) — среднеквадратической ошибки от соответствующего сигнала S() , гдесреднее взято и по длине ряда, и по числу реализаций:1 1 ∑︁ ̂︀()MSE =‖S − S() ‖2 .=1В качестве критерия остановки STOP в алгоритме Cadzow(quadratic) было использовано условие‖X −X−1 ‖W‖X0 ‖W< 10−6 .
Заметим, что в этом случае результатыалгоритмов Cadzow(quadratic) и Cadzow(box) из раздела 4.6.1 совпали (так какалгоритмы поиска весов дали одинаковый ответ). В качестве длины окна былавзята половина длины ряда = 20, следуя рекомендации в [61].124Для начала, мы промоделировали = 10000 реализаций каждого ряда,и к одними и тем же реализациям мы применили метод Cadzow(quadratic) приразличных , использующихся при решении задач поиска весов (4.7) и (4.8).Напомним, что чем ближе к 0, тем ближе матрица весов R в (4) к вырожденной, но тем лучше может быть задаваемая ими аппроксимация весов W в (3).Для MSE были построены 95-процентные выборочные доверительные интервалы («Simulation»), а так же вычислена теоретическая оценка MSE линеаризованного алгоритма Кэдзоу («Theory») c помощью замечания 4.7.4 и границаРао-Крамера с помощью предложения 2.5.2 («CRB»).
Более того, было вычислено среднее число шагов, затраченное алгоритмом Cadzow(quadratic) («Numberof iterations»), и на том же графике изображена величина (4.32) («Convertedeigenvalue»), объясняющая число шагов.Для начала, мы рассмотрим сигнал S(1) и шум c 1 = 0, то есть случайбелого шума, и с двумя разными уровнями шума = 1.5 и = 0.5. Результатыпредставлены на рисунках 5.5 и 5.6.SimulationTheoryCRBMSE of Cadzow signal estimate ~ alpha, b1 = 0 sigma= 0.5SimulationTheoryCRB●●0.0340.32MSE of Cadzow signal estimate ~ alpha, b1 = 0 sigma= 1.50.30●●●●●0.032●●●0.28●●●0.26●●0.030MSEMSE●●●●●●0.028●●●●●0.24●●●●●0.026●●●●●●●●●●●0.22●0.01●●0.02●●●●●●●●●А●●●0.050.10alpha0.200.501.000.01●●0.02Б●0.050.100.200.501.00alphaРис.
5.5. Сравнение для оценки сигнала S(1) с помощью алгоритма Cadzow(quadratic)для белого шума при различных и (A) = 1.5 (Б) = 0.5.Прокомментируем результаты на рисунках 5.5 и 5.6. Заметим, что под1251012100880606●●●−1/log(max non−unit eigenvalue)4●●4●Number of iterationsConverted eigenvalue●40●Mean number of iterations660●40Mean number of iterations880●Average speed of Cadzow algorithm, b1 = 0 sigma = 0.510Number of iterationsConverted eigenvalue●−1/log(max non−unit eigenvalue)100Average speed of Cadzow algorithm, b1 = 0 sigma = 1.5●●●●●0.01А0.020.050.10alpha0.20●2●2020●2●●●●●0.50●●1.000.01Б0.020.050.10●0.20●●●0.50●●1.00alphaРис.
5.6. Сравнение числа итераций алгоритма Cadzow(quadratic) для сигнала S(1) для белогошума при различных и (A) = 1.5 (Б) = 0.5.твердился теоретический результат леммы 4.7.2 для оценок ошибок алгоритмаКэдзоу, однако заметим, что в случае (Б) на рисунке 5.5, в котором соотношениесигнал/шум больше, чем в случае (А) ( = 0.5 против = 1.5, что должно датьувеличение MSE примерно в 9 раз), практический эксперимент оказался ближек теории (теоретические значение ближе к середине доверительных интервалов), что естественно, так как сходимость ошибок первого порядка возникаетпри соотношении сигнал/шум, стремящемся к бесконечности.