Модели экспоненциального сглаживания Брауна, взвешенного среднего, Холта – Уинтерса и подобные

Модель экспоненциального сглаживания [14, 17, 18] обладает лучшей адаптивностью при моделировании, прогнозировании коротких нестационарных временных рядов по сравнению с другими известными методами. Из-за своей простоты реализации и удобства применения он широко используется при прогнозировании и моделировании процессов в радиоэлектронике [4, 7], системах диагностики и контроля, прогнозирования состояния различных технических объектов [27, 28].

Согласно метода экспоненциального сглаживания вводится понятие экспоненциальной средней p-го порядка для временного ряда [29,30]:

, (3.7)

,

, ,

где и – параметры сглаживания ( ); .

Итерационные формулы (3.7) можно записать через отсчеты временного ряда и экспоненциальных средних в виде:

, где t=2,..,M (3.8)

В качестве прогнозирующего выбирают полином p-го порядка, в частности, полином первого порядка [29, 30]:

, (3.9)

где , – коэффициенты полинома; – упреждение прогноза.

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

Определим основные формулы для адаптивного определения коэффициентов , при прогнозе, для этого в выражение (3.8) подставим (3.9):

(3.10)

(3.11)

Решив систему линейных уравнений (3.10), (3.11) относительно неизвестных коэффициентов , , получим [31]:

, (3.12)

. (3.13)

В большинстве источников по методам адаптивного прогнозирования и экспоненциального сглаживания [14, 17, 18, 29], судя по ссылкам и в первоисточнике (работа Р.Г. Брауна и Р.Ф. Мейера [30]), формулы (3.12) и (3.13) приводятся в упрощенном виде:

, (3.14)

. (3.15)

Такое упрощение формул (3.12) и (3.13) можно использовать, если прогнозируемый временной ряд содержит большое число отсчетов . В общем же случае формулы (3.14), (3.15) являются приближенными и легко показать, что они совпадают с точными формулами (3.12), (3.13) только в пределе:

и .

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

Использование новых формул (3.12), (3.13) при прогнозировании процессов позволяет говорить о новом уточненном методе экспоненциального сглаживания. Общая линейная прогнозирующая модель в этом случае примет вид [31]:

,

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

,

или минимизации иной метрики вектора .

Абсолютная погрешность прогнозирования на один шаг определена по формуле .

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

- монотонный характер изменения без резких кратковременных скачков;

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

- процесс изменения является нестационарным в статистическом смысле, а временной ряд процесса относится к коротким временным рядам с числом отсчетов до 100.

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

Недостатками этого метода являются:

- невозможность прямого учета в модели действия того или иного влияющего фактора (как, например, в ARX-модели);

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

Одна из попыток вывода общей аналитической формулы для прогнозирующих полиномов любого порядка метода экспоненциального сглаживания предпринята автором в работе [76].

Кроме перечисленных, к моделям временных рядов также относят модели Холта-Винтера, простейшего скользящего среднего и ряд других. Но они в той или иной степени являются упрощением или модификацией изложенных моделей, подробное их описание можно найти в литературе [14, 16 – 18, 32].

3.3. Спектральные статистические модели (канонические разложения, разложение Карунена – Лоэва, модели компонентного и факторного анализов)

При моделировании объектов электропотребления для соблюдения требуемой адекватности модели, возникает необходимость использования многомерных математических моделей, которые позволяют выявлять и учитывать явные и неявные взаимосвязи и закономерности, объективно существующие между параметрами – отсчетами графиков нагрузки [3, 23].

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

Сложные технические системы, в том числе объекты радиоэлектроники за определенные интервалы времени могут быть всесторонне охарактеризованы при помощи наборов параметров (признаков) или, например, векторов дискретных реализаций V1, V2,…, VN этих параметров.

При описании стохастической системы многомерной моделью используются методы многомерного статического анализа [12, 26], в которых набор параметров объекта для j-го момента времени представляется в виде координат (компонент) случайного вектора Vj=[V1j, V2j ,…, Vnj] T.

К основным задачам, решаемым многомерным статистическим анализом [12, 26] относятся: устранение коррелированности координат между собой векторов V1, V2,…, VN ; снижение размерности пространства параметров Ln, модели объекта, т.е. представление многомерного случайного процесса V×=[V1×,V2×,…,Vn×]T через иной многомерный случайный процесс F× =[f1×, f2×,…, fm×]T меньшей размерности, m<n.

Однако глубоко проработаны эти вопросы в многомерном анализе только для случая стационарных случайных многомерных процессов, имеющих многомерный нормальный (гауссовский) закон распределения плотности вероятностей [11,12, 26]:

, (3.16)

где – ковариационная матрица многомерного процесса V; – алгебраическое дополнение элемента матрицы . Наиболее важная особенность n-мерного нормального закона распределения состоит в том, что он полностью определяется средними значениями (математическими ожиданиями) отдельных координат и ковариационной матрицей . Другая важная особенность нормального закона распределения (3.16) заключается в существовании многомерной центральной предельной теоремы [11, 25], которая утверждает, что распределение векторной суммы большого числа взаимно независимых n-мерных случайных величин стремится при достаточно общих условиях к нормальному n-мерному закону.

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

. (3.17)

Этот факт (3.17) называют взаимной независимостью координат в вероятностном смысле. Таким образом, каждая координата случайного многомерного процесса F может моделироваться одна независимо от другой, в случае же процесса Vнезависимое моделирование координат в общем случае недопустимо.

Перечисленные задачи многомерного анализа, в частности, могут решаться с помощью канонического разложения многомерного случайного вектора [10, 33]. Случайный вектор V с коррелированными в общем случае координатами в реализациях V1, V2,…,VN и корреляционной матрицей

выражается через случайные вектора с некоррелированными координатами в реализациях следующим образом [10,33,34]:

(3.18)

где – вектор математического ожидания компонент вектора V; Ui – векторы координатных функций канонического разложения.

Преобразование (3.18) позволяет перейти от случайного многомерного процесса Vj с симметричной корреляционной матрицей KV к процессу Fj с диагональной корреляционной матрицей . Так как приведение симметричной положительно определенной матрицы KV к диагональной форме не является однозначным [35 – 37], то из этого следует, что для любого случайного вектора с конечным моментом второго порядка существует бесконечное множество канонических разложений [10, 33]. Если число координат вектора F равно n, то равенство (3.18) является точным и m = n [34].

Методика моделирования случайных многомерных процессов, подобная изложенной, широко используется в компонентном и факторном анализах [38 – 43]. В этих случаях при построении математических моделей объектов по экспериментальным данным Х учитывают ошибки измерений ξ:

(3.19)

Задача построения факторной модели (3.19) во многом подобна задаче канонического разложения (3.18). Однако в отличие от методов канонического разложения, когда известна корреляционная матрица KV моделируемого процесса, в данном случае по результатам наблюдения можно получить только оценку корреляционной матрицы KХ многомерного процесса X, содержащего ошибку. Если задаться, что случайные процессы Х и ξ некоррелированы, то корреляционную матрицу KХ можно представить следующим образом:

.

В случае, если задаться, что ошибки измерения некоррелированы, то матрица Кξ диагональная. Следовательно, недиагональные элементы матриц КХ, КV с одинаковыми индексами равны

.

Таким образом, недиагональные элементы матрицы КХ можно принять равными соответствующим элементам матрицы KV, но диагональные элементы матриц КХ и KV не равны из-за ошибки измерения. Отсюда задача построения факторной модели (3.19) сводится к нахождению канонического разложения вектора Vj при неизвестных диагональных элементах его корреляционной матрицы [38 – 40].

Существуют два основных подхода построения факторной модели многомерного случайного процесса. Первый (компонентный) подход заключается в том, что каноническое разложение многомерного процесса Х заменяется каноническим разложением процесса V. Определяются первые m<n слагаемых (главных компонент) в сумме из (3.19), а все оставшиеся n – m члены суммы, обеспечивающие точное воспроизведение, переносятся в ошибку ξ. В данном случае естественно, что компоненты изменившегося вектора ошибки ξ станут коррелированными. Однако выбор m слагаемых в (3.19) осуществляют так, что минимизируется евклидова норма ковариационной матрицы ошибки воспроизведения процесса, т.е.:

, где .

Данный подход реализован в компонентных методах моделирования, в частности, в методе главных компонент (МГК) [38 – 42], в дискретном разложении Карунена – Лоэва (ДРКЛ) [34, 43, 44]. Данный подход также используется как аппарат в задачах математической редукции измерений [77], т.е. математической обработки результатов измерений таким образом, чтобы свести к минимуму искажающее влияние помех и непосредственно процесса измерения или прибора (редукция – приближение измерений к объекту).

Другой подход, относящийся к методологии факторного анализа, состоит в замене неизвестных факторных дисперсий или общностей (диагональных элементов матрицы КV) произвольными или подобранными по некоторому методу положительными числами [38 – 42]. В этом случае m координатных функций в сумме модели (3.19) выбирается из условия минимизации ковариации в векторе ошибок [45]:

Оба этих подхода позволяют сократить размерность модели (3.19) или моделируемого пространства, что равнозначно сокращению числа координат получаемого многомерного процесса F. Число координат процесса F сокращают до m (m<n) и это сокращение, в том или ином смысле, минимально ухудшает точность моделирования многомерного процесса.

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

Положительным свойством компонентного подхода (МГК, ДРКЛ) является возможность оценки в процессе получения координатных функций вклада каждой из них в точность воспроизведения (моделирования) по (3.19) многомерного процесса V (оценка информативности координат) [45].

Компонентный подход может использоваться на практике при решении задач моделирования как в детерминированной, так и статистической постановках. В детерминированной постановке он позволяет получить наилучшую проекцию совокупности точек наблюдений в пространство меньшей размерности [34, 40]. В статистической постановке он дает возможность выделить обобщенные случайные величины, имеющие максимально возможную дисперсию [12, 41, 42]. Одним из преимуществ этого подхода является то, что он является линейным и аддитивным.

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

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

Широкое распространение факторные и компонентные подходы получили в теории распознавания образов, в которой они используется для сжатия признакового пространства, т.е. выбора небольшого числа наиболее информативных признаков, характеризующих исследуемый процесс или объект [21, 34].

При решении технико-экономических задач эти подходы применяются при анализе производственной деятельности предприятий или отраслей промышленности, а также для разработки математических моделей управления производственным процессом [17, 41, 42].

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

В химии этот подход применяется для описания наборов данных, относящихся к изучению состава и физико-химических свойств веществ и соединений, для анализа данных химии окружающей среды (экология) и в биомедицинских исследованиях [48 – 52].

Недостатком моделей многомерной статистики, компонентных и факторных моделей является то, что в общем случае они были теоретически обоснованы для случая многомерных стационарных процессов, имеющих нормальный многомерный закон распределения. Только в этом случае правомерно и обосновано применение ковариационных и корреляционных матриц, канонических разложений и т.п. [12, 25, 26, 38 – 42]. Иное применение данных моделей для других частных случаев нестационарных и нелинейных процессов, как правило, никак не обосновывается [78].

3.4. Обобщенный фильтр Винера (метод прогнозирования Заде – Рагаззини)

Помимо рассмотренных моделей для моделирования и прогнозирования нагрузки применяются методы оптимальной фильтрации [53 – 55], в частности модификации алгоритмов фильтрации Винера (метод Заде – Рагаззини) [56].

Данные алгоритмы фильтрации базируются на уравнении Винера – Хопфа [53, 54, 57]:

, (3.20)

где – взаимно-корреляционная функция входа xвх(t) и выхода y(t) объекта;

– автокорреляционная функция входного сигнала x(t). Интервал TR обычно выбирают как интервал затухания корреляционных функций, а интервал Tз выбирается как интервал затухания импульсной переходной функции (ИПФ) k(t).

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

В случае дискретных прогнозируемых сигналов, в частности электропотребления, уравнение Винера – Хопфа (3.20) перепишется в виде [57, 58]:

, . (3.21)

При этом используется дискретное представления корреляционных функций и ИПФ [9, 14, 57 – 59] или их оценка.

Уравнение Винера – Хопфа составляет один из важнейших результатов теории фильтрации Колмогорова – Винера [60]. Показано, что интегральное уравнение первого рода (3.20) или дискретное (3.21) позволяют определять оптимальные ИПФ фильтра, обеспечивающего воспроизведение полезного сигнала из некоррелированного шума с минимальной среднеквадратической ошибкой (СКВО) . Если при этом сам полезный сигнал и помеха представляют собой стационарные, эргодические центрированные случайные функции [59].

Однако те обстоятельства, что процессы xвх(t) и y(t) должны отвечать условиям стационарности, эргодичности и иметь нулевое математическое ожидание, резко уменьшают применимость выражений (3.20) и (3.21).

В ряде случаев остаточная составляющая может быть представлена в виде суммы [55, 56]: , где – стационарный процесс с нулевым средним и корреляционной функцией вида [55]:

(3.22)

а составляющая – полином первой степени с изменяющимися во времени коэффициентами: , . Положительная константа T определяет длину интервала предыстории, на котором функция может представляться в виде указанной суммы.

Для указанного частного случая и построен метод прогнозирования Заде – Рагаззини, который является обобщением классической теории Винера для нестационарных процессов с полиномиальными нестационарными компонентами. Прогнозирующая математическая модель оптимального среднеквадратичного линейного фильтра Заде – Рагаззини для оценивания будущих значений процесса, наблюдаемого на интервале , в непрерывной форме имеет вид:

, (3.23)

где λ – время упреждения; – ИПФ линейного фильтра.

Когда задано значение l – степень полинома , то определяется путем решения системы из l+2 уравнений [55, 56]:

(3.24)

Первое выражение в системе (3.24) является обобщением уравнения Винера – Хопфа, уравнение (3.23) можно интерпретировать как интеграл свертки [55, 57, 58, 61]. Константы имеют смысл множителей Лагранжа при решении оптимизационной задачи с ограничениями: найти ИПФ , минимизирующую дисперсию ошибки предсказания процесса и удовлетворяющей условию .

Модели Заде – Рагаззини присущи ошибки, обусловленные: неполным соответствием корреляционной функции выражению (3.22) у реального процесса электропотребления; нестационарные компоненты процесса электропотребления могут описываться полиномами лишь приближенно. Все перечисленное зачастую может привести к катастрофическому росту погрешности и неприменимости модели [55].

Использование уравнений Винера-Хопфа для определения ИПФ очень часто приводит к решению некорректной плохо обусловленной задачи [57]. Это означает, что при незначительных изменениях входных данных, получаемая в результате решения (3.22), (3.24) ИПФ может измениться сколь угодно сильно. Но так как входные данные неизбежно подвержены искажению в результате действия помех, то очень важно чтобы эта задача была корректной или сводилась к ней. Для регуляризации изложенной задачи используют методы на основе минимизации функций невязки и сглаживания решения [57, 62]. Это обстоятельство затрудняет прямое использование описанной модели.

Кроме того, в работах показано [6], что адаптивная регрессионная ARMA-модель, по сравнению с рассмотренной, в большинстве случаев дает лучший результат по погрешности в силу своей адаптивности, хотя в этой модели явно не учитывается нестационарность процесса, как в рассмотренной модели. Поэтому рассмотренная модель в настоящее время имеет ограниченное применение.

3.5. Прогнозирующий фильтр Калмана (фильтр Калмана – Бьюси)

В настоящее время все шире начинают применяться при моделировании и прогнозировании графиков нагрузки методы и модели адаптивной цифровой фильтрации, модели оценивания состояния, в частности, фильтрация Калмана [1, 2, 4, 15, 19, 63].

Алгоритм оптимальной фильтрации (идентификации) для динамической системы электропотребления представляет собой рекуррентную реализацию алгоритма МНК, оценивания вектора параметров модели по имеющейся предыстории (дискретный оптимальный фильтр Калмана-Бьюси). Рекуррентные алгоритмы позволяют адаптивно с минимальными вычислительными затратами обновлять параметры модели. Общий вид модели прогнозирующего фильтра Калмана может быть представлен [1, 5, 9, 14, 15] как:

(3.25)

, , (3.26)

где – вектор измерений входных и выходных величин (векторов) процесса в j – 1момент времени; – значение i-го выходного параметра процесса в j-й момент времени; – ошибка оценивания параметра , случайная стационарная величина, распределенная по нормальному закону. Моделируемое значение может представлять как остаточную составляющую графика нагрузки или непосредственно весь график нагрузки .

В качестве выходной величины в (3.26) могут рассматриваться либо значение электропотребления в j-й интервал времени в отдельном i-ом узле сети [1], либо точка j-го суточного графика нагрузки в i-й момент времени [2, 4, 5, 8].

В развернутой форме последнее рекуррентное выражение (3.25) записывают в виде [1, 5, 9, 14, 15]:

(3.27)

,

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

Данная рекуррентная модель Калмана, представленная в форме алгоритма оптимальной фильтрации (3.25), (3.26) имеет ряд существенных недостатков [9, 64, 65]. В частности, алгоритм рекуррентной идентификации является не асимптотически устойчивым по Ляпунову, это означает, что он становится неустойчивым при постоянстве измеренных сигналов , т.е. при неинформативных измерениях происходит накопление ошибок оценки параметров модели. Это обстоятельство определяет только эпизодическое применение данного алгоритма при возникновении в объекте динамических режимов.

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

Перечисленные проблемы оптимальной фильтрации характерны и для процессов в энергетических системах [1, 2, 5]. Определение параметров фильтра требует оценки статистических характеристик моделируемого сигнала и помехи. Кроме того, полагают, что эти характеристики неизменны для моделируемого процесса (что не всегда верно). Учет изменчивости характеристик процесса и устранение эффекта “расходимости” увеличивает число вычислений [5, 64, 65], во время идентификации и требует итерационной настройки фильтра. Это, в свою очередь, в некоторых случаях приводит к расходимости фильтра вследствие ошибок вычислений.

Другая особенность прогнозирующего фильтра Калмана (3.25) – (3.27) состоит в том, что при моделировании процесса электропотребления не реализуется концепция многомерного моделирования графиков нагрузки. В выражениях фильтра Калмана, если скалярная величина представляет значение процесса в j-й интервал времени в отдельном i-м узле сети, то реализуется моделирование и прогнозирование «продольного среза» изменения потребления в одном из узлов во времени, причем в модели не учитывается «поперечная» зависимость отдельных потребителей (узлов) сети. Если же скалярная величина представляет собой точку j-го суточного графика нагрузки в i-й момент времени, то при ее прогнозировании учитывается «поперечное» изменение нагрузки внутри суток, но не учитывается «продольная» зависимость изменения этой i-й точки от графика к графику.

В работе [66] фильтр Калмана сравнивается с ARIMA-моделью (в обобщенном виде (3.5) при прогнозировании процесса и показано, что результаты прогнозирования очень близки между собой, это позволяет говорить о взаимозаменяемости моделей. Кроме того, показано [79], что калмановская фильтрация повышает точность прогнозирования при использовании ее совместно с моделями типа экспоненциального сглаживания.

Все вышесказанное ограничивает применимость фильтра Калмана при прогнозировании нагрузки. Основное применение оптимальных фильтров рассмотренного типа – это выделение и моделирование полезного сигнала на уровне шумов (в том числе оценивание состояния в различных точках электрической сети в реальном масштабе времени [2, 5], первичная обработка быстроизменяющейся информации и т.п. [19, 65]), при этом процессы линейные и имеют гауссовский характер [65] с известными характеристиками, а входные возмущения и шумы не коррелированы между собой.