Замечания к использованию главных компонент в математическом моделировании

МАТЕМАТИКА

DOI: 10.18721/JPM.11307 УДК 519.24

ЗАМЕЧАНИЯ К ИСПОЛьЗОВАНИЮ ГЛАВНЬ1Х КОМПОНЕНТ В МАТЕМАТИЧЕСКОМ МОДЕЛИРОВАНИИ

Ю.А. Пичугин

Санкт-Петербургский государственный университет аэрокосмического приборостроения,

Санкт-Петербург, Российская Федерация

В статье рассматриваются вопросы, связанные с использованием метода главных компонент (ГК) в математическом моделировании: априорная оценка дисперсии погрешности регрессии на ГК для случаев большой и малой выборки; оценки минимального риска; оценка структурного подобия; проблема пропущенных данных и прогноз нестационарных временных рядов. Во всех случаях автор предлагает свои решения рассматриваемых вопросов или исправляет ранее допущенные неточности.

Ключевые слова: главная компонента, априорная оценка, дисперсия погрешности регрессии, малый объем выборки

Ссылка при цитировании: Пичугин Ю.А. Замечания к использованию главных компонент в математическом моделировании // Научно-технические ведомости СПбГПУ. Физико-математические науки. 2018. Т. 11. № 3. С. 74-89. DOI: 10.18721/JPM.11307

NOTES ON USING THE PRINCIPAL COMPONENTS IN THE MATHEMATICAL SIMULATION

Yu.A. Pichugin

Saint-Petersburg State University of Aerospace Instrumentation, St. Petersburg, Russian Federation

The paper discusses the issues related to the use of principal components analysis (PCA) in mathematical simulation. The paper significantly expands the range of the solved problems using PCA. In particular, the solutions of the following three tasks are given: (i) structural similarity and homogeneity estimation for random Gaussian vectors; (ii) recovery of missing data; (iii) the forecast of non-stationary time series based on the caterpillar method, which is a generalization of PCA for non-stationary time series. To solve the problems, to restore missing data and to predict the data, the author offers an unbiased estimation of the variance of the error of the regression on the PCs base for the cases of large and small samples. All the main statements are formulated in the form of theorems proved by the author.

Key words: principal component, variance of the regression error, small sample volume

Citation: Yu.A. Pichugin, Notes on using the principal components in the mathematical simulation, St. Petersburg Polytechnical State University Journal. Physics and Mathematics. 11 (3) (2018) 74-89. DOI: 10.18721/JPM.11307

Введение

Анализ, или метод главных компонент (далее сокращения для главных компонент — PCs, а для метода — PCA (Principal Component Analysis)) представляет собой хорошо известный аппарат математической статистики. Этот метод был предложен К. Пирсоном в 1901 году [1], и суть PCA состоит в следующем.

Если должным образом повернуть (посредством ортогонального преобразования) систему координат я-мерного пространства так, чтобы оси координат совпадали с главными осями эллипсоида рассеяния, то компоненты нормально распределенного я-мерного центрированного вектора будут некоррелированными и, в силу нормального закона распределения, независимыми.

В алгебраическом смысле это есть не что иное, как приведение ковариационной матрицы к диагональному виду путем ортогонального преобразования, а квадратичной формы в экспоненте функции плотности многомерного нормального распределения — к каноническому виду. Хорошо известное преобразование Карунена — Лоева [2, 3] есть, по сути, именно это преобразование координат. Переход же к независимым переменным PCs позволяет, как правило, существенно сократить размерность исследуемой задачи с минимальной потерей информации.

В связи с этим в литературе PCs нередко выводятся как решение оптимизационной задачи, хотя все их оптимальные свойства достаточно хорошо видны из самого спектра ковариационной матрицы (спектр показывает, какова доля отбрасываемой дисперсии, см. далее).

В анализе временных рядов PCA известен как анализ сингулярного спектра (SSA — Singular Spectrum Analysis), где посредством данного метода решается проблема избыточности классического спектрального анализа [4 — 6]. Особенность SSA состоит в том, что размерность вектора в этом случае равна N, а размерность матрицы взаимных ковариаций равна N * N (N - длина исследуемого временного ряда). При этом элементы ковариационной ма-

трицы вычисляются особым способом, когда делитель, независимо от величины сдвига, а соответственно и от числа слагаемых, равен N. Такие оценки, очевидно, относятся к классу смещенных оценок, но именно они не приводят в спектральном анализе к искажению (в сторону завышения, как отмечают Г. Дженкинс и Д. Ваттс [7]) длины волны. Проблема большой размерности ковариационной матрицы в реализации SSA, как показано в работе [8], легко решается применением итераций фон Мизеса, так как все строки ковариационной матрицы временного ряда могут быть получены из первой строки посредством сдвига, дублирования и перестановки элементов. Собственные значения и собственные векторы ковариационной матрицы получаются последовательностью простых итераций без вращения матрицы размерности N * N.

Альтернативой SSA является «метод гусеницы» [9], а также метод, предлагаемый нами далее в настоящей работе (см. разделы «Проблема относительно малой выборки» и «Прогноз нестационарных временных рядов»), где схема прогноза строится на основе метода гусеницы.

Следует отметить, что существует множество методов, близких к PCA, например, в методе независимых компонент (ICA) последние могут подчиняться не только распределению Гаусса, но распределениям Стьюдента, Коши, Дирихле. Отметим, что метод независимых компонент известен также как анализ указанных компонент (ICA — Independent Component Analysis).

Обобщением PCA является метод главных кривых и многообразий. В последнее время PCA широко используется для визуализации и графического представления многомерных данных (рассматривается проекция выборки на плоскость первых двух главных осей [10, 11]). При этом требование нормального распределения исходных данных не возникает.

Здесь мы имеем большое разнообразие достаточно близких по сути методов, таких как многомерное шкалирование, нелинейный маппинг, поиск наилучшей проекции, а также методы нейросетевых задач, такие как метод «узкого горла», самоорганизую-

щиеся карты Кохонена и т. п. Следует также отметить, что графическое представление многомерных данных проекцией на плоскость первых двух главных осей PCs позволяет получить достаточно хорошее начальное приближение разделения выборки в решении задачи классификации в работе [12].

Цель настоящей работы — расширение спектра задач, решаемых на основе метода главных компонет.

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

Краткое описание математического аппарата PCA

Предполагается, что вектор y имеет размерность m (dimy = m) и подчиняется многомерному нормальному распределению, т. е. y ~ N(0у, ^).

Пусть P - ортогональная матрица, такая, что

PTV^P = Л = diag(X:, Х2,..., %m) и X, >Х2 > ... > X„,

где Т - символ (оператор) транспонирования.

Напомним, что столбцы матрицы P есть собственные векторы матрицы Vy, а совокупность собственных чисел {X:, X2, ..., Xm} называется спектром этой матрицы. В математической статистике столбцы матрицы P называются базисом главных компонент, а главными компонентами называются компоненты вектора PT (у - 0у). В приложениях параметры распределения N(0y, V), как правило, неизвестны. При наличии выборки {yj, j = 1,2,..., п} мы можем вычислить несмещенные оценки неизвестных параметров:

0 y = 1Z y jj

j=1

Vy =

n -

1Z (y j -

1 j=1

0 )(У; - '0y )T.

В этом случае в качестве P берем ортогональную матрицу Г, которая приводит к диагональному виду оценку V , т. е.

PTVP = Л = diag

и X1 > X2 > ... > Xm

Очевидно, что матрицы Vy и Vy в принципе не равны, следовательно, не равны и матрицы Р и Р.

Задача снижения размерности, как и другие задачи, рассмотренные в данной работе, непосредственно связана с проверкой следующей гипотезы:

Н: X, > Х2 > ... > X, > Х,+1 = Хк+2 = ... = Хм.

Принятие этой гипотезы позволяет рассматривать вектор меньшей размерности

Р(! )(у - 0у),

где матрица Р(к) содержит только первые к столбцов матрицы Р.

Первым известным тестом для проверки гипотезы Н (см. формулу (3)) был тест Бартлетта [13, 14] (см., например, монографию [15] или справочник [16]). Однако, если первое цитирование теста Бартлетта [15] предполагать безупречным, то во втором обнаруживаются сразу две неточности (знак и множитель).

Действительно, в работе [15] х2-стати-стика для проверки гипотезы Н выражается как

Yn = n' \ (m - к)ln

1

где

n' = n - к — 6

m - к г=к+1

" !n I П Ъ

ZX ,j-

2( m - к) +1 +

тогда как в справочном издании [16] эта статистика следует выражению

Yn = (n - 1) \(m - k)ln

1

m - k :

,1+

=k+1

Оба издания указывают одинаковое число степеней свободы n-В работе [15]:

П = 1 (m - k + 2)(m - k - 1); В работе [16]:

П = 1 (m - k + 1)(m - k) - 1В настоящее время широкую известность приобрел тест, называемый правилом «сломанной трости» (англ. Broken stick model, см, например, статью [17])- Согласно этому тесту, следует выбрать максимальное значение k, при котором выполняется неравенство

X k 1/ k + 1/(k + 1) + ...+1/ m

trV

m

где ^ — след матрицы.

Наиболее простым средством определения k является визуальный анализ графического представления спектра ковариационной матрицы в порядке убывания. В этом случае в качестве k берется значение, которое предшествует смене относительно быстрого убывания собственных значений на относительно медленное (плавное), что, в сущности, повторяет метод сломанной трости (на интуитивном, неформализованном уровне). Во многих справочных и учебных изданиях предлагается определять k из соотношения

/-ч - . /-ч 1

/ . 1100% - K%,

где K — заранее установленный процент общей дисперсии.

Выбор метода определения k, в конечном итоге, зависит от характера решаемой задачи. Приступим к изложению результатов.

Очевидно, что при снижении размерности возникает ошибка, дисперсия которой

должна быть как-то связана с отбрасываемой частью выборочного спектра. Решению этой задачи посвящены два следующих раздела статьи-

Априорная оценка дисперсии погрешности регрессии на PCs

Регрессией вектора y на PCs (на компоненты вектора z) называется связь, выраженная уравнением вида

y = + V +

Согласно общим принципам классического регрессионного анализа, предполагается, что

г ~ N

где 0 — нулевой вектор, а I — единичная матрица соответствующей размерности-

Предположение (5) означает, что

Е

и cov(s, sj) = 0 (i ф j),

где E, var и cov — операторы математического ожидания, дисперсии и ковариации, соответственно-

Пусть в описанных выше условиях рассматривается некоторая реализация вектора y. Это может быть одна из тех реализаций, которые использовались при вычислении оценок параметров (Qy, Vy), т. е. yj, (/ = 1,2,...,«), или одна из последующих yn+l, (l > 1), что не принципиально. Поэтому нижний индекс пока опускаем. В практических задачах моделирования вместо уравнения (4) мы имеем регрессионную модель вида

у = Qy + Р

Прежде всего, заметим, что, согласно классическому регрессионному анализу, уравнение (6) следовало бы называть регрессией на базис PCs, так как сами PCs (компоненты вектора z) на этапе применения построенной модели (6) являются определяемыми параметрами. Регрессия именно на PCs фигурирует в доказательстве следующей теоремы.

теорема 1. В условиях модели (6), где элементы модели вычислены по формулам

(1), (2), при верной (принятой) гипотезе Н

(см. выражение (3)) и предположении (5),

2

априорная несмещенная оценка а выражается как

а2 =

n -1

— £ х,.

(т - к)(п - к - 1) г=к+1

Доказательство. В регрессионной модели (4), которая связывает три нормально распределенных центрированных вектора (у - 0^, г и г), перейдем к ковариациям. Тогда имеем следующее равенство:

Vy = Р№, Z + Vs ,

где

V, = а21, V = ^(а?,а2,..., ак),

а2 = уаг(гг), (/ = 1,2,...,к).

Умножим равенство (8) слева на Р(, а справа - на Р. В результате имеем:

diag(Xl, Х2 , ... , Хк , Xk+1,..., Хт ) =

= diag

+ diag(а2, а2,..., а2),

т. е.

Х. = а2 + а2 (/ = 1,2,...,к),

Х(. = а2 (/ = к + 1, к + 2,..., т).

Здесь следует особо отметить, что не только равенство (9) справедливо лишь при верной гипотезе Н (см. (3)), но и предположение (5) в случае большой выборки также возможно лишь при верной гипотезе Н. Это непосредственно видно из равенства (9), где Х(. = а2. Однако эта логика нарушается в случае малой выборки (см. далее). Из равенства (9) также следует, что

а2 =

1

( - £) i 1

Рассмотрим классическую модель регрессии [18]:

XJ = Ро + e1 Zj ,1 + Р2 Zj ,2 + ... ... + Р,-1 Zj,,-1 +Sj (j =1,2,..., n).

Для произвольно выбранного значения 1

подставим значения выборочных PCs, вычисленные по формуле

(Zj,1,Zj,2,..., Zj,2 = (Уj - 0y)TР(к), где y j — j-я реализация y в исходной выборке.

В последнем случае мы оставляем прежние обозначения. Это PCs исходной выборки, но индексы (номера) компоненты и реализации (l и j) меняем местами, приводя в соответствие стандартам регрессионного анализа, т. е. модели (11). В данном контексте эти значения PCs

{Zлl, Zj, 2,..., Zj; к} считаем известными и модель (11) здесь действительно есть регрессия на PCs. Тогда уравнение (11) будет соответствовать i-й строке матричного равенства (6).

Из условия

Ъ = 1 ¿Zj,, = 0 (l=1,2,..., к

j=1

следует, что

во = х,

1 n „ 1 n

где х = ~Zxj, т.е. ро = - Z ytj (см. выше).

j=1

j=1

PCA обеспечивает минимальность остаточной суммы квадратов, как и метод наименьших квадратов (OLS — Ordinary Least Squares), следовательно, вектор-строка OLS-оценок

(в 1, в 2,..., в ,-1) совпадает с i-й строкой матрицы рк), т.е.

(Р1, Р2,..., Р,-1) = [Р(*)], где [ Р(к)] — i-я строка матрицы рк).

Это утверждение нетрудно проверить непосредственным вычислением. Определим следующие матрицы исходных данных:

Y = (y1 -0,y2 -0,...,yn -0),

z = Рк ).

Тогда аналог модели (6) для этих матриц можно записать в виде

Y = ]р(к)+Е,

где E — (m х n)-матрица всех остатков регрессии.

Соответственно, аналог (11) запишется в виде

¥г = Т Г Р(£,+Е Г. Далее нам нужно проверить равенство

Р(£) = ( г1 ТYГ,

которое и означает, что матрица Р(к), а, соответственно, и все ее строки — суть ОЬ8-оценки.

Подставляя выражение для матрицы Т (см. выше), умножая справа на Р(к)Р

и

расставляя дополнительные скобки, имеем очевидное тождество:

( Р ) =

к )г

= (РГ У^ГР Г1(PГYY:ГP )Р

(к)'

так как Р(£) Р(Л) = I. _ _

Корректность умножения на Рда Р(£) основана на том, что ранг матриц при этом не снижается:

гапк(Р(к)Р(к)) = гаикР(к) = к.

Из принятых предположений и равенства (10) следует, что остаточная сумма квадратов S2 регрессионной модели (11) выражается как

п -1

52 =

( - к) г=!+1 В соответствии с теорией линейной регрессии [18], оценка

52 п-1

ст2 =

Е ^ I

п - р (т - к)(п - к - 1) ,•1

является несмещенной (напомним, что р = к + 1).

Теорема 1 доказана.

Отметим, что оценка (7) и краткий набросок доказательства теоремы 1, содержащий основную идею (см. формулу (9)), были предложены ранее автором настоящей статьи в работе [19]. В ряде задач, касающихся оценки информативности (см. [20]), требуется именно смещенная оценка ст2. В этом случае необходимо брать оценку

ст 2 = п -1 Е £ I,

(т - к)п I^ 1 которая непосредственно следует из формулы (10). К вопросу об использовании

априорной оценки мы вернемся далее (см. раздел «Восстановление пропущенных данных»).

Проблема относительно малой выборки

Во многих задачах нередко возникает ситуация, когда объем выборки п меньше размерности вектора т (п <т). В этом случае

1 = 4+2 = ... = 4 = 0,

что никак не соответствует равенству (9), а, следовательно, и (10). Это обстоятельство, однако, не исключает возможности проверки гипотезы

Н1: А-1 > ^ > ... > Хк > Хк+1 =

= ^к+2 = ... = ^ п

для дальнейшего рассмотрения модели (6).

Если считать вектор 6у и матрицу Р(к) известными, то для какой-либо из последующих реализаций вектора у, например (п + /)-й (I > 1), несмещенная апостериорная оценка ст2 в регрессионной модели (6) следует выражению

1 т

т-ЧЕ

т — к I=1

ст2 =

где у — 1-я компонента априорной оценки вектора 0у (/ = 1,2, ...,/и), а число оцениваемых параметров равно к.

Если же считать известной лишь матрицу Р(к), тоэтаоценка будет следующей:

1 т

—,-Т Е &, п+ I Уп + 1

т - к -11=1

ст2 =

[Р(к)]^п +1 ) ,

1

где уп+1 =—Е^,п+/; у, п+1 — |-я компонента

т £1 , вектора уя+,.

В обоих случаях, отвечающих формулам (13) и (14), в этих формулах, как и выше, [Р(к)]I — 1-я строка матрицы Р(к) ,а

£ я+1 = % )(Уя+/ -6 у ).

Однако в выражении (13) вектор 0у вычислен по формуле (1) с использованием всех реализаций исходной выборки, а в (14) мы заменяем все компоненты 0у

средним значением по компонентам новой реалгоации уя+, (0 = 0 = ...=0т = уп+,, см выше). Однако вектор 0у используется при вычислении Vy и Р, что необходимо учесть в априорной оценке а2. В связи с этим рассмотрим иной метод вычисления элементов модели (6).

Определим иначе матрицу У. Пусть

V = (У1 - 01, у 2 - 02,..., у п - 0п),

где компоненты каждого вектора 0. равны между собой и равны среднему по компонентам реализации у. т. е.

1 т

01. =02. = ... = 0 . = у. = —У V. .

1] 2] т] '] .Г,,.

п -1

т 1=1

(] =1,2,..., п).

Вычислим оценку V, по формуле

^^у =

Далее произведем все перечисленные выше действия: вычислим Р и {Х1, Х2, ..., Хп}; проверим гипотезу Н1 (см. (12), п < т) и определим матрицу Р(к). При этом могут несколько измениться не только V и Р, но и значение к.

Теорема 2. В условиях модели (6), где элементы модели вычислены по схеме (15), при верной (принятой) гипотезе Н1 (см. (12)) и предположении (5), априорная оценка а2, равная

а2 =

п -1

(т - к - 1) п 1=к+1

У Х,,

является несмещенной.

Доказательство. В силу линейности оператора математического ожидания, среднее значение несмещенной оценки (14), вычисленное по исходной выборке (1=1,2,...,и) и есть несмещенная априорная оценка, равная

1

1 п

=У

п т~1 т - к - 1

(у, - 0]]Т

(I - Р(

Р0(у, - 0]] =

1

(т - к -1) п 1

(т - к -1) п

1г(У((I - Р(к)Р(к))У) = (1гУ(У - \гЪТ Ъ) =

УХ,.

(т - к -1) п !1

(здесь, как и прежде, Z = Р(к) У).

Теорема 2 доказана.

Заметим, что оценка (16) не требует совпадения математических ожиданий оценок средних (0.), а требует лишь совпадения математических ожиданий оценок (14). Это означает, что выражение (16) подходит для случая временных рядов, содержащих тренды. В свою очередь, это может быть иным основанием (не только относительная малость объема выборки) к применению формул (15) и (16). Сравнивая (7) и (16), замечаем также, что в случае малой выборки (п < т) некорректное использование (7) будет давать завышенное значение оценки а2. Смещенным вариантом оценки (16) будет оценка

п -1

а2 =

тп

УХ,.

Если принять во внимание тот факт, что индекс ] в математическом моделировании нередко соответствует некоторому временному отсчету, то при достаточной стационарности и несущественном различии средних значений компонент вектора у, вычисленных первым способом (см. оценки (1)), можно применять оценку (16), не прибегая к описанной выше коррекции модели (6). Возникающие при этом ошибки и неточности будут несущественны. И, наоборот, в нестационарном случае следует использовать выше изложенный метод оценивания ^. При этом подразумевается пересчет 0у (см. выше) по реализации у, непосредственно фигурирующей в модели (6).

Рассмотренный в этом разделе метод оценивания элементов модели (6) наиболее эффективен, когда метод гусеницы применяется к прогнозу нестационарных временных рядов с явно выраженным трендом (см. раздел «Прогноз нестационарных временных рядов»).

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

невырожденного распределения и большого объема выборки мы можем смотреть на главные компоненты как на регрессоры, а на элементы матрицы Р(к), представляющей собой базис главных компонент, как на оцениваемые параметры, что мы видим в доказательстве теоремы 1. Случай же малой выборки (преобладание размерности над объемом), наоборот, приводит нас к необходимости рассматривать элементы матрицы Р(к) как регрессоры, а главные компоненты — как оцениваемые параметры (теорема 2). При этом следует осуществлять центрирование вычитанием не среднего значения по реализациям для каждой компоненты, а вычитанием среднего по компонентам для каждой реализации (см. выше). Суть проблемы состоит в том, что при малой выборке равенство (9) не имеет выборочного аналога, так как выборочный спектр не полон и логика теоремы 1 рушится. Поэтому приходится менять точку зрения (обзора ситуации).

И, наконец, отметим, что случай принципиально вырожденного распределения, когда выборочный спектр неполон и большой объем выборки, приводит нас к вычислительной схеме данного раздела и оценке (16).

Оценки минимального риска

При выборе модели (6) для использования в каких-либо прикладных задачах обычно применяется формула:

У = 6 у + Р

где г = Р((у- 6у).

В справочном издании [21] предлагается применять формулу

у = 6 у + Р

где С = g2,..., gk).

При этом значения §,. определяются из условия минимума квадратичного риска

Я2 = )2 (I =1,2,...,к).

В этой ситуации компоненты вектора а называются оценками минимального риска. С учетом того, что

- §Л = (! - §1) - §1 (£ - )'

имеем равенство

Я2 = ^ - 1)2 + §2а|,

где Е(|- - г,-) = ст|.

Приравнивая к нулю производную Я 2 по §., получаем равенство

1

или ё

1 + ст |/л2

Подставляя вместо оценку минимального риска , получаем простое квадратное уравнение вида

£2 - ё, +82 = 0, §2 =а|/£2 со следующей окончательной формулой для ё; :

§ = 2 11

В справочном издании [21] авторами предлагается следующий алгоритм вычисления величины :

_1 = 2

1 § 2 1

— - §,, §2 < —;

4 ' 4'

1

= 0, если 82 > —.

Автором настящей статьи в работе [22] предложен альтернативный алгоритм:

1 12 ^ 1 § = 24~8г', 8г " 4;

1 1 §2 < 1

§. = -, если ^ < §2 <1;

ё =

1

1 + 8;

-, если 82 > 1.

Теорема 3. В предположении, что различие между величинами ст|/^2 и ст2./^2 пренебрежимо мало, предлагаемый алгоритм (21) обеспечивает меньшее значение квадратичного риска Я2, по сравнению с алгоритмом (20), для случая 82 > -1.

Доказательство. При условии §2 >>

алгоритм (20) дает величину Я2 = z2. Если 1 2

— < §; < 1, подставляем в (18) значение 4 1

= —. Учитывая, что ст< I2, получаем,

что

* 2= + ст 2 )< z2< Z

Если же 82 > 1, то, согласно алгоритму (21), значение g. таково, что выполняется равенство (19). Подставляя (19) в (18), получаем, что

R2 = gf(z2+aР) - 2giZf + Z2 = = giZf - 2gtzf + zf = zf(1 - gt) < Z2.

Теорема 3 доказана.

Полученный в данном разделе результат можно отнести к уточнению деталей. При исследовании различных многомерных процессов и явлений, в особенности природных, не всегда использование PCs подчинено задаче снижения размерности, а может быть направлено на исследование внутренней структуры явления, что показано в следующем разделе данной работы.

Структурное подобие и однородность

Предположим, что кроме вектора y исследуется и вектор x такой же размерности:

dim x = dim y = m.

Предположим, что имеется выборка реализаций этого вектора. Объемы выборок реализаций векторов x и y могут быть различны. Используя формулы (1) для вектора x, вычислим оценки параметров распределения 0 x и Vx. Пусть ортогональная матрица Q такова, что выполняется равенство

QTVQ = diag

Р, >Д2 > ... >£m.

В этой ситуации мы имеем два множества статистических характеристик PCs обоих векторов:

{Xi, Pi }m=1 и {Д,, q,. }m=1,

причем q,., p,. есть i-е столбцы матриц Q и р, соответственно, т. е.

Q = (q1, q2,.., qra),

P = (p1, p2,..., Pm ).

Оценкой коэффициента структурного

подобия векторов x и y назовем величину

¿\/ДЛ

qT Pi

M trVy

Этот коэффициент показывает, насколько структуры колебаний исследуемых векторов согласуются в относительных долях дисперсии. В некоторых случаях 5' целесообразно вычислять по формуле

ZJ

Рф

X;

qT(i) Pi

s, = max -

xy

ф

M trVy

где ф(i) есть перестановка индексов, т. е. варьируется порядок статистических характеристик PCs одного из векторов (здесь это x ).

Необходимость применения формулы (25) может возникнуть в случае обнаружения близких собственных значений хотя бы в спектре одного из векторов. Если верна (тестирована) гипотеза

Нx: Д1 > Д2 > ... > Pi > Р1+1 = Pi+2 = ... = Pm,

то можно рассматривать отфильтрованный коэффициент структурного подобия в виде

р f =_i=1_

г,—p— ,

JZ P. Z Xi

V i=1 i=1

где p = min(k, I),

или относительный коэффициент структурного подобия вида

ZT^ |qT p.-1

. _i=1_

M trVy

Оценки коэффициентов и 5гу в случаях, когда соседние собственные значения мало различимы, целесообразно находить по формулам, аналогичным формуле (25).

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

i=1

ластей и районов. В микроэлектронном производстве, когда в каждой ячейке кристаллической пластины производится несколько типов микроэлектронных приборов (см., например, работу [23]), коэффициент структурного подобия удобен для оценки, насколько погрешности изготовления различных приборов связаны между собой и зависит ли эта погрешность от положения ячейки на кристаллической пластине. Такие вопросы могут возникать и в процессе настройки оборудования. В анализе временных рядов базой для применения коэффициента структурного подобия является упомянутый во введении анализ сингулярного спектра (SSA). Первые попытки построения коэффициента структурного подобия были предприняты автором данной работы именно при анализе временных рядов [24].

Если в качестве выборки вектора x используется выборка значений того же вектора у, а второй набор характеристик PCs {Д,, q m вычисляется по оценке корреляционной матрицы R , т. е.

QTRyQ = diagfai, Д2, ..., Дm),

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

Возникает естественный вопрос о возможности проверки гипотезы о равенстве коэффициента структурного подобия нулю, т. е. Hs: sxy = 0. Прежде всего, отметим, что на практике мы всегда имеем лишь оценку s , а истинное значение sxy мы имели

xy' xy

бы лишь в случае использования в наших вычислениях матриц Vy, Vx, P и Q, что возможно лишь гипотетически.

Рассмотрим два составных вектора:

Ys = (V^р г,<ДГP p m )T =

= (Y1, Y2, ..., YM ) ;

Xs = (Л/Д7qf q 2,...^Л/ДГq m)T =

= (X1, X2, ... , XM ) ,

где M = m x m.

Нумерацию подвекторов

л/ДД-

составного вектора X 3, при необходимости, можно установить в соответствии с формулой (25). Это может нарушить лишь условие (23), что не имеет для нас значения. Знаки некоторых столбцов матрицы

0 =

нужно поменять на противоположные так, чтобы все произведения 5р( или р.. были положительными. Эта смена знаков не нарушит равенства (22). Тогда, если все эти условия выполнены, из ортогональности матриц Р и 0 следует, что

Sxy

XT Ys

I X,Y,

VXT X s YT Ys

Е х21 ъ2

I 1 =1 1 =1

Рассмотрим два регрессионных уравнения, связывающие компоненты составных векторов У5 и X 5:

= р Х1 + е,.;

у.= Ь0 + Ь X. + в.(1 =1,2,..., М),

где предполагается, что остатки (погрешности) кождой из регрессий подчиняются нормальному закону распределения, взаимно независимы и имеют одинаковую дисперсию.

Теорема 4. Если принята (не отвергнута) гипотеза Не: Е(е) = 0 (или эквивалентная гипотеза Н0: Ь0 = 0), то, при верной гипотезе Н8: я = 0, величина

i

1 - s2

xy

t,

т. е. подчиняется распределению Стьюдента с числом степеней свободы М - 1.

Способы проверки гипотез Н0: Ь0 = 0 и Не :Е(е) = 0 хорошо известны. Поэтому рассмотрим простое доказательство, которое мало чем отличается от хорошо известного для обычного коэффициента корреляции.

Доказательство. Из общей теории метода наименьших квадратов, предположения Б(е) = 0 (см. гипотезу Н6) и формулы (26) следует, что ОЬ8-оценка параметра регрессии (27) имеет вид

м

= 5х

У х,г,

в=^м—

Ух2

.=1

уаг(р) =

У 12

У х,2

а

У X2

где а2 = уаг(е 1), а

У х,. 8,. = 0.

Из равенства (30) следует, что несме-

щенная оценка а выражается как

а2 =

1

м -

1 ! м

- У (1, -рх,. )2 =

1 ,=1

м-

71У^2-в2 У X?

1 V '=1 '=1

Очевидно, что равенство 5 = 0 эквивалентно в = 0 (см. (29)). Следовательно, принимая во внимание выражение (31), имеем:

Р<

м

Ух2

1

м ,_

У?; 24м-1

а

м м

-в2Ух,2

4

1 - 52

Теорема 4 доказана.

Как было отмечено выше, смена всех знаков матрицы ^ не нарушает формулы (22), но меняет знак 5 и знак (28) на противоположный. Это вполне соответствует симметрии распределения Стьюдента. Нетрудно заметить, что гипотеза Н5 может быть проверена и для отфильтрованного коэффициента структурного подобия . Различие сводится лишь к тому, что в этом случае М = т х р, а из = 0 следует и

= 0. Невозможность проверить гипотезу Н5 (гипотеза Н8 отвергнута) не является принципиальным препятствием для вычисления 5. В большинстве приложений мы скорее заинтересованы в том, чтобы отвергнуть гипотезу Н5, т. е. установить структурное подобие. Исключение составляет случай, когда мы собираемся установить неоднородность некоторого поля или временного ряда вышеуказанным способом.

Перейдем к рассмотрению задач практического характера, где непосредственно используются результаты разделов «Априорная оценка дисперсии погрешности регрессии на РСз», «Проблема относительно малой выборки» и «Оценки минимального риска».

Восстановление пропущенных данных

Рассмотрим ситуацию, когда распределение невырожденное и объем выборки достаточно превышает размерность. Для уменьшения количества индексов и прочих знаков перепишем уравнение (6) в виде

у = у + Бг + г,

где у = 0у, Б = Р(к).

Пусть уп+1 (I > 1) — какая-либо очередная реализация вектора у, которая имеет и измерений и V пропусков (т = и + V).

Тогда, при соответствующей нумерации, имеем следующие разбиения на блоки:

Уп+1 = (У( , у2 ) , У = (У( , у2 ) ,

Б = (БТ, Б2Т )Т,

где блоки у1, у1 и Б1 соответствуют измеренным компонентам уп+/, а блоки у2, у2 и Б2 — пропущенным данным.

Вычислим оценку у2 по формуле

у2= у2+Б2С

где матрица в (см. формулу (17)) имеет размерность к х к, а ее компоненты вычисляются по алгоритму (21).

При этом в алгоритме

а| = а2[( Б( Б1)-1], „

где [...],] — оператор взятия элемента матрицы с указанными номерами строки (;) и столбца (]); ^ — ,-я компонента вектора г

,=1

=1

,=1

[=1

г

м-1

ОЬ8-оценок:

2 = (Ц1Ц)-1^(у: - у:) (/ = 1,2,..., к).

Из-за отсутствующих строк столбцы матрицы Е1 не ортогональны, матрица ЕТЕ1 не диагональна, а компоненты вектора 2 зависимы между собой:

V = а2(ЦТ Ц)-1.

Апостериорная оценка ст2 в этом случае имеет вид (см. выражение (13)):

1

ст2 =

"(У 1 - У:) X

и - к

х (I - ВДТ)-1 Щ - у,).

Однако именно в этой ситуации естественно воспользоваться более надежной априорной оценкой (7) как вычисленной по выборке существенно большего объема. Эта оценка используется во всех следующих формулах данного раздела.

Здесь следует отметить, что у2 есть, в сущности, центр условного распределения у2. Действительно, безусловное или априорное распределение у2 в некотором приближении (при достаточно большом п) есть N (у 2,У2), где У2 получается из Уу удалением строк и столбцов, которые не имеют отношения к у2. Условное или апостериорное распределение у2, также лишь в некотором приближении, есть N (у 2, ст21), что непосредственно следует из соотношения (5) и уравнений (32) и (33).

Для пропущенных в количестве V компонент вектора у можно указать относительно точное распределение Стьюдента:

У, - У1

ст. 1 + 1 + ГгС( Т

'„,

где число степеней свободы п = п - к -1, если используется оценка (7) и п = т - к -1, если по какой-либо причине используется оценка (16); % — 1-я строка Е (/ = и + 1, и + 2,..., т; т. е. строка % относится к блоку Е2).

По формуле (34) и выбранному уровню доверия 1 - а (а — уровень значимости) для этих значений получаем границы (1 - а)% -х стьюдентовских доверительных интервалов:

у, ± С/2ст^ 1 + п + Е^Е1 ,

где 'па/2 — квантиль распределения Стьюден-та с учетом двустороннего характера критерия (берется а /2 в силу симметрии распределения Стьюдента), (/ = и + 1, и + 2,..., т).

Целесообразность применения модели (32) и формулы (33) тем выше, чем существеннее неравенство

Ш1 +1 + %С(Е1 Ц)-^ <

1 п __

т. е. чем меньше ширина доверительного интервала, построенного по формуле (35), по сравнению с соответствующей шириной, где доверительный интервал построен по выборочным оценкам параметров распределения (Уу ), т. е. меньше неопределенность относительно восстанавливаемых данных.

Прогноз нестационарных временных рядов

Модель (32) и формула (33) могут использоваться для прогноза нестационарных временных рядов [25, 26]. Здесь модель (32) строится методом гусеницы, когда у есть скользящий отрезок (гусеница) временного ряда {у }^=1, а Шшу = т есть длина гусеницы.

При вычислении оценки ковариационной матрицы выборка строится пошаговым сдвигом, т. е.

У1 = (У1> У2>...» Ут )Т , У2 = (У2 , У3,..., Ут+1 )Т ,

Уп = (УN -т+1' УN -т+2' •" ' УN )Т ,

где п = N - т + 1, N — общая длина временного ряда.

Оценку ковариационной матрицы следует вычислять по алгоритму вычисления оценки (15). При определении длины гусеницы окончательное значение (здесь предварительные расчеты будут весьма полезны) т следует брать не меньше, чем период

волны, несущей наибольшую долю дисперсии, т. е. соответствующей £1. Стремление удовлетворить указанному требованию может приводить к ситуации малого объема выборки (п < т). В разделе «Проблема относительно малой выборки» было отмечено, что применение оценок (15) и (16) может быть обусловлено двумя причинами: малый объем выборки или нестационарность.

Авторы метода гусеницы [9] придерживаются классической схемы вычислений оценки ковариационной матрицы (см. разделы «Краткое описание математического аппарата РСА» и «Проблема относительно малой выборки») и ошибочно определяют число степеней свободы по минимальному размеру выборочной матрицы исходных данных (выбор п = п - к -1 или п = т - к -1 (см. выше) определяется через шт(т, п)). В методе гусеницы строки и столбцы матрицы исходных данных имеют одну и ту же природу, в связи с чем двойственность, которая изложена в конце раздела «Проблема относительно малой выборки», становится очевиднее. Однако в работе [9] метод гусеницы использовался преимущественно для фильтрации временных рядов, и данный вопрос о числе степеней свободы не стоял столь принципиально, как в случае прогноза по рассматриваемой здесь схеме.

При прогнозировании (в формуле (33)) у1 есть вектор последних и значений временного ряда:

У = (У N-и+1' Ум - и+2' "• ' Ум) , (и < т,, так как т = и + V, см. «Восстановление пропущенных данных»), а у2 — вектор прогнозируемых значений:

у2 = Ом+1' Ум+2'...' Ум+V ) .

Из вышеизложенного (см. раздел «Проблема относительно малой выборки») следует, что в качестве априорной оценки ст2 необходимо использовать оценку (16). Подлинный смысл прогнозирования заключается не столько в вычислении значений у2, сколько в построении достаточно узкой доверительной полосы (см. выше) для компонент у2. При этом в формулах (35) и (36) п = т - к -1, а вместо 1/п будет

стоять значение 1/т, если равные между собой компоненты (см. раздел «Проблема относительно малой выборки») вектора у вычислять как среднее по последним т значениям временного ряда, или 1/и, если вычислять среднее по последним и значениям (по компонентам), что вполне допустимо. В любом случае окончательный выбор параметров схемы прогноза определяется неравенством (36).

Заключение

В результате проведенного анализа существующих методов решения задач на основе главных компонет и предлагаемой модификации методов, можно сформулировать следующие основные итоги.

1. Получены оценки дисперсии погрешности регрессии на главные компоненты для случаев большой и малой (относительно размерности задачи) выборки и доказана несмещенность этих оценок. Полученные оценки являются важной деталью в схемах восстановления пропущенных данных и прогноза нестационарных рядов, предложенных автором настоящей работы. При этом несмещенность оценок служит непременным условием в построении доверительных интервалов для восстановленных или прогнозируемых значений (см. далее).

2. Теоретически обоснованно уточнены ранее известные оценки минимального риска, которые также использованы в перечисленных выше практических задачах.

3. Введен в рассмотрение коэффициент структурного подобия и теоретически обоснована статистика для проверки гипотезы о равенстве этого коэффициента нулю.

4. Предложены схемы восстановления пропущенных данных и прогноза нестационарных рядов. Указаны критерии применимости и доверительные интервалы для восстанавливаемых или прогнозируемых значений.

В заключение отметим, что в построении статистических моделей нам приходится выбирать, какие элементы модели наделить свойством статистической устойчивости и включить в модель, а какие нет. От адекватности этого выбора во многом зависит успех применения построенных моделей на прак-

тике, возможно даже в большей степени, чем от точности применяемых формул. В настоящей работе предполагается, что базис главных компонент как раз и является наиболее статистически устойчивой частью модели.

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

СПИСОК ЛИТЕРАТУРЫ

1. Pearson K. On lines and planes of closest fit to systems of points in space // The London, Edinburgh, and Dublin Philosophical Magazine and Journal of Science. 1901. Ser. 6. Vol. 2. No. 11. Pp. 559-572.

2. Karhunen K. Uber lineare Methoden in der Wahrscheinlichkeitsrechnung //Ann. Acad. Sci. Fennicae. Ser. A. I. Math.-Phys. 1947. Vol. 1947. No. 37. Pp. 1-79.

3. Лоев М. Теория вероятностей. М.: Изд-во иностр. лит-ры, 1962. 719 с.

4. Broomhead D.S., King G.P. Extracting qualitative dynamics from experimental data // Physica. D: Nonlinear Phenomena. 1986. Vol. 20. No. 2-3. Pp. 217-236.

5. Broomhead D.S., King G.P. On the qualitative analysis of experimental dynamical systems // Nonlinear Phenomena and Chaos (Malvern physics series). 1st edition. Sarkar S. (Ed.). Bristol: CRC Press, 1986. Pp. 113-144.

6. Ghil M., Vautard R. Interdecadal oscillations and the warming trend in global temperature time series // Nature. 1991. Vol. 350. No. 6316. Pp. 324-327.

7. Дженкинс Г., Ваттс Д. Спектральный анализ и его приложения. Вып.2. М.: Мир, 1972. 287 c.

8. Пичугин Ю.А. Итерационный анализ сингулярного спектра в оценке естественных цикличностей метеорологических наблюдений // Метеорология и гидрология. 2001. № 10. С. 34-39.

9. Главные компоненты временных рядов: метод «Гусеница». Под ред. Д.Л. Данилова и А.А. Жиглявского. СПб.: Пресском, 1997. 308 с.

10. The transform and data compression handbook. Eds. K. Rao, P. Yip. Boca Raton (USA): CRC Press, 2001.

11. Muresan D.D., Parks T.W. Adaptive principal components and image denoising // Proceedings of IEEE International Conference on Image Processing (ICIP). 14-17 Sept. 2003. Vol. 1. Pp. I-101— I-104.

12. Пичугин Ю.А. О классификации летних режимов погоды в Санкт-Петербурге // Метеорология и гидрология. 2000. № 5. С. 31-39.

13. Bartlett M.S. The effect of standardization on a х2 approximation in factor analysis // Biometrika.

1951. Vol. 38. No. 3-4. Pp. 337-344.

14. Bartlett M.S. A note on the multiplying factor for various 2 approximations // J. Roy. Statist. Soc. 1954. Vol. B16. Pp. 296-298.

15. Лоули Д., Максвелл А. Факторный анализ как статистический метод. М.: Мир, 1967. 144 с.

16. Айвазян С.А., Бухштабер В.М., Енюков И.С., Мешалкин Л.Д. Прикладная статистика. Классификация и снижение размерности. М.: Финансы и статистика, 1989. 607 с.

17. Jackson D. Stopping rules in principal components analysis: A comparison of heuristical and statistical approaches // Ecology. 1993. Vol. 74. No. 8. Pp. 2204-2214.

18. Себер Дж. Линейный регрессионный анализ. М.: Мир, 1980. 456 c.

19. Пичугин Ю.А. К проблеме статистического контроля данных наблюдений за приземной температурой на отдаленных станциях // Метеорология и гидрология. 2000. № 10. С. 18-24.

20. Пичугин Ю.А. Экологический мониторинг и методы многомерной математической статистики // Астраханский вестник экологического образования. 2012. № 2. С. 101-105.

21. Айвазян С.А., Бухштабер В.М., Енюков И.С., Мешалкин Л.Д. Прикладная статистика. Исследование зависимостей. М.: Финансы и статистика, 1985. 487 с.

22. Пичугин Ю.А. Учет сезонных эффектов в задачах прогноза и контроля данных о приземной температуре воздуха // Метеорология и гидрология, 1996. № 4. C. 52-64.

23. Михальчук А.С., Пичугин Ю.А. Дисперсионный анализ погрешностей технологических процессов микроэлектроники // Моделирование и ситуационное управление качеством сложных систем. Сб. докл. СПб.: ГУАП, 2017. С. 35-38.

24. Пичугин Ю.А. Естественные составляющие годового хода приземной температуры воздуха // Метеорология и гидрология. 1994. № 12. С. 34-43.

25. Пичугин Ю.А., Малафеев О.А. Оптимизация и прогноз в динамической модели управления портфелем ценных бумаг // Матер. секц. заседаний симп. «Нобелевские лауреаты по экономике и российские экономические школы».

СПб., 2003 г. СПб.: СПбГУ, 2003. С. 183 -185.

26. Пичугин Ю.А. Главные компоненты многомерных временных рядов: анализ и прогноз

// Сб. докл. XIII Междунар. научн. конф. по мягким вычислениям. СПб., 2010 г. Т. 1. СПб.: СПбГЭТУ «ЛЭТИ», 2010. С. 160-163.

Статья поступила в редакцию 26.03.2018, принята к публикации 21.06.2018.

СВЕДЕНИЯ ОБ АВТОРЕ

ПИЧУГИН Юрий Александрович — доктор физико-математических наук, профессор Института инноватики и базовой магистерской подготовки Санкт-Петербургского государственного университета аэрокосмического приборостроения, Санкт-Петербург, Российская Федерация.

190000, Российская Федерация, г. Санкт-Петербург, Большая Морская ул., 61. yury-pichugin@mail.ru

REFERENCES

[1] K. Pearson, On lines and planes of closest fit to systems of points in space, The London, Edinburgh, and Dublin Philosophical Magazine and Journal of Science, Ser. 6. 2 (11) (1901) 559-572.

[2] K. Karhunen, Über lineare Methoden in der Wahrscheinlichkeitsrechnung, Ann. Acad. Sci. Fennicae, Ser. A, I. Math.-Phys. 1947 (37) (1947) 1-79.

[3] M. Loeve, Probability theory, Vol. II, 4th ed., Graduate texts in mathematics, Springer-Verlag. 46 (1978).

[4] D.S. Broomhead, G.P. King, Extracting qualitative dynamics from experimental data, Physica. D: Nonlinear Phenomena. 20 (2-3) (1986) 217-236.

[5] D.S. Broomhead, G.P. King, On the

qualitative analysis of experimental dynamical systems, Nonlinear Phenomena and Chaos, S. Sarkar. (Ed.), CRC Press, Bristol (1986) 113-144.

[6] M. Ghil, R. Vautard, Interdecadal oscillations and the warming trend in global temperature time series, Nature. 350 (6316) (1991) 324-327.

[7] G.M. Jenrins, D.G. Watts, Spectral analysis and its aplications, Holden-Day, San Fracisco — Cambridge - London - Amsterdam, 1969.

[8] Yu.A. Pichugin, Iterative singular-spectrum analysis in estimating natural cyclicities in meteorological observation data, Meteorology and Hydrology. 10 (2001) 34-39.

[9] Glavnyye komponenty vremennykh ryadov: metod "Gusenitsa" [Principal components of time series: Caterpillar method], D.L. Danilov, A.A. Zhiglyavskiy (Eds.), SPbSU, St. Petersburg, 1997.

[10] The transform and data compression handbook, K. Rao, P. Yip (Eds.), CRC Press LLC, Boca Raton, USA, 2001.

[11] D.D. Muresan, T.W. Parks, Adaptive principal components and image denoising, Proceedings of IEEE International Conference on Image Processing (ICIP), 14-17 Sept. 1 (2003) I-101— I-104.

[12] Yu.A. Pichugin, Classification of summer

weather regions in St. Petersburg, Meteorology and Hydrology. 5 (2000) 31-39.

[13] M.S. Bartlett, The effect of standardization on a x2 approximation in factor analysis, Biometrika. 38 (3-4) (1951) 337-344.

[14] M.S. Bartlett, A note on the multiplying factor for various 2 approximations, J. Roy. Statist. Soc. B16 (1954) 296-298.

[15] D.N. Lawley, A.E. Maxwell, Factor analysis as a statistical method, Butterworths, London, 1963.

[16] S.A. Ayvazyan, V.M. Bukhshtaber, I.S. Enyukov, L.D. Meshalkin, Prikladnaya statistika. Klassifikatsiya i snizheniye razmernosti [Applied Statistics. Classification and dimension reduction], Finansy i statistika, Moscow, 1989.

[17] D. Jackson, Stopping rules in principal components analysis: A comparison of heuristical and statistical approaches, Ecology. 74 (8) (1993) 2204-2214.

[18] G.A.F. Seber, Linear regression analysis, John Wiley & Sons, New York, London, Sydney, Toronto (1977).

[19] Yu.A. Pichugin, The problem of statistical control of observation data on surface temperature at distant stations, Meteorology and Hydrology. 10 (2000) 18-24.

[20] Yu.A. Pichugin, Ekologicheskiy monitoring i metody mnogomernoy matematicheskoy statistiki [Environmental quality monitoring and multivariate mathematical statistics], Astrakhanskiy vestnik ekologicheskogo obrazovaniya. (2) (2012) 101-105.

[21] S.A. Ayvazyan, V.M. Bukhshtaber, I.S. Yenyukov, L.D. Meshalkin, Prikladnaya statistika. Issledovaniye zavisimostey [Applied Statistics. Relation studies]. Finansy i statistika, Moscow, 1985.

[22] Yu.A. Pichugin, Consideration of seasonal effects in problem of SAT forecasting and data control, Meteorology and Hydrology. 4 (1996) 52-64.

[23] A.S. Mikhalchuk, Yu.A. Pichugin, Dispe-rsionnyy analiz pogreshnostey tekhnologicheskikh

protsessov mikroelektroniki [The variance analysis of errors in the microelectronics technological processes], In collection of papers: Modelirovaniye i situatsionnoye upravleniye kachestvom slozhnykh sistem: sbornik dokladov [Simulation and quality control of complicated systems], SUAI, St. Petersburg (2017) 35-38.

[24] Yu.A. Pichugin, Empirical components of annual march of surface temperature, Meteorology and Hydrology. 12 (1994) 34-43.

[25] Yu.A. Pichugin, O.A. Malafeyev, Optimizatsiya i prognoz v dinamicheskoy modeli upravleniya portfelem tsennykh bumag [Optimization and prediction in the dynamic model of investment

Received 26.03.2018, accepted 21.06.2018.

portfolio governance], In: Materialy sektsionnykh zasedaniy simpoziuma «Nobelevskiye laureaty po ekonomike i rossiyskiye ekonomicheskiye shkoly» [In: Proceedings of symp. "Nobel Prize winners in economics and Russian economic schools of sciences"] , SPbSU, St. Petersburg (2003) 183-185.

[28] Yu.A. Pichugin, Glavnyye komponenty mnogomernykh vremennykh ryadov: analiz i prognoz [Principal components of multivariate time series: Analysis and prediction], In: Collection of papers of "The 13th International Youth Scientific Conf. on Soft Computing", Vol. 1, The 1st Electrotechnical University «LETI», St. Petersburg (2010) 160-163.

THE AUTHOR

PICHUGIN Yury A.

Saint-Petersburg State University of Aerospace Instrumentation

61 Bolshaya Morskaya St., St. Petersburg, 190000, Russian Federation

yury-pichugin@mail.ru

© Санкт-Петербургский политехнический университет Петра Великого, 2018