Анализ динамики поступлений в бюджет от московского рынка наружной рекламы
Рубрика: Математическое моделирование в экономике и управлении
Статья в выпуске: 4, 2013 года.
Бесплатный доступ
В работе представлена методика анализа, моделирования и прогнозирования нестационарных временных рядов (НВР) поступлений в бюджет г. Москвы от объектов наружной рекламы. Показана возможность оценки основных статистик НВР по малому объёму исходных данных. Проанализированы способы выбора параметров адаптивных прогнозных моделей НВР.
Похожие статьи в разделе Общее машиностроение. Ядерная технология. Электротехника. Технология машиностроения
Короткий адрес: https://sciup.org/148160166
IDS: 148160166 | УДК: 621.031:681.332.5
The analysis of the dynamics of revenues from the Moscow market of outdoor advertising
This paper presents a methodology for analysis, modeling and forecasting of nonstationary time series (NTS) revenues from Moscow of outdoor advertising. Shows the possibility of estimation for main statistics NTS in a small sample size. Analyzed the methods of selecting the parameters for adaptive predictive models of NTS.
Текст научной статьи Анализ динамики поступлений в бюджет от московского рынка наружной рекламы
нятия решений (СППР) [1]. Информационнопоисковые, оперативно-аналитические, а также интеллектуальные СППР обеспечивают, соответственно, возможность поиска необходимых данных, их группировку и обобщение, а также поиск функциональных и логических закономерностей в накопленных данных, построение моделей и правил, которые объясняют найденные закономерности или прогнозируют развитие некоторых процессов.
Подсистема хранения информации в СППР использует современные системы управления базами данных (СУБД) и концепцию хранилищ данных. Предложенная Эдгаром Коддом модель многомерных данных в виде так называемого информационного гиперкуба и сформулированные им правила формирования реляционных СУБД являются фундаментальной основой для проектирования указанной выше структурной компоненты СППР.
В данной работе исследуется многомерная модель данных в виде помесячных поступлений в бюджет г. Москвы от объектов наружной рекламы (ОНР) в период с 2008 по 2011 год.
Подсистема интеллектуального анализа данных использует современные методы и алгоритмы Data Mining. К типичным задачам этого раздела теории искусственного интеллекта относят такие, как формирование описательных и прогнозных моделей. Первые улучшают понимание человеком анализируемых данных с помощью поиска ассоциативных правил и алгоритмов кластеризации, т.е. группировки однородных объектов. Вторые призваны идентифицировать модели данных по набору обучающих пар «вход – выход» и выполнить прогнозирование результатов для новых наборов наблюдений. К такому виду относят задачи классификации объектов по их характеристикам и поиска регрессионных зависимостей.
В данной работе представлена методика анализа, моделирования и прогнозирования нестационарных временных рядов (НВР) помесячных поступлений в московский бюджет от ОНР по типам рекламы и городу в целом. Важной особенностью предлагаемой методики является ее адекватность экспериментальным измерениям в условиях малого объема выборки. Этот результат достигается за счет применения интеллектуальных моделей локально взвешенной полиномиальной регрессии и робастного оценивания основных статистик НВР.
Современная методология адекватного прогнозирования НВР основана на рациональном сочетании экспертных моделей данных, а также методов адаптивного анализа и моделирования временных рядов [2]. Экспертный подход к анализу исходной информации позволяет учесть ее неопределенность. Плодотворным в этом смысле является, на наш взгляд, применение систем, основанных на обработке знаний экспертов с помощью алгоритмов нечеткого логического вывода.
Важной особенностью адаптивного направления анализа НВР является учет ожиданий участников инвестиционного процесса, а также возможности возникновения и влияния случайных факторов. Исследование различных сценариев развития событий позволяет естественным образом учесть указанные выше особенности прогнозирования НВР. Существенное значение в этом смысле приобретает управление результатами прогноза с помощью параметров эконометрических моделей НВР. В данной работе представлен анализ различных сценариев прогноза поступлений в бюджет г. Москвы от ОНР на основе гибких адаптивных моделей динамики временных рядов.
Модель многомерных данных
Трехмерная модель исходных данных представляет собой информационный «куб» [1] со следующими измерениями:
-
1) «время» в виде последовательности месяцев n = 1,2,..., N ( N = 48 ) в период с 2008 по 2011 год включительно;
-
2) «округ» в виде упорядоченной последовательности административных округов г. Москвы m = 1,2,..., M ( M = 11 ) , представленных в таблице 1 по состоянию на 2011 год;
-
3) «тип ОНР» в виде упорядоченной последовательности типов наружной рекламы k = 1,2,..., K ( K = 30 ) , представленных в таблице 2.
Меры, расположенные на пересечении осей измерений, представляют собой отчисления в московский бюджет от ОНР в миллионах рублей. СУБД поддерживает стандартные операции среза, вращения, консолидации и детализации над многомерной моделью данных [1].
Таблица 1
Административные округа
|
n 2 |
Наименование |
|
1 |
Общегородской заказ (ОГЗ) |
|
2 |
Восточный административный округ (ВАО) |
|
3 |
Западный административный округ (ЗАО) |
|
4 |
Зеленоградский административный округ (ЗелАО) |
|
5 |
Северный административный округ (САО) |
|
6 |
Северо-восточный административный округ (СВАО) |
|
7 |
Северо-западный административный округ (СЗАО) |
|
8 |
Центральный административный округ (ЦАО) |
|
9 |
Южный административный округ (ЮАО) |
|
10 |
Юго-восточный административный округ (ЮВАО) |
|
11 |
Юго-западный административный округ (ЮЗАО) |
ВЕСТНИК 2013 № 4
В качестве примеров на рис. 1 и 2 представлены графики НВР отчислений в городской бюджет по месяцам. Рис. 1 иллюстрирует данные, полученные в результате консолидации по округам и среза для ОНР «Щит отдельно стоящий». Рис. 2 демонстрирует динамику поступлений, консолидированных как по округам, так и по типам рекламы.
ВЕСТНИК 2013 № 4
Рис. 1. Консолидированные по округам поступления от ОНР «Щит отдельно стоящий»
Рис. 2. Консолидированные поступления в бюджет от ОНР
Таблица 2
Типы наружной рекламы
|
n 3 |
Наименование |
|
1 |
Витраж |
|
2 |
Выносная щитовая конструкция (штендер) |
|
3 |
Информационная конструкция предприятия и организации по обслуживанию населения на здании |
|
4 |
Кронштейн на здании (площадью одной стороны более 1 кв. м и высотой более 1,5 м) |
|
5 |
Кронштейн на здании (площадью одной стороны менее 1 кв. м и высотой менее 1,5 м) |
|
6 |
Кронштейн на несветовой опоре |
|
7 |
Кронштейн на световой опоре |
|
8 |
Крышная установка |
|
9 |
Маркиза |
В процессе анализа данных и поиска решений часто возникает необходимость в построении зависимостей между различными показателями. Число таких показателей может варьироваться в широких пределах. Таблицы реляционных баз данных не удовлетворяют требованиям объединения, просмотра и анализа информации с точки зрения множественности измерений. Оперативный многомерный анализ выполняют OLAP (On-Line Analytical Processing) – системы, предложенные Эдгаром Коддом в 1993 году. Основное назначение OLAP-технологии состоит в поддержке произвольных запросов аналитиков. Цель OLAP-анализа – проверка возникающих гипотез. Пример оперативного анализа структуры поступлений в московский бюджет от ОНР иллюстрирует рис. 3.
Визуальный анализ диаграмм рассеяния помесячных поступлений, сформированный OLAP-системой для 30-ти типов рекламы, позволяет аналитику ранжировать рекламные объекты по уровню их доходности. Вполне очевидной является следующая последовательность шести ОНР, упорядоченных по убыванию значимых суммарных поступлений за период с 2008 по 2011 год:
Рис. 3. Диаграммы рассеяния поступлений в бюджет от ОНР
-
1) щит отдельно стоящий (27) –
6301,21481116 млн руб.;
-
2) транспарант-перетяжка (23) –
1404,36650349 млн руб.;
-
3) настенное панно площадью более 2 кв. м (10) – 421,79276631 млн руб.;
-
4) крышная установка (8) – 321,92105990 млн руб.;
-
5) кронштейн на световой опоре (7) – 228,86775148 млн руб;
-
6) тумба отдельно стоящая (24) –
177,54042206 млн руб.
Графики временных рядов поступлений в московский бюджет от этих типов рекламы представлены в логарифмическом масштабе на рис. 4. Поступления для первых трех типов (с номерами 27, 23 и 10 из таблицы 2) измерены по левой вертикальной шкале. Остальные – по правой шкале.
Рис. 4. Помесячные поступления от ОНР со значимой доходностью
Главная задача аналитика – это формирование гипотез на основе своих знаний и опыта. Гигабайты и терабайты накопленных данных содержат значительный объем «скрытых» знаний, которые человек не в состоянии исследовать самостоятельно. Существует высокая вероятность пропустить значимые гипотезы.
Методологию обнаружения «скрытых» знаний в больших объемах накопленной информации называют добычей данных (Data Mining). Data Mining – это алгоритмическое обнаружение в исходных данных «скрытых» знаний, которые ранее не были известны, нетривиальны, практически полезны и доступны для интерпретации человеком.
Рекуррентные статистики временного ряда
Рассмотрим одну из задач, решаемых с помощью подсистемы интеллектуального анализа данных СППР, а именно: задачу анализа и моделирования НВР x ( n ) , n = 1,2, . , N в виде консолидированных помесячных поступлений в бюджет г. Москвы от ОНР. При формировании эконометрической модели указанного временного ряда необходимо учитывать относительно малый объем выборки, содержащей N = 48 отсчетов, и наличие аномальных значений (см. рис. 2). Надежной методической основой для решения указанной выше задачи является применение рекуррентного сглаживания НВР [3] с помощью процедур локально взвешенной полиномиальной регрессии [4; 5]. В частности, для оценивания основных статистик НВР удобными, на наш взгляд, являются модели:
-
– ядерной робастной регрессии (Robust Locally Weighted Regression and Smoothing Scatterplots) [4];
-
– наименьших квадратов, взвешенных «расстоянием» (Distance Weighted Least Squares – DWLS) [5];
– экспоненциально взвешенной робастной регрессии (Exponential Weighted Regression – EWR) [6];
– комбинация моделей DWLS и EWR (DEWR) [3].
В алгоритме рекуррентного сглаживания параметрические модели начального ( i = 1 ) и остаточных ( i = 2,3,... ) трендов у, ( n - т ) для отсчета относительного времени т аппроксимируют локально постоянной y , ( n - т ) = с 0 ' ) ( т ) , линейной у , ( n - т ) = с ( 0 ' ) ( т ) + т с 1 ' ) ( т ) или квадратичной у , ( n - т ) = с 0 ' ) ( т ) + т с 1 ' ) ( т ) + т 2 с 2 ) ( т ) зависимостью. Относительное время τ будем вычислять по схеме прогнозирования «вперед»,
ВЕСТНИК 2013 № 4
ВЕСТНИК 2013 № 4
т.е. τ = n – t . Здесь t – временн о й отсчет в пределах интервала сглаживания n - s i N + 1 < t < n , для которого оценивают значение начального или остаточного тренда, а 0,2 < s i < 0,8 - доля наблюдений, по которым выполняют сглаживание [4] на i -ой итерации удаления трендов [3]. В этом случае 0 < т < s i N - 1 и начало отсчета т = 0 соответствуют правой границе t = n интервала сглаживания. Альтернативный выбор т = n - t - ( s i N - 1 ) Д отвечает «центрированной» оценке тренда, поскольку в этом случае |т |< ( s i N - 1 ) /2 и начало отсчета т = 0 соответствуют середине t = n - ( s i N - 1 ) Д интервала сглаживания.
В дальнейшем для определенности будем рассматривать квадратичную аппроксимацию, которая, в большинстве случаев, позволяет хорошо оценивать тренды, а также их темп изменения и кривизну, в особенности в окрестностях локальных экстремумов исходного и остаточного временных рядов.
Вектор Ci (т) = { c(0i) (т),c 1 ) (т),c2i) (т)} оптимальных коэффициентов полинома на i-ой итерации сглаживания выбирают из условия ми- нимума
Ci ( т )| opt = Arg min “ Z w 1 i ) ( т - j ) w 2 i ) ^ с , (т ) [ j = 0
^{ ri ( n -j)} ri2 ( n -j)' (1)
взвешенного квадрата текущего остатка [4; 5]
r ( n ) = Г - 1 ( n )- У, ( n ) , i = 1, 2, — .
Здесь J, = siN; w 1 i)(т) и w(21){ ri (n - т)} -локальный и робастный веса остатка [7]. Исходный временной ряд r0 ( n ) = x ( n) инициализирует алгоритм рекуррентного оценивания. Рациональными критериями завершения итераций являются сходимость средней энергии текущих трендов к нулевому уровню ei =
N
Z w 2 i ) { r ( n ) } y 2 ( n )
n = 1
N
N Z w (2 i){ ri ( n )}
n = 1
< £ 1 , i = 1, 2, ...
или незначительное изменение средней энергии остаточных трендов за одну итерацию
I ei+1 - ei |< £ 2, i = 2, 3, —, где £1 > 0 и £ 2 > 0 - выбранные пользовате- лем достаточно малые уровни значимости критериев. Критерий окончания итераций i < Imax страхует неудачный выбор допусков. Финаль- ные модели остаточного ряда и тренда рассчитывают по формулам
r ( n ) = г, ( n ) , y ( n ) = x ( n ) - r ( n ) , где I < I max - номер заключительной итерации.
Локальный вес w 1 i ) ( т - j ) текущего остатка r i ( n - j ) обратно пропорционален временному интервалу между отсчетами т и j относительного времени. Поэтому значимый вклад в оценку параметров C i ( т ) полиномиальной регрессии y i ( n - т ) вносят лишь те наблюдения r i - 1 ( n - j ) , для которых отсчеты j близки к моменту времени т по критерию веса w ( i ) ( т - j ) .
С этой точки зрения масштаб s i весовой функции удобно интерпретировать как параметр сглаживания данных. Вычислительные эксперименты показали, что для обеспечения быстрой сходимости процедуры рекуррентного оценива- ния в соответствии с указанными выше критериями текущий масштаб si рационально увели- чивать по мере увеличения номера итерации, т.е.
si-1 + As , si-1 < ^max max , si-1 max
I . (2)
Локальный вес в модели Lowess w(i)(т) = = Ker{(т) / si}! si определяется выбором положительной четной функции ядра Ker(и), заданной на стандартном интервале [–1; 1] и интегрируемой с единицей J Ker (и) du = 1. Модели наиболее популярных ядер представлены в монографии [8, с. 140]. В частности, Клевеланд (Cleveland) рекомендует применять трижды кубическое ядро [4]
Ker ( и ) = <
70 ( 1 -| и |3) 3^81, | и | < 1
, | и | > 1.
Начальный масштаб ядра s 1 выбирают с помощью процедуры скользящей проверки [8].
Локальный вес в модели DWLS удобно, на наш взгляд, вычислять с помощью формулы w 1 ) ( т ) = s i N Д v ( s i ) т 2 + s i N } , где v ( s ) - скорость убывания веса в виде полинома, степень и коэффициенты которого выбирают эмпирически. В наших вычислительных экспериментах выбран полином пятой степени следующего вида
v ( s ) = 3,9301 - 33,9621 s + 122,1792 s 2 -
- 221,3187 s 3 + 198,7453 s 4 - 70,3125 s 5.
График этого полинома представлен на рис. 5. В модели DWLS, как известно, отсутствует механизм робастного взвешивания остатков. Поэтому все веса w ( 2 i ) { r i ( n - j ) } , i = 1,2, _ , I полагают равными единице.
В алгоритмах, реализующих модели Lowess и DEWR, предусмотрен дополнительный (внут-
Рис. 5. Зависимость скорости убывания локального веса от параметра сглаживания модели DWLS
ренний) цикл робастного оценивания. Для фиксированной итерации i внешнего цикла локального взвешивания вес остатков w2) {ri ( n - j)} на первой итерации робастного сглаживания по- лагают равным единице, поскольку отсутствует информация об остатке ri ( n - j). На последующих итерациях вес обратно пропорционален остатку. Такая зависимость позволяет исключить влияние аномальных значений временного ряда на параметры C, (т) полиномиальной регрессии y , (n - т) [6, 7]. Важно также отметить, что на второй и последующих итерациях робастного оценивания реализуют процедуру скользящего контроля, т.е. из суммы квадратов остатков в выражении (1) для целевой функции исключают отсчет j = т . Сглаживание продолжают до тех пор, пока робастные веса не стабилизируются. На практике процесс сходимости ограничивают двумя-тремя итерациями.
В модели Lowess Клевеланд рекомендует применять робастное взвешивание остатков, предложенное Хубером (Huber). Весовая функция Хубера имеет вид
Г2 (n - j), |r( n - j )|> R, 1 R^ n - j )|-R,2, |r( n - j )|> R, , w 2){ r( n - j )} = *
где R , = 1,345 RMS, ; RMS , - робастная мера масштаба остаточного НВР на - ой итерации локально взвешенного сглаживания. В качестве такой меры рационально применять абсолютное медианное отклонение [9]. Различные варианты выбора робастных весовых функций представлены в работах [6; 7].
В комбинированной модели DEWR оптимальные параметры полиномиальной регрессии выбирают из условия минимума функции потерь [3]
Г
C , ( т )| opt = Arg imin i a - N1+ X ) ( т ) E w ( , ) ^
Ci(т) [
^(т-j)w21){r,(n-j)} -.
Здесь w 2 ) { г ( n - j ) } = exp { - 2 r , 2( n - j )/ j ( 2 a , 2 ( n - j ) ) } - робастная весовая функция Л.Д. Мешалкина [6]; a , ( n - j ) - экспоненциально взвешенная (ЭВ) оценка дисперсии остаточного НВР r , ( n - j ) на текущем временном интервале сглаживания 0 < j < s , N - 1; X > 0 -параметр устойчивости ЭВ-регрессии к аномальным значениям данных.
Решение задач квадратичной оптимизации (1) и (3) приводит к системам линейных и нелинейных уравнений, соответственно [3]. Канонический вид этих систем позволяет решать их методом последовательных приближений. В качестве начальных оценок для тренда и дисперсии исходного и остаточных НВР выбирают стандартные оценки, полученные с помощью процедуры простой скользящей средней
1Ji y,( n )= E r-i(n - j),
J , j = 0
1J aA n )= E rA n- j), ^1,2,-, I.
J j = 0
Важным этапом формирования оценок основных статистик НВР является обоснование выбора параметров локально взвешенной полиномиальной регрессии.
ВЕСТНИК 2013 № 4
Выбор параметров локально взвешенной полиномиальной регрессии
Как упоминалось выше, начальное значение параметра сглаживания s 1 обычно выбирают с помощью функции скользящей проверки (CrossValidation). Альтернативный критерий оптимальности параметра s1 , реализующий идею регуляризации модели локально взвешенной полиномиальной регрессии, основан на введении штрафа за ее сложность. В частности, рациональной, на наш взгляд, является целевая функция si| opt = Arg min (DRF);
DRF = DR + ^ DF , (4)
где
NN
DR =^E r 2 ( n ) ; DF =t^E(a 2 y ( n ) ) ;
N n = 1 N 2 n = 3
ВЕСТНИК 2013 № 4
γ > 0 – вес штрафа за кривизну финального тренда
A2 y ( n ) = y ( n)- 2 y ( n -1) + y ( n - 2), n = 3,4,..., N.
Уместно отметить, что целевая функция (4) аналогична критерию оптимальности цифрового фильтра Ходрика – Прескотта [3]. Кроме того, компоненту DF в выражении (4) удобно интерпретировать как дисперсию ошибки прогноза в среднем. Действительно, вторую разность A 2 y ( n ) финального тренда можно представить в виде ошибки
А 2 y ( n ) = £ ( n ) = y ( n ) — y 1 ( n — 1 ) наивного локально линейного прогноза
У 1 ( n - 1 ) = у ( n - 1 ) + 8 1 ( n - 1 ) тренда на один шаг времени. Здесь 8 1 ( n - 1 ) = = y ( n - 1 ) - y ( n - 2 ) - текущий темп измене -ния тренда. Иными словами, смысл критерия (4) состоит в выборе оптимального значения параметра s 1 модели локально взвешенной полиномиальной регрессии, обеспечивающего разумный компромисс между дисперсией DR ошибки r ( n ) аппроксимации данных и дисперсией DF ошибки 8 ( n ) локально линейного прогноза в среднем.
Рис. 6. Зависимость целевой функции от параметра сглаживания модели Lowess: 1 – дисперсия ошибки аппроксимации данных;
2 – дисперсия ошибки локально линейного прогноза тренда; 3 – целевая функция для γ = 10;
4 – целевая функция для γ = 50
Рис. 6 и 7 иллюстрируют зависимости целевой функции (4) от начального значения s 1 параметра сглаживания в рекуррентной процедуре оценки тренда консолидированных поступлений в московский бюджет от ОНР (см. рис. 2) для моделей Lowess и DWLS. Целевая функция, изме-
Рис. 7. Зависимость целевой функции от параметра сглаживания модели DWLS: 1 – дисперсия ошибки аппроксимации данных;
2 – дисперсия ошибки локально линейного прогноза тренда; 3 – целевая функция для γ = 100
ренная в масштабе левой вертикальной шкалы, соответствует различным значениям штрафа γ за ошибку прогноза.
В модели Lowess выполнялись три итерации робастного взвешивания. Параметры правила (2) управления локальным взвешиванием для указанных моделей сглаживания сведены в таблицы 3 и 4. Четвертый и пятый столбцы таблиц содержат значения средней энергии начального и финального остаточных трендов.
На рис. 8 представлены оценки трендов консолидированных поступлений в бюджет г. Москвы от ОНР, полученные с помощью моделей Lowess и DWLS с параметрами s 1 = 0,4 и s 1 = 0,25, соответственно.
Таблица 3
Параметры управления рекуррентным сглаживанием в модели Lowess для оценки тренда
|
s 1 |
S max |
I |
e 1 |
e i |
|
0,2 |
0,6 |
10 |
26,184581 |
0,038549 |
|
0,25 |
0,65 |
13 |
25,781354 |
0,062346 |
|
0,3 |
0,7 |
15 |
25,887734 |
0,079519 |
|
0,35 |
0,75 |
10 |
25,694735 |
0,086252 |
|
0,4 |
0,8 |
17 |
25,185234 |
0,111986 |
|
0,425 |
0,8 |
11 |
25,066667 |
0,125209 |
|
0,45 |
0,8 |
14 |
24,984593 |
0,123597 |
|
0,5 |
0,8 |
14 |
24,749617 |
0,131539 |
|
0,55 |
0,8 |
14 |
24,716589 |
0,144613 |
|
0,6 |
0,8 |
10 |
24,722549 |
0,161205 |
Таблица 4
Параметры управления рекуррентным сглаживанием в модели DWLS для оценки тренда
|
s 1 |
S max |
I |
e 1 |
eI |
|
0,2 |
0,6 |
20 |
27,857751 |
0,004018 |
|
0,25 |
0,65 |
20 |
27,867429 |
0,004717 |
|
0,3 |
0,7 |
20 |
27,881724 |
0,005683 |
|
0,35 |
0,75 |
20 |
27,890537 |
0,006215 |
|
0,4 |
0,8 |
20 |
27,896649 |
0,004679 |
|
0,45 |
0,8 |
20 |
27,902536 |
0,003619 |
|
0,5 |
0,8 |
20 |
27,909339 |
0,002367 |
|
0,55 |
0,8 |
20 |
27,917260 |
0,001814 |
|
0,6 |
0,8 |
20 |
27,924725 |
0,002382 |
Рис. 9. Зависимость целевой функции от параметра сглаживания модели Lowess для γ = 100
Рис. 8. Финальные тренды для различных моделей сглаживания:
1 – консолидированные поступления в бюджет от ОНР ;
2 – модель Lowess; 3 – модель DWLS;
4 – робастные веса Хубера
Очевидно, что экспертные ожидания тенденции развития московского рынка наружной рекламы плохо согласуются с результатами, формируемыми моделью Lowess. Такой результат обусловлен весьма малыми величинами робастных весов Хубера, измеренных на рис. 8 по правой вертикальной шкале, для 36-го, 42-го и 47-го отсчетов НВР. Иными словами, вес штрафа за ошибку прогноза в среднем 10 ≤ γ ≤ 50 дает несколько заниженное оптимальное значение s 1 = 0,4 параметра сглаживания для модели Lowess.
Достаточно хорошо согласованные оценки финального тренда обеспечивает выбор веса штрафа γ = 100 для модели Lowess и оптимальное значение s 1 = 0, 425 . Зависимость целевой функции (4) от параметра сглаживания s 1 для этого случая демонстрирует рис. 9.
Оценки трендов, соответствующие указанному выше выбору параметров моделей сглаживания Lowess и DWLS, представлены на рис. 10.
Рис. 10. Оптимальные оценки тренда для различных моделей сглаживания:
1 – консолидированные поступления в бюджет от ОНР;
2 – модель Lowess; 3 – модель DWLS;
4 – робастные веса Хубера
ВЕСТНИК 2013 № 4
Рекуррентные робастные оценки тренда
Наряду с моделью Lowess на практике широкое распространение получила модель EWR – экспоненциально взвешенной робастной регрессии Л.Д. Мешалкина. Ее преимущества обусловлены, прежде всего, гибким механизмом управления робастным взвешиванием и высоким прогнозным потенциалом [6].
Анализ зависимости целевой функции (4) от параметра λ устойчивости статистик к аномальным значениям НВР показывает (рис. 11), что в широком диапазоне величин 10 ≤ γ ≤ 1600
Рис. 11. Зависимость целевой функции от параметра λ устойчивости модели EWR: 1 – дисперсия ошибки аппроксимации данных;
2 – дисперсия ошибки локально линейного прогноза тренда
ВЕСТНИК 2013 № 4
штраф DRF ( λ ) практически не отличается от значений дисперсии DR ( 2 ) остаточного ряда. Поэтому для временного ряда консолидированных поступлений в московский бюджет от ОНР предпочтительным является величина 2 = 0,1.
Процесс стабилизации робастных весов w 2 i ) { r i ( n ) } , 1 - n - N по критерию расстояния в N -мерном пространстве Евклида для этого случая иллюстрирует рис. 12.
Высокий потенциал модели экспоненциально взвешенной робастной регрессии для прогнозирования тренда НВР демонстрирует рис. 13. В зависимости от выбора значения 2 темп изменения тренда за один шаг времени и робастные веса, измеренные в масштабе правой вертикальной шкалы, варьируются в широких диапазонах для правой границы сегмента данных.
Рис. 12. Сходимость процедуры робастного взвешивания модели EWR для λ = 0,1
Рис. 13. Оценки тренда для модели сглаживания EWR: 1 – консолидированные поступления в бюджет от ОНР;
2 – тренд для λ = 0,1; 3 – тренд для λ = 0,25;
4 – тренд для λ = 0,5
Иными словами, модель EWR обеспечивает возможность управления результатами прогноза в среднем с помощью адаптации параметра 2 по мере поступления новых данных. Очевидно, что комбинация локального и робастного взвешивания, реализованная в модели DEWR [3], и соответственно адаптация по параметрам s 1 и 2 может усиливать указанное выше свойство.
а) λ = 0,1
б) λ = 0,25
в) λ = 0,5
Рис. 14. Зависимость целевой функции от параметра сглаживания модели DEWR: 1 – дисперсия ошибки аппроксимации данных; 2 – дисперсия ошибки локально линейного прогноза тренда; 3 – целевая функция для γ = 50;
4 – целевая функция для γ = 100
от выбора параметра устойчивости статистик к аномальным значениям в диапазоне значений 0,1 < X < 0,5.
Рис. 16. Оценки остаточного НВР для модели DEWR с параметром сглаживания s 1 = 0,35:
1 – λ = 0,1; 2 – λ = 0,25; 3 – λ = 0,5
Рациональный выбор параметров s 1 и λ комбинированной модели сглаживания по критерию (4) обосновывает зависимости, представленные на рис. 14. Соответствующие оценки тренда консолидированных поступлений в бюджет г. Москвы от ОНР иллюстрирует рис. 15.
Рис. 15. Оценки тренда для модели DEWR с параметром сглаживания s 1 = 0,35:
1 – консолидированные поступления в бюджет от ОНР ;
2 – тренд для λ = 0,1; 3 – тренд для λ = 0,25;
4 – тренд для λ = 0,5
Наш выбор квазиоптимальных параметров модели DEWR состоит в следующем: s 1 = 0,35 и X = 0,1. В этом случае результаты моделирования хорошо согласуются с экспертными ожиданиями темпа изменения тренда для московского рынка наружной рекламы. Более того, оценки остаточного НВР r ( n ), 1 ≤ n ≤ N (рис. 16) для его правого сегмента практически не зависят
Рекуррентные робастные оценки волатильности
Масштаб данных рационально оценивать с помощью алгоритма рекуррентного сглаживания квадрата остаточного НВР r ( n ) [3]. В этой процедуре начальную ( i = 1) и остаточные ( i = 2,3,. . ) меры волатильности v 2 ( n ) удобно аппроксимировать следующими локальными регрессионными моделями: v 2 ( n - т ) = g (01 ) ( т ) — постоянной; v 2 ( n - т ) = g ( 0 1 ) ( т ) + т g 1 1 ) ( т ) - линейной; v 2 ( n - т ) = g 0* ) ( т ) + т д 1 1 ) ( т ) + т 2 g 2 1 ) ( т ) - квадратичной. Здесь 0 < т < J, - 1 - относительное время, вычисленное по схеме прогнозирования «вперед». В дальнейшем для определенности будем рассматривать квадратичную аппроксимацию.
Векгор G 1 (т ) = { g 0 1 ) т ) , g 1 1 ' ( г ) , g 2 1 ) ( т ) } оптимальных коэффициентов полинома на i -ой итерации сглаживания выбирают из условия минимума
Г J - 1
G 1 ( т ) I opt = Arg Im in 1 £ w ( 1 1 ) ( т - j ) w ( 2 1 ) ^ 5 1 (т ) [ j = 0
^ { R1 2 ( n - j ) } R 2 ( n - j ) ' взвешенного квадрата остатка
R 1 ( n ) = u 1 - 1 ( n )- v 1 2 ( n ) , i = 1, 2, - .
Здесь w 11 1 ) ( т ) и w ( 21 ) { R 2 ( n - т ) } - локальный и робастный веса. Квадрат остаточного НВР u 0( n ) = r 2( n ) инициализирует алгоритм рекуррентного оценивания. На каждой итерации локального сглаживания выполняют нормировку u1 ( n ) = u ^1 ( n )/ v 2 ( n ) , что обеспечивает схо-
ВЕСТНИК 2013 № 4
димость текущей меры волатильности v2 ( n ) , 1 ≤ n ≤ N к единичному уровню в метрике Евклида. Соответственно стандартные критерии завершения итераций имеют вид:
E i =
N 2
X w ( 2 ) { R i ( n ) }{ v ( n ) - 1 } n = 1
N
NЕ w2°{ Ri2 ( n)}
n = 1
< 5 1 ,
i = 1, 2, (5)
или E+1 - Ei |< 5 2, i = 2, 3, ..., где 51 > 0 и 52 > 0 - выбранные пользователем достаточно малые уровни значимости критериев. В отличие от модели тренда, финальную оценку меры волатильности НВР рассчитывают по мультипликативной формуле
I
v( n) = л/Пv2 ( n ).
V i = 1
Очевидно, что центрированная и нормиро- ванная структурные компоненты
- ( n) = r (n) / v ( n)
ВЕСТНИК 2013 № 4
представляет собой временной ряд, стационарный по критериям тренда и волатильности.
Оптимальное начальное значение параметра сглаживания s1 в рассмотренной выше рекуррентной процедуре оценки волатильности НВР будем выбирать в соответствии с критерием регуляризации, аналогичным (4). Компоненты штрафа DRF за сложность локально взвешенной полиномиальной регрессии в этом случае приоб- ретают вид
N
DR = ^ X- - 2 ( n ); DF =
N n = 1
Рис. 17. Зависимости компонентов штрафа от параметра сглаживания модели Lowess: 1 – дисперсия ошибки аппроксимации квадрата остаточного НВР; 2 – дисперсия ошибки локально линейного прогноза меры волатильности
Рис. 18. Зависимости штрафа от параметра сглаживания в модели Lowess для оценки меры волатильности: 1 – γ = 10; 2 – γ = 50; 3 – γ = 100
Зависимости этих составляющих и целевой функции от параметра сглаживания s 1 для модели Lowess представлены на рис. 17 и 18, соответственно.
Штраф DRF на рис. 18 для значения веса γ = 10 измерен по левой вертикальной шкале, а для весов γ = 50 и 100 – по правой шкале. В вычислительном эксперименте выполнялись три итерации робастного взвешивания Хубера. Параметры правила (2) управления локальным сглаживанием сведены в таблицу 5.
Таблица 5
Параметры управления рекуррентным сглаживанием в модели Lowess для оценки волатильности
N - 2
N
Z ( A 2 v ( n ) ) .
n = 3
|
s 1 |
S max |
I |
E 1 |
EI |
|
0,2 |
0,6 |
10 |
1367,6288 |
0,000953 |
|
0,25 |
0,65 |
13 |
1377,9393 |
0,000512 |
|
0,3 |
0,7 |
15 |
1377,4365 |
0,000570 |
|
0,35 |
0,75 |
10 |
1387,0383 |
0,000388 |
|
0,4 |
0,8 |
8 |
1375,1016 |
0,000521 |
|
0,45 |
0,8 |
22 |
1373,9722 |
0,001165 |
|
0,5 |
0,8 |
23 |
1370,7831 |
0,001605 |
|
0,55 |
0,8 |
20 |
1351,5053 |
0,001988 |
|
0,6 |
0,8 |
20 |
1323,5300 |
0,002489 |
Четвертый и пятый столбцы таблицы содержат значения критерия (5) для начальной v 2 ( n ) и финальной v 2 ( n ) остаточной меры волатильности. Процесс сходимости остаточных мер волатильности v i 2 ( n ) , 1 < n < N к единичному уровню в метрике Евклида по итерациям i алгоритма рекуррентного сглаживания демонстрирует рис. 19. Наибольшее абсолютное отклонение финального НВР v 2 ( n ) , 1 < n < N от единицы не превышает величин 0,007182 и 0,029608 для значений параметра сглаживания s 1 = 0,35 и 0,55.
а) s 1 = 0,35
Рис. 20. Меры волатильности при сглаживании с помощью модели Lowess: 1 – r ( n ) для s 1 = 0,425; 2 – ν ( n ) для s 1 = 0,35;
3 – ν ( n ) для s 1 = 0,55
а) s 1 = 0,55
Рис. 19. Сходимость остаточных мер волатильности к единичному уровню: 2222
V 2 ( П ) ; V 3 ( П ) ; J V 5 ( П ) ; V 7 ( n )
ВЕСТНИК 2013 № 4
Рис. 21. Оценки центрированного и нормированных рядов при сглаживании с помощью модели Lowess:
1 – s 1 = 0,35; 2 – s 1 = 0,55
Финальные оценки меры волатильности v ( n ) , 1 < n < N консолидированных поступлений в московский бюджет от ОНР для указанных выше параметров модели сглаживания Lowess представлены на рис. 20. Робастные веса Хубера измерены по правой вертикальной шкале. Соответствующие оценки центрированного и нормированного временных рядов z ( n ) , 1 < n < N иллюстрирует рис. 21. Более предпочтительным, на наш взгляд, является выбор s 1 = 0,425 и 0,35 при анализе тренда и масштаба данных.
Альтернативные DEWR-оценки меры волатильности консолидированных поступлений в бюджет г. Москвы от ОНР демонстрирует рис. 22. В вычислительном эксперименте выполнялось рекуррентное сглаживание квадрата центрированного ряда r ( n ) , 1 < n < N , полученного с помощью модели DEWR с параметрами s 1 = 0,35 и X = 0,1 (рис. 16). Предпочтительными, по нашему мнению, являются следующие параме-
Рис. 22. Меры волатильности при сглаживании с помощью модели DEWR:
1 – r ( n ) для s 1 = 0,35 и λ = 0,1; 2 – ν ( n ) для λ = 0,1;
3 – ν ( n ) для λ = 0,25
тры управления локальным сглаживанием ряда r 2 ( n ) : 5 1 = 0,35, A 5 = 0,05; S max = 0,75, I = 48. Наибольшее количество итераций робастного взвешивания Л.Д. Мешалкина не превышало 10.
Сходимость по критерию (5) текущей меры волатильности v 2 ( n ) , 1 < n < N к единичному уровню в метрике Евклида за двадцать итераций 2 < i < 20 иллюстрирует рис. 23. Наибольшее абсолютное отклонение финального НВР v 4 8 ( n ) , 1 < n < N от единицы не превышает величин 0,006273 и 0,000856 для значений параметра устойчивости статистик X 1 = 0,1 и 0,25. Соответствующие DEWR-оценки временного ряда z ( n ), 1 < n < N представлены на рис. 24. В итоге приемлемый, на наш взгляд, выбор параметров DEWR модели сглаживания при анализе меры волатильности консолидированных поступлений в московский бюджет от ОНР состоит в следующем: 5 1 = 0,35 и X = 0,1.
ВЕСТНИК 2013 № 4
Рис. 23. Сходимость рекуррентной процедуры оценки волатильности НВР с помощью модели DEWR для s 1 = 0,35 и λ = 0,1
Рис. 24. Оценки центрированного и нормированного рядов z ( n ) при сглаживании с помощью модели DEWR для s 1 = 0,35: 1 – λ = 0,1; 2 – λ = 0,25
Кластеризация доходности ОНР
Дальнейшее исследование динамики поступлений от рекламы связано с анализом статистик центрированного и нормированного временных рядов z ( n ) , 1 < n < N . В частности, практический интерес представляет изучение корреляционных зависимостей этого ряда. Анализ формы диаграмм рассеяния для значений корреляционных лагов от 1 до 6 показал, что данные не обладают сколько-нибудь выраженным генеральным направлением. Указанное свойство иллюстрирует рис. 25, на котором представлены двумерная и трехмерная диаграммы для лагов двух и шести месяцев. Стандартные выборочные оценки коэффициента корреляции для этих лагов составили величины - 0,229 ± 0,138 и 0,207 ± 0,132, соответственно.
а)
б)
Рис. 25. Диаграммы рассеяния временного ряда z ( n ): а) лаг два месяца; б) лаги два и шесть месяцев
Иными словами, отсчеты временного ряда z ( n ) , 1 < n < N практически не коррелированы. Обнаруженное свойство позволяет рассчитать доверительные интервалы и соответствующие им вероятности пребывания ряда в заданном диапазоне его значений. Удобной методической основой такого анализа является гистограмма, сглаженная сдвигом (Average Shifted Histogram – ASH) [8] 1 m - 1
*-)= Ns , 5/ (M-;
z k = z min + k A z ; k = 0, •••, K ;
K = ( z max - z min )/ A z •
Здесь 5 = 2 IQ J 3[N - робастная оценка ши- рины разрядных интервалов (bins) Фридмана –
Дьякониса [8]; IQ – интерквартильный диапазон данных; Az = 5 / m и m - ширина «суженных» интервалов (narrow bins) и их количество; Sl -количество наблюдений, попавших в l-ый «су- женный» интервал (Sl = 0, если l < 0 или l > K).
Окно данных w( j ) выбирают из условия ^mmw( j ) = m. В этом случае гистограмма интегрируема с единицей. Такой нормировке удовлетворяет обобщенное окно вида m-1
w ( j ) = m Ker ( j / m ) ^ Ker ( l / m ) .
l = 1 - m
Расчеты выполнялись для центрированных и нормированных временных рядов, полученных с помощью моделей Lowess (рис. 21) и DEWR (рис. 24) с параметрами s 1 = 0,35 и X = 0,1. Пунктирными линиями на указанных рисунках отмечены уровни поддержки ZS и сопротивления ZR , выбранные по критерию робастных весов. В вычислительном эксперименте гистограммы сглаживались окном данных, сформированным на основе трижды взвешенного ядра Епанечникова [8]. Оценки распределений представлены на рис. 26.
б) модель сглаживания DEWR
Рис. 26. ASH-оценки распределений отсчетов временного ряда z ( n ):
1 – плотность распределения вероятностей (ПРВ);
2 – интегральное распределение
= P { z < Z S } и pR = P { z > Z R } выбросов ряда z ( n ) , 1 < n < N за пределы уровней поддержки и сопротивления сведены в таблицу 6. Последний столбец таблицы содержит вероятности
pSR = P { ZS < z < ZR } = 1 — pS - pR для «типичных» значений доходности ОНР.
Таблица 6
Уровни поддержки и сопротивления центрированных и нормированных консолидированных поступлений от ОНР
|
Модель |
Z S |
Z R |
pS |
pR |
p SR |
|
Lowess |
–1,208046 |
2,231238 |
0,0233 |
0,059472 |
0,917228 |
|
DEWR |
–1,214729 |
2,265476 |
0,0046 |
0,026147 |
0,969253 |
ВЕСТНИК 2013 № 4
Важно также отметить, что ASH-оценки демонстрируют наличие кластерной структуры данных, а именно наличие трех классов «низкого» – Low, «среднего» – Mid и «высокого» – High уровней консолидированных поступлений от ОНР. Дискриминантные границы этих классов доходности ОНР, сформированные в виде абсцисс ZLM и ZMH локальных минимумов гистограмм, а также вероятности p L = P { z < ZLM } , p M = P { ZLM < z < ZMH } , pH = P { z > ZMH } таких доходностей сведены в таблицу 7.
Таблица 7
Границы классов доходности ОНР
|
Модель |
ZLM |
7 MH |
pL |
pM |
pH |
|
Lowess |
0,653179 |
1,967250 |
0,604889 |
0,313794 |
0,081317 |
|
DEWR |
0,679055 |
1,533509 |
0,679603 |
0,171313 |
0,149084 |
Рассмотренная выше методика интеллектуального анализа НВР позволяет формировать прогнозные модели доходности московского рынка наружной рекламы.
ВЕСТНИК 2013 № 4
Адаптивная модель прогнозирования доходности ОНР
Будем предполагать, что изменения основных статистик для поступлений от ОНР г. Москвы достаточно медленные. Такого рода персистентность данных позволяет описывать тренд y ( n) и меру волатильности v( n ) для доходности ОНР с помощью адаптивной модели Хольта (Holt) [10]. В рамках этой модели локально линейные прогнозы y 1 ( n -1) = a 1 ( n -1) + a 2 ( n -1), v1 ( n -1) = b 1 ( n -1) + b 2 ( n -1) указанных выше статистик на один шаг времени рационально рассчитывать в терминах ошибок аппроксимации
n ( n ) = y ( n )- y 1 ( n - 1 ) ,
д( n ) = v (n)-V1 ( n -1), 2 < n < N.
Адаптацию средних значений a 1 ( n ) , b 1 ( n ) и темпов изменения a 2 ( n ) , b 2 ( n ) статистик по мере поступления новых оценок y ( n ) и v ( n ) реализует прогнозная форма процедуры экспоненциального сглаживания
a 1 (n) = a 1 (n -1) + a 2 ( n -1) + а1 n( n), a 2 (n) = a 2 ( n -1) + а1 а 2 n( n),
b 1 ( n) = b 1 (n -1) + b 2 ( n -1) + в1 Ц (n),
b2 ( n) = b2 ( n -1) + в1 в2 ц( n).
Здесь а 1 , в 1 и а 2 , в 2 - параметры сглаживания статистик и их первых разностей. Эти параметры удовлетворяют известным ограничениям 0 < а 1 < а 2 < 1 и 0 < в 1 < в 2 < 1 [10]. Относи-
а) модель сглаживания Lowess
б) модель сглаживания DEWR
Рис. 27. Относительные ошибки аппроксимации статистик: 1 – тренда; 2 – меры волатильности тельные ошибки аппроксимации
Et ( n) = 100 n( n)/y ( n),
Ev ( n ) = 100 ^v ( n), 2 < n < N тренда и меры волатильности, полученные с помощью процедур рекуррентного сглаживания Lowess и DEWR, демонстрирует рис. 27.
Сравнение рисунков свидетельствует о том, что модель Хольта аппроксимирует DEWR-оценки статистик с меньшими ошибками. Соответствующие оптимальные значения параметров модели Хольта сведены в таблицы 8 и 9.
Таблица 8
Параметры модели Хольта для аппроксимации тренда
|
Модель |
a 1(1) |
a 2(1) |
α 1 |
α 2 |
|
Lowess |
210,1 |
1,025 |
0,736 |
1 |
|
DEWR |
205,3 |
4,495 |
0,772 |
1 |
Таблица 9
Параметры модели Хольта для аппроксимации волатильности
|
Модель |
b 1(1) |
b 2(1) |
β 1 |
β 2 |
|
Lowess |
20,61 |
9,712 |
0,681 |
1 |
|
DEWR |
53,81 |
3,173 |
0,396 |
1 |
Актуальные линейные прогнозы тренда и меры волатильности НВР на заданное количество m шагов времени приобретают вид
ym ( N ) = a 1 ( N) + ma 2 ( N) ,
vm ( N )= b 1 ( N ) + mb 2 ( N) .
Соответствующие прогнозы
x m ( n ), x m ( n ), x mM ( n ), x MH ( n )
характерных уровней доходности ОНР вычисляют по формуле
x (Л N ) = У m ( N ) + Z ( . ) v m ( N ) , где под символом ( * ) подразумевают один из четырех индексов S , R , LM , MH уровней доходности ОНР, представленных в таблицах 6 и 7.
Важно отметить, что каждому уровню x m ( N ) , x m ( N ) или x mM ( N ) , x MH ( N ) отвечают свои вероятности pS , pR , pSR или pL , pM , pH попадания доходности в соответствующие доверительные интервалы. В этом смысле прогнозные оценки уместно интерпретировать как интервальные. Результаты такого рода прогнозирования на три месяца ( m = 3) поступлений в бюджет г. Москвы от ОНР с помощью моделей Lowess и DEWR иллюстрируют рис. 28.
Сравнение рисунков свидетельствует о том, что модель Lowess формирует несколько завышенные оценки тренда и меры волатильности для правого сегмента данных.
а) модель сглаживания Lowess
б) модель сглаживания DEWR
Рис. 28. Прогноз рынка наружной рекламы:
1 – консолидированные поступления в бюджет; 2 – тренд;
3 – граница «низкой/средней» доходности;
4 – граница «средней/высокой» доходности;
5 – линии поддержки и сопротивления
Заключение
В работе представлена методика интеллектуального анализа данных на примере исследования основных статистик и прогнозирования НВР в виде поступлений в бюджет г. Москвы от ОНР. Продемонстрирована возможность и эффективность оценивания тренда и меры волатильности рынка наружной рекламы с помощью моделей локально взвешенной робастной регрессии и процедуры рекуррентного сглаживания выборочных данных малого объема. Исходя из гипотезы медленного изменения статистик НВР, сформированы интервальные краткосрочные прогнозы динамики поступлений в рамках адаптивных эконометрических методов.
Список литературы Анализ динамики поступлений в бюджет от московского рынка наружной рекламы
- Барсегян А.А., М.С. Куприянов, В.В. Степаненко, И.И. Холод. Технологии анализа данных: Data Mining, Visual Mining, Text Mining, OLAP. -2-е изд., перераб. и доп. -СПб.: БХВ -Петербург, 2007. -384 с.
- Лукашин Ю.П., Рахлин Л.И. Современные направления статистического анализа взаимосвязей и зависимостей/отв. ред. Ю. П. Лукашин. -М.: ИМЭМО РАН, 2012. -54 с.
- Лабунец Л.В., Лабунец Н.Л., Чижов М.Ю. Рекуррентные статистики нестационарных временных рядов//Радиотехника и электроника. -2011. -Т. 56. -№ 12. -С. 1468-1489.
- Cleveland W.S. Robust Locally Weighted Regression and Smoothing Scatterplots//Journal of the American Statistical Association. -1979, V. 74, N 368, P. 829-836.
- McLain D.H. Drawing contours from arbitrary data points//The Computer Journal. -1974, V. 17, N 4, P. 318-324.
- Шурыгин А.М. Прикладная стохастика: робастность, оценивание, прогноз. -М.: Финансы и статистика, 2000. -224 с.
- Воронцов К.В. Лекции по алгоритмам восстановления регрессии. -2009.
- Scott D.W. Multivariate Density Estimation: Theory, Practice, and Visualization. -N.-Y.: John Wiley & Sons, Inc, 1992. -317 c.
- Шестаков О.В. О скорости сходимости распределения выборочного абсолютного медианного отклонения к нормальному закону//Информатика и её применения. -2011. -Т. 5. -Вып. 5. -С. 74-79.
- Лукашин Ю.П. Адаптивные методы краткосрочного прогнозирования временных рядов: учебное пособие. -М.: Финансы и статистика, 2003. -416 с.