Общий подход к идентификации экономической динамики на основе моделей авторегрессии Текст научной статьи по специальности «Математика»

УДК 519.6:311

ОБЩИЙ ПОДХОД К ИДЕНТИФИКАЦИИ ЭКОНОМИЧЕСКОЙ ДИНАМИКИ НА ОСНОВЕ МОДЕЛЕЙ АВТОРЕГРЕССИИ

© 2004 В. К. Семёнычев

Самарский государственный аэрокосмический университет

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

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

В отсчетах Ук(Ук= У(Тк), Тк = кА -моменты получения с периодом А отсчетов, (к = 0, 1, 2,..., Ы) анализируемых временных рядов вместе с неслучайными компонентами неизбежно присутствует эволюционный стохастический компонент ^к, общим предположением о котором является нормальность закона его распределения Ыо, нулевое математическое ожидание, постоянство дисперсии, некоррелированность отсчетов [5, 6]:

^к = Ыо (0, а), Б[^к] = а2, М[£к£к + 1] = 0,

(1)

где М - оператор математического ожидания, Б - оператор дисперсии, а - среднеквадрати-ческое отклонение Ъ,к, N - объем выборки временного ряда.

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

чительными погрешностями или даже принципиально невозможна.

Например, при широко распространенном в практике экономических исследований случае моделирования временного ряда экс-понентой

Ук=А1 ехр( - а1 Тк) (2)

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

Помеха при предположении (1) центрирована и, в силу этого, знакопеременна. Более того, нелинейное логарифмическое преобразование отсчетов Е, к нарушает предположение (1), и среднеквадратические оценки параметров Л1 и а1 модели (2) окажутся в результате этого смещенными и неэффективными [3 - 5]. Указанные обстоятельства зачастую не принимают в расчет.

Появление слагаемого А0 в модели (2) для обобщенной экспоненциальной функции

Ук = А0 + А1ехр( - а1 Тк) (3)

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

Известные способы оценки параметров гармоники (не сводимой каким-либо преоб-

разованием к линеинои по параметрам модели), обычно описывающей сезонные и циклические компоненты

YK = A1Cos(& Тк + ф)

(4)

ражении соотношением YZ - ^Z"1 Y(Z) = Ао,

(5)

где Х1=ехр(-а\А), 2 = ехр(р) ,р - комплексная переменная.

Вернувшись в область оригиналов, будем иметь

первого порядка YK= X.YK- 1 + К

(6)

сложны, предполагают анализ десятков и сотен отсчетов на длительности более одного периода гармоники [5, 6].

Еще сложнее осуществляется идентификация известными методами сочетаний полиномиальных, экспоненциальных, гиперболических и тригонометрических функций при аддитивной помехе ^к[3 - 5].

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

2. Демонстрация подхода на экспоненциальной модели. Для выражения (2) выполним 2-преобразование (его называют также преобразованием Лорана) [6, 7]

да А0 У{2)= Е Ао ехр( - а\Ак) * = - ,

к=0 \-X\Z'1

которое будет представлено в области изоб-

Видим, что из (6) по любым двум отсчетам Кки Ук- 1 (т. е. на участке в 10^15% от длительности переходного процесса модели (2), оцениваемого обычно « 3 < а1) можно определить параметр а затем и а1 по формуле

aj= - 1/А LnXy

(7)

При наличии более двух отсчетов можно оценить параметр А по методу наименьших квадратов, реализуя необходимое усло -вие экстремума:

N

А°1 = arg min Х{ Yk - А Yk-1 }2,

Ai к=1

которое приводит к следующему действию над отсчетами

N

N

А°1 = Z{ YKYK-i}{Z{ Y Vi}}"1

к=1

к=1

(8)

Можно утверждать, что полученные оценки (и а°1) при выполнении предположения (1) являются оптимальными, т. к. использованное для их определения соотношение (6) является линейным по отношению к отсчетам Yк, Yк-1 и £к.

Параметр А можно также оценить оптимально в среднеквадратическом смысле, обеспечив выполнение на линейной по отношению к отсчетам Yк выборке условие

N

А°1 = arg min !{ Yk- (A°i)k А1}2 А1 к=0

Yk=A1 Yk-1 + Аобк, где

Yk = 0, при к < 0, Sk =

1, к = 0 - дискретный

аналог дельта -

функции 0, к Ф 0.

При k > 1 и аддитивной помехе ^кдля

которое ведет к следующей операции над N + 1 отсчетами (при рассчитанном ранее значении А,°):

N

N

А°1= EYK(A°i)K{E(A°i)2K}-1.

к=0 к=0

Значения a° и Л°, можно использо-

(2) будет справедлива модель авторегрессии вать для расчета помехозащищенных «состо-

явшихся» или, что очень важно, «будущих» прогнозных значений У к при тех или иных значениях «к», подставляя их в модель (2).

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

3. Идентификация обобщенной экспоненциальной функции. Аналогично можно получить и для модели (3) (при к >2 и аддитивной помехе) модель авторегрессии второго порядка, записанную в первых разностях отсчетов:

AYK = \AYK- 1 +

(9)

где АУк = Ук- Ук-\, АУк- 1 = Ук-1 - Ук-2, параметр ^ определен в (5).

Можно рассчитать параметры модели (3) из условий

годовом сезонном цикле. При большем объеме используемой выборки можно осуществить среднеквадратическое приближение

N

v°1 = arg min Z{ Yk - v1 Yk-1 + Yk-2}2 . (13) v1 k=2

Тогда помехоустойчивая оценка частоты из (12) и (13) будет равна

ш° = 1/Д Arccos(v°1/ 2),

а оценки амплитуды и фазы определяются условием

N

A°2, A°3 =arg min !{ Yk - A2Cos®°kA+ A2,A3 k=0

+ A3Sinra °kA}2

и соотношениями

а°1 = (а°22 + а°з2)1/2,

N

X°i = arg min E(AYk- ^iAYk-1}2, (10)

N

A°o, A\° = arg min !{ Yk- Ao - (X°{f A1}2.

ao,ai k=0

(11)

Условие (10) приводит к соотношению (8), но в первых разностях отсчетов, а (11) -к соответствующей системе двух линейных алгебраических уравнений (СЛАУ второго порядка) относительно А°0 и А° Минимально необходимое число отсчетов для идентификации модели (3) равно трем.

4. Идентификация гармонической компоненты. Для модели (4) при к > 2, у1 = 2Cos&А и аддитивной помехе будет справедлива авторегрессия

Yk= v1 Yk-1 -Yk-2 + £к.

(12)

Значение у1 можно определить из (12) через любые три отсчета: Ук, Ук-\ и Ук-2, расположенных на доле периода гармоники: например, через три месячных отсчета на

ф = Arctg(A°3/A°2),

в которых введение новых параметров А2 = = А1Со^ф и А3 = А^тф демонстрирует прием определения ф.

5. Идентификация многокомпонентных моделей временных рядов.

5.1. К ним относится, например, сумма линейного тренда и гармоники (аддитивной сезонной или циклической компоненты) [1 -

3]:

Yk = А1Тк + А2 + A3Cos(& Тк + ф), (14)

которой будет соответствовать при k > 4 и аддитивной помехе следующая модель авторегрессии:

YK= (YK-1 + YK-3 - 2YK-2)| + + 2 YK-2 - YK-4 + (15)

где | = 2(Cos&A + 1).

Очевидны среднеквадратическая оценка из (15)

N

|°= arg minL{ Yk-|( Yk-1+Yk-3 - 2 Yk-2)-

| K=4 - 2 YK-2 + YK-4}2

и последующим расчет через частоты

ш°=1/ AArcCos(^°/2 - 1).

Оставшиеся параметры модели (14) определят условия

N

А°рЛ°2,Л°4,Л°5= Mg min !{ Yk - А1кА -А2-

А1,А2,А4,А5 к=0

- А4^ш°кА + А^шш°кА }2,

где Аг А2 - параметры модели (14), А4 = = А3Cosф и А5 = А^шф. После решения соответствующей СЛАУ четвертого порядка относительно А° А°2, А°4, Л°5 вычислим

А°з = (А°42 + А°52)1/2,

ф° = Arctg(А°5/А°4).

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

5.2. Часто требуется идентифицировать линейный тренд с мультипликативной сезонной или циклической компонентой [1 - 3]

Yk = (А1Тк+ А2) Cos(o Тк+ ф).

N

{А°у}= arg min Z{ Yk - А3кАcosш°Ак -А]- к=0

-АлCosш°Ак + А5кАSinш°Ак + А^тш°Ак}2,

где] = 3,4,5,6; Л3= Л1^'ф, Л4 = А2Cosф, Л5 = = Л1Sinф, Л6 = Л^шф.

Тогда

Л°1 = (Л°32 + Л°5)1/2, Л°2 = (Л°42 + Л°62)1/2, ф° = °б/Л°4).

Параметры тренда Л1, Л2 и параметр ф мультипликативной компоненты определяет СЛАУ четвертого порядка, а частоту мультипликативной компоненты - СЛАУ второго порядка. Минимальное число отсчетов для идентификации модели (16) равно пяти.

5.3. Известно [1, 4, 8] моделирование экономических процессов выражением

Yк=А1ехр(а1 Тк) +

+А2ехр(-а2 TK)Cos(ш Тк +

(18)

(16)

Поставим в соответствие модели (16) авторегрессию отсчетов при k > 4 при аддитивной помехе

Yk=Vj( Yk-1 + Yk-3) - v2 Yk-2 + Yk-4 + ¿¿к, (17)

где v1 = 4CosшА, v2 = 2 + v^/4.

Найдем из (16) и (17) помехозащищен-ные оценки параметров гармоники и линейного тренда:

N

v°1,v°2 = arg min 2{ Yk - vi( Yk-1 + Yk-3) +

v1,v2 k=4

+v2 Yk-2 - Yk-4}2,

ш° = (1/А)Агс Cos(v°1/4),

которому при k> 3 и аддитивной помехе, как показано в [8], соответствует уравнение регрессии третьего порядка

Yk= ^Yk-1 - Yk-2 + ^3Yk-3 + ¿к. (19)

Из (19) следуют среднеквадратические приближения (при числе отсчетов больше шести)

N

Я,°19 = arg min E{Yk- ^Yk-1 +

2Д з к=3

+l2 Yk-2 - Yk-3}2,

где \ = 2V3V2 + Vl2' ^2 = V2 + 2V3V2V12,

= v2vb vl = ехр(а1А), v2 = ex^( -а2А), v3 = CosшА.

Решив соответствующую СЛАУ третьего порядка относительно A,°l? ^°3, вычислим последовательно v°l, v°2 = (Ä,°3/ v°l), v°3 = (^°2 - v°2)/ (2v°lv°2).

Искомые динамические параметры модели (18) с учетом обозначений в (19) могут быть рассчитаны следующим образом:

a°i = 1/ÄLnv°i, а°2 = ( - 1/Ä)Lnv°2, ю°= (1/Ä)ArcCosv°3,

N

А°1, А°3, А°4 = arg min £{ Yk - А1ехр(а°1Ак) -AiA з,А 4 k=0

- А3ехр( - a°2Äk)Cosro°Äk+

+ А4ехр( - a°2Äk)Sinro°Äk}2,

где А3 = А2^5ф, А4 = А^шф,

А°2 = (А°32 + А°42)1/4, ф° = Arctg(Ä°jÄ°).

5.4. Известны [1, 4 - 6, 8] примеры моделирования экономической динамики рядами отсчетов

YK=Axexp( - a/к) + А2ехр( - а/к), (20)

Ук=ехр( - a3 Тк)(А3 Тк + A4), (21)

которым адекватна при k > 2 и аддитивной помехе авторегрессия второго порядка

YK= Л1 YK-1 - Л2 YK-2 +

где Л1 = ехр( - ajÄ) + ехр( - a2Ä) для модели (20) и Л1 = 2ехр( - a3Ä) для модели (21); Л2 = ехр( - (aj + a2)Ä) для модели (20) и Л2 = ехр( - 2a3Ä)/2 для модели (21).

Условием отнесения анализируемого временного ряда к модели суммы двух экспонент будет следующая система неравенств: 0 < Л1 < 2, 0 < Л2 < 0,25—^. Условиями принятия модели произведения экспоненты на линейную форму аргумента являются соотношения: 0< Л1 < 2, Л2 = 0,25—^. При N > 4 можно получить помехозащищенные оценки параметров

N

Л°ь —°2= arg min £(Yk - Л1 Yk-1 + Л2 Yk-2}2,

—1,—2 k=2

которые определятся из СЛАУ второго порядка и через которые можно рассчитать

1

1 ° ТО 2

Л 1 Л 1

1/2

a°i,2=--Ln ( — + ( — ± Л°2) ),

Ä

2

4

a3° = -—Ln —— 3 Ä 2

Отметим возможность идентификации параметров на доле существенного изменения ординат модели, за которую можно при-н , т ь ^3/ш1п(а1,а2}. Видно, что помехозащищенные оценки параметров определят усло -вия

N

{А°1, А°2} = arg min £{ Yк-A1еxp( -a°1ÄK) A1,A2 k=0

- A2еxp(-a°2Äк)}2 для (20),

N

{А°з, A°4} = arg min £{ Yk - (—°j/2)k x

A3,A4 K=0

x(A3äK+A4)}2

для (21)

и соответствующие СЛАУ второго порядка.

5.5. Практический интерес представляет и временной ряд, состоящий из произведения экспоненциального тренда и гармоники

Yk = А1ехр( - a1 Tk)Cos(&Tk + ф).

(22)

Для (22) при к > 2 и аддитивной помехе будем иметь следующую модель авторегрессии:

YK= Л, YK-1 - — YK-2 +

(23)

где Л1 = —0ехр( - a1Ä), Л2 = ехр( - 2a1Ä), Л0 =2CosraÄ.

Метод наименьших квадратов (при N > 4) дает оценки

N

Л°1, Л°2 = arg min£{ Yk - —1 Yk-1 + —2 Yk-2}2, — 1,— 2 к=2

N

a°2, a°3 = arg min £{ Yk -a2, a3 k=0

-a2 ехр( - a°1 Äk)Cos(°Äk +

+ a3 ехр( - a°1 ÄK)Sin( °Äk}2,

где A2 = A1Cosф, A3 = А^^ф.

Идентифицируемые параметры будут YkK = 02{ Yk-1 (К-l) + Yk-3(K-3)} -равны - 02 Yk-2(K-2) - Yk-4(K-4) + ¿к,

a°i= -

1 1

-LnА°2, ш°= — ArcCos(-

N

2Д

А

A°i = (А°22 + А°з2)1/2, ф° = Arctg(A°3/A°2).

5.6. Сумма экспоненциального тренда и гармонического компонента

Yk = A1 ехр( - aj Тк) + A2Cos(ra Тк + ф) (24)

определит следующую модель авторегрессии при k > 2 при аддитивной помехе:

Yk= n Yk-1 - n2 Yk-2 + n3 Yk-3 + ¿к, (25)

где п1 = ехр( - a2A) + 2CosraA, П2 = 1 + ехр( -ajA) + 2CosraA, n3 = ехр( -a2A). Предложенный подход дает оценки

N

П01,Л°2,Л03= arg minZ{ YK-^ YK-1 +

П1,П2,Пз к=0

2(А°2)

), 0°1, 0°2 = arg min Е{ YkK -01{ Yk-1 (К-1)

!/2

+ n2 yk-2 - n3 yk-3}2

(26)

N

a°ia°3a°4 = arg min 2{ Yk-alехр(a°lДк) a1a3a4 k=0

Aзехр(-a°2Дк)Cosю°Дк + A4COSq°Ak}2

(27)

где A3 = A2Cosф, A4= A2Si^, a°l = l/A Ln n

о

1/ A Arcсos{(n02 - l)/2n°},

2 (AO32 + A042)l/2, фо = Arctg(A°JAO3). (28)

ш

А°2 = (А°32 + А° 2)1/2

Соотношения (25) - (28) для модели (24) обладают внешним сходством с формулами для модели (18), но диапазон значений коэф -фициентов уравнений регресии (19) и (25) различен, что может служить признаком классификации моделей.

5.7. Для произведения гиперболического тренда и гармоники

Ук= (А1 + А2/ Тк)Со<ш Тк+ ф) (29)

целесообразен прием представления его в

виде УкТк=(АТк+А2)Сда(ш7К+ ф), а при к > 4 и аддитивной помехе - использования модели авторегрессии и формул

+

0l,02 К=4

+ Yk-3(K-3)} + 02 yk-2(K-2) + yk-4(K-4)}2

A°3, A°

N

A°5, A°6 = arg min Z{ Yk AK ■

A3A4AA k=0

- A3Aк Cosra°Aк + A4Aк Sinra°Aк -

- A5Cosa°AK+ A5 Sinra°AK}2,

где 0=4CosoA, 02 = 16 Cos2raA+2, A3 = A2Co^, A4 = A2Si^, A5 = A2Cosф, A6 = A2Si^, й° = l/A ArcCos{(0°2/4), A°6= (A°32 + A°42)l/2, A°2= (A°52 + A°62)l/2, ф°= Arctg(A'°4/ A3).

5. 8. Сумму гиперболического тренда и гармоники

Yk = Al + A2/ Тк+A3Cos(®TK+ ф)

(30)

целесообразно представить в виде Yk Тк = = A1Tk+A2 + A3TkCos(&Tk+ ф), для которо -го при k > 6 и аддитивной помехе получим

YkK = pl{ Yk-1 (К-l) + Yk-5(K-5)} -

- р2{ YK-2(K-2) + YK-4(K-4)}+

+ p3 YK-3(K-3)} - YK-6(K-6)) + ¿К,

N

p°l,p°2,p°3 = arg minE{ YkK - pl{ Yk-1 (K-l)+

pi,p2,p3 K=6

+ Yk-5(K-5)} + p2{ Yk-2(K-2) + Yk-4(K-4)}-

- p3 Yk-3(K-3)} + Yk-6(K-6)}2,

N

A°lA02A°lX,A0i = arg min!{ YkAK - A1AK -

AiA2AfA5 k=0

- A4AK Cos&°AK + A 5AK Sinra°AK}2,

где = l/A ArcCos{(p°l/4 - l/2),

pl =2(2CosoA + l),

p2 = (2CosoA + l)(2CosoA + 3),

p2 = 4(2Cos2oA + CosraA + l),

A4 = A3Cos^, A5 = A3Siиф, A°3 = (A°42 + A°52)l/2,

ф°= Arctg (A°5/A4).

о

3

2

Особенностью двух последних моделей авторегрессий является их «нестационарность», обусловленная присутствием в них произведений отсчетов на номера этих отсчетов. Для (29) N > 6, а для (30) N > 8.

5.9. Логистическая динамика наиболее часто моделируется выражением (моделью Верхулста или Перла-Рида)

Л

Ук =-, (31)

1 + В ехр( - СТк)

где А,В,СеЯ( и, как правило, положительны). Перейдем в (31) к значениям Хк =1/Ук:

2к = Е ехр( - СТк) + О + %к,

для которых, очевидно, будут справедливы модель авторегрессии вида (6) и соответствующие соотношения для расчета: Е°=В/А, О° = 1/А, а через них, соответственно, и

А°, В°.

В отличие от известных методов идентификации данный метод не требует априорного знания параметра «насыщения» Л.

6. Заключение. На примерах 13-ти рассмотренных моделей, широко применяемых в эконометрической практике, показана возможность получения оптимальных оценок параметров на основе авторегрессии отсчетов путем решения соответствующих СЛАУ невысокого (до четвертого) порядка.

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

Проведена и находится в режиме опытной эксплуатации программная реализация

предложенного подхода идентификации моделей (2), (4), (14), (16), (22), (24) в составе автоматизированных информационных систем управления администрации: «АИС - Го -род» г. Новокуйбышевска и «АИС - Регион» Самарской области, подтвердившая справедливость и эффективность предложенного подхода.

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

Список литературы

1. Гранберг А. Г. Динамические модели народного хозяйства. - М.: Экономика. 1985.

2. Статистическое моделирование и прогнозирование./ Под ред. А. Г. Гранберга. - М.: Финансы и статистика. 1990.

3. Айвазян С. А. Прикладная статистика. Основы эконометрики. - М.: 2001. ЮНИ-ТИ - ДАНА. 2001.

4. Кобринский Н. Е., Кузьмин В. И. Точ -ность экономико-математических моделей. -М.: Финансы и статистика. 1981.

5. Эконометрика/ Под ред. И. И. Елисеевой. - М.: Финансы м статистика. 2002.

6. Колемаев В. А. Математическая экономика. - М.: ЮНИТИ-ДАНА. 2002.

7. Деч Г. Руководство к практическому применению преобразования Лапласа и Х-преобразования. - М.: Наука.1971.

8. Семёнычев В. К. Эконометрическое моделирование расширенного воспроизводства на основе авторегрессии. Сб. «Вестник Самарского аэрокосмического университета. №4. Самара. 2003. - С. 26-33.

9. Грицан В. Н. Эконометрика. - М.: Дашков и К. 2001.

GENERAL APPROACH TO THE ECONOMIC DYNAMICS IDENTIFICATION ON THE BASISOF AUTOREGRESSION MODELS

© 2004 V. K. Semenychev

Samara State Aerospace University

General approach to the identification of economic dynamics models containing exponential and harmonic components is proposed.

The approach is based on autoregression of economic index readings.