Федеральное государственное автономное образовательное учрежение высшего образования



Pdf көрінісі
бет25/26
Дата22.01.2017
өлшемі5,2 Mb.
#2428
1   ...   18   19   20   21   22   23   24   25   26
ГЛАВА

 

ПОСТРОЕНИЕ ЭКОНОМЕТРИЧЕСКИХ 

МОДЕЛЕЙ И ИХ КОРРЕКЦИЯ 

 

Мотина Валентина Григорьевна

старший преподаватель 

 

6.1.



 

Эконометрические модели 

 

Применяемые  в  эконометрике  регрессионные  методы  и 



методы  анализа  временных  рядов  в  большинстве  своем 

основываются 

на 

разложении 



наблюдаемых 

значений 

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

составляющие.  Данные  в  эконометрике  могут  относиться  к 

различным  типам:  перекрестная  выборка  –  это  форма 

организации  данных,  в  которых  наблюдения  одной  и  той  же 

характеристики  фиксируются  одномоментно  или  собираются  за 

один временной промежуток для разных объектов. Во временных 

рядах  проводится  фиксация  характеристик  одного  объекта  в 

последовательные  моменты  времени.  В  панельных  данных 

фиксация проводится по обоим измерениям – в разные моменты 

времени для одного и того же набора объектов. 

Хотя данные разных типов отличаются по своим свойствам, 

применяется много общих методов. 

Так,  в  методе  декомпозиции  структура  временного  ряда 

разлагается  на  составные  части,  которые  суммируются 

(аддитивная 

форма 


модели) 

или 


перемножаются 

(мультипликативная  форма).  Применяется  также  смешанная 

мультипликативно-аддитивная  форма.  Обычно  в  методе 

декомпозиции  выделяют  4  структурно  образующих  элемента 

временного ряда: 

t

T

  –  тренд  –  устойчивое  изменение  среднего  уровня 

показателя в течение длительного времени; 

t

S

  –  сезонная  компонента  –  устойчивые  периодические 

колебания  уровней  ряда  (суточные,  недельные,  месячные, 

квартальные и т.д.); 



 

406 


t

C

  –  циклическая  компонента  –  довольно  длительные, 

обычно несколько лет, периоды подъема и спада; 

t

  –  остаточная  компонента,  образующаяся  под  влиянием 



различных неучтенных причин. 

В  аддитивной  форме  модель  декомпозиции  может 

выглядеть так: 

 

t



t

t

t

t

S

C

T

y





(1) 


Единственным  обязательным  элементом  временного  ряда 

является остаточная компонента. 

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

"белым  шумом",  то  есть  представляет  собой  стационарный 

случайный  процесс  второго  порядка  (постоянное  среднее: 

)

(



i

M

= 0 в  любом наблюдении,  и  постоянная дисперсия 



)

(

i



D

=



)

(

j



D

=



2

  для  любых  наблюдений  с  номерами 



i

  и  )  при 

отсутствии  автокорреляции  и  отсутствии  корреляции  между 

необъясненными  остатками  и  исходным  временным  рядом.  Эти 

условия 

называют 

условия 

Гаусса-Маркова. 

Остатки 

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

некоторому закону, обычно  нормальному, и предполагается, что 

параметры закона распределения не изменяются. 

В более общем виде аддитивную модель временного ряда с 

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

 

t

t

t

X

t

f

y



)

,



(

(2) 



где: 

t

y

 – наблюдаемые значения временного ряда; 

)

,

(



t

X

t

f

 – детерминированная составляющая, 



t

X

  –  значения  детерминированных  переменных,  факторов, 

оказывающих  влияние  на  величину  f  в  момент  или  период 

времени t

В формуле (2) переменные 

t

X

  могут  отсутствовать.  Кривая 

с  единственным  аргументом  "время"  называется  кривая  роста 

(Рис.  5).  Во  временных  рядах  реальных  экономических 

показателей, 

представленных, 

например, 

на 


сайте 

Территориального органа Федеральной службы государственной 

статистики  по  Республике  Крым  [1],  большинство  зависимостей 

нелинейные. 



 

407 


 

Рис. 5. Наблюдаемые значения временного ряда и модель 

детерминированной составляющей. 

 

Для перекрестных данных модель имеет вид: 



 

i

i

i

X

f

y



)

(



(3) 


где 

i

 – номер наблюдения. 

Для  перекрестных  данных  наиболее  распространенная 

форма  модели  –  линейная,  в  особенности,  если  модель 

многофакторная. Применяются и нелинейные модели, например, 

производственная  функция  Кобба-Дугласа,  но  для  оценки 

параметров  нелинейные  модели  обычно  линеаризуются  заменой 

переменных или логарифмированием. 

Для  моделирования  сложных  систем  используют  системы 

одновременных или регрессионных уравнений. 

Наиболее  широко  для  оценки  параметров  моделей 

используется 

метод 

наименьших 



квадратов 

(МНК), 


минимизирующий  сумму  квадратов  остатков.  Основной 

недостаток МНК – придание большого веса наблюдениям, далеко 

отклоняющимся  от  основной  тенденции.  При  моделировании 

перекрестных  данных  резко  выделяющиеся  наблюдения  как  по 

оси  Ox,  так  и  по  оси  Oy  (выбросы)  часто  удаляют  из  выборки, 

иногда выделяют в отдельную группу. По возможности в набор X 

добавляют  фактор,  который  объясняет  причину  выбросов.  Этот 

фактор  может  быть  и  фиктивной  переменной,  принимающей, 

например,  значения  1  в  праздничные  дни  и  0  в  остальные  дни. 

Для  временных  рядов  МНК  придает  большой  вес  начальным  и 

последним  наблюдениям.  Для  упрощения  вида  функции  f 

начальные наблюдения можно удалить (рис.Рис. 5). 



t

y

Модель


Наблюдение

y

t

t

 

408 


В  результате  получаем  оценки  параметров  функции  f

Значения,  полученные  подстановкой  в  модель  значений  t  и  X

обозначают 

t

y

  или 



i

y

,  а  отклонения  модели  от  наблюдаемых 



значений – соответственно 

t

e

 или 


i

e

Для  моделей,  содержащих  константу  (свободный  член) 



выполняется  правило  сложения  дисперсий:  общая  сумма 

квадратов  равна  сумме  объясненной  суммы  квадратов  и  суммы 

квадратов остатков: 

 











n

t

t

t

n

t

t

n

t

t

y

y

y

y

y

y

1

2



1

2

1



2

)

(



)

(

)



(



(4) 


где: 

n

 – число наблюдений. 

Поскольку нет уверенности в правильности выбора функции 

f  (это  касается  и  выбора  факторов,  входящих  в  набор  X), 

полученный  результат 



y

  можно  считать  достаточно  "хорошей" 



моделью неизвестной функции 

f

, если модель обладает высокой 

точностью.  

Основным 

показателем 

точности 

модели 

является 



остаточная дисперсия: 

 

k



n

e

Se

n

t

t



1

2



2

(5) 



где: 

k

 – число параметров модели, 



t

t

t

y

y

e



Величину 



2

Se

Se 

 

называют 



среднеквадратической 

(стандартной) ошибкой модели. 

Рассчитывают  также  среднюю  относительную  ошибку 

модели 


%

100




y



Se

V

. Для реальных экономических показателей 

рекомендуется, чтобы значение 

V

 не превышало 20%. 

Для 

модели 


множественной 

линейной 

регрессии 









)

(

)



2

(

2



)

1

(



1

...


m

m

x

x

x

y

 среднеквадратические ошибки 

значений 

j

a

 – оценок коэффициентов 



j



m

j

...


1

, – определяют 



по формуле: 

 

jj



j

Z

Se

a

s



2

2

)



(

(6) 



где: 

2

Se

 – см. формулу (5); 

jj

Z

 – диагональный элемент 

матрицы 

 


1



X



X

T

Для проверки  гипотез о статистической значимости оценок 



и  определения  доверительных  интервалов  для  истинных 

 

409 


значений  коэффициентов  регрессии  используется  утверждение: 

величина 

 

)

(



)

(

j



j

j

j

a

s

a

a

t





m



j

...


1



(7) 

подчиняется распределению Стьюдента с числом степеней 

свободы 

1



 m

n

Кроме  того,  проверяется  выполнение  условий  Гаусса-



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

t

e

  белым 


шумом,  т.е.  оценкой 

t

.  Если  да,  то  функцию 



f

  можно  считать 

моделью детерминированной составляющей. 

Априорно  проверить  выполнение  условий  Гаусса-Маркова, 

как правило, невозможно. Проверку выполняют, получив оценку 

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

то  есть  фактическую  ошибку  модели.  Иногда  случайная 

компонента  может  быть  больше,  иногда  меньше  по  абсолютной 

величине,  однако  не  должно  быть  причины,  порождающей 

большую  ошибку  в  одних  наблюдениях  и  малую  в  других, 

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

исследуемого процесса описана моделью регрессии.  

Полученные  результаты  моделирования  используются  для 

анализа  характеристик  ряда,  например,  особых  точек, 

показателей  динамики  временного  ряда  и  т.п.,  а  также  для 

получения прогноза. Точечный прогноз получают подстановкой в 

выражение для 

y

 будущих значений 



t

 и ожидаемых значений 



p

X

,  а  характеристики  остатков  используют  для  получения 

интервального прогноза. 

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

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

на  треть  размаха  (по  каждому  фактору).  Формулы  средней 

квадратической ошибки прогноза для линейной модели: (8) – для 

множественной  регрессии,  (9)  –  для  парной  регрессии.  Аналог 

формулы  (9)  существует  для  полиномов  2  или  3  степени,  для 

других зависимостей формул нет. 

 

 




T



p

T

p

X

y

X

X

X

X

Se

s

p





1

1



(8) 



 





2



2

)

(



)

(

)



(

1

1



x

x

x

x

n

Se

s

i

p

x

y

p



(9) 

 

410 


Границы интервального прогноза с заданной доверительной 

вероятностью (обычно 0,95) получают по формуле: 

 

кр

)



(

)

(



t

s

p

p

x

y

x

y





(10) 

где: 


кр

t

 – критическая точка распределения Стьюдента. 

При большом объеме выборки вклад слагаемых под корнем 

незначительный,  и  тогда  формулу  (10)  можно  записать  в 

упрощенном виде: 

 

кр



)

(

t



Se

p

x

y



(11) 



 

6.2.

 

Методы анализа и устранения 

гетероскедастичности 

 

Одним  из  условий  Гаусса-Маркова,  которое  часто 



нарушается  при  моделировании  как  временных  рядов,  так  и 

перекрестных 

данных, 

является 

гомоскедастичность 

– 

постоянство  дисперсии  остатков.  Невыполнение  этого  условия 



называется гетероскедастичность. 

Гетероскедастичность  делает  невозможным  определение 

доверительного 

интервала 

для 

истинных 



значений 

коэффициентов  линейной  регрессии  по  формуле  (7).  Возможно, 

стандартные  ошибки  будут  занижены,  а  следовательно, 

t

-

статистика  –  завышена,  и  будет  получено  неправильное 



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

регрессии. 

При 

условии 


гомоскедастичности 

коэффициенты 

линейной  регрессии,  найденные  обычным  МНК  (методом 

наименьших  квадратов),  имеют  наиболее  низкую  дисперсию 

среди  всех  несмещенных  оценок,  являющихся  линейными 

функциями  от  наблюдений 



Y

,  т.е.  эффективны.  При 

гетероскедастичности  оценки  МНК  неэффективны.  Затруднено 

также  построение  доверительных  границ  для  прогноза:  в  основе 

экстраполяции лежит предположение, что дисперсия отклонений 

останется без изменений и на периоде упреждения.  

На  практике  гетероскедастичность  может  наблюдаться 

как  в  перекрестных  данных  (например,  зависимость 

потребления  некоторого  продукта  от  дохода  у  разных 

домохозяйств),  так  и  во  временных  рядах,  имеющих 

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


 

411 


Гетероскедастичность  можно  наблюдать  на  графическом 

представлении зависимости остатков или их квадратов от 



t

, от 


x

от 



y

.  Для  проверки  используются  различные  статистические 



критерии  (тесты).  Например,  ранговой  корреляции  Спирмена, 

Уайта, Парка, Глейзера, Голдфельда-Квандта. 

Для 

теста 


Уайта 

оценивается 

(обычным 

МНК) 


вспомогательная  регрессия  квадратов  остатков 

2

i



e

  на  все 

регрессоры 

)

j



x

, их квадраты и попарные произведения. В случае 

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

незначимой.  Для  проверки  этой  гипотезы  используется 

статистика 

2

R



,  где 


2

R

  –  коэффициент  детерминации 

вспомогательной 

регрессии. 

При 

отсутствии 



гетероскедастичности  данная  статистика  имеет  асимптотическое 

распределение 

)

1

(



2



N

,  где 


N

  –  количество  параметров 

вспомогательной  регрессии.  Можно  использовать  для  проверки 

также F-статистику. 

В  условиях  непостоянства  дисперсии  рекомендуется 

применять  другие  методы  оценки  коэффициентов  регрессии, 

например,  обобщенный  МНК  (ОМНК)  для  модели  с 

гетороскедастичностью, или взвешенный МНК. 

В ОМНК оценку параметров модели регрессии получают по 

формуле [2]: 

 





Y

X

X

X

A

T

T





1

1



(12) 

где: 


 – ковариационная матрица остатков. 

 

При полном отсутствии автокорреляции остатков матрица 



 

– диагональная. Элементы главной диагонали – значения 



2

i



i

 – 


номер  наблюдения.  В  таком  случае  ОМНК  называется 

взвешенным МНК. 

Каждому  значению 

i

e

  присваивается  "вес",  равный 



i

1



Идея  метода  заключается  в  компенсации  недостатка  МНК, 

заключающегося в придании большого веса выбросам. Придание 

меньшего  веса  более  рассеянным  данным  обеспечивает 

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

приводит  к  получению  наиболее  эффективных  оценок 

параметров модели. 


 

412 


Поскольку  значения 

i

  обычно  неизвестны,  в  качестве 



2

i

 



используют оценки 

2

i



e

, полученные, например, по квадратичной 



модели  в  ходе  применения  теста  Уайта.  Возможны  и  другие 

формы  зависимости 

2

i

e

  от  регрессоров 

)

j



x

например,  линейная 

или экспоненциальная. 

Принимая 

2

i

i

e

 



, получают новые переменные 



i

i

i

y

y



*



i

j

i

j

i

x

x



)

(

*



)

(



.  Оценки  параметров  модели  регрессии  для 

новых переменных и есть оценки взвешенного МНК. 

Известны  и  другие  рекомендации  [3]:  сделать 

предположение 

относительно 

характера 

гетероскедастичности,  например:  "

2

i

 

пропорциональны 



значениям 

x

" или "


2

i

 пропорциональны значениям 



2

x

". 


При  первом  предположении  уравнение  парной  линейной 

регрессии 









x

y

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

части на 

i

x

к виду (13) (уравнение не содержит свободного члена). 

 

v

x

x

x

y





(13) 


Для 

отклонений 



v

 

выполняется 



условие 

гомоскедастичности.  

Для  многофакторной  модели  можно  делить  не  на 

x

,  а  на 



y



При  втором  предположении  выполняют  деление  на 

x

  и 


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

v

 тоже гомоскедастичны: 

 

v

x

x

y





(14) 


Пример 

построения 

модели 

в 

условиях 

гетероскедастичности.  Применим  приведенные  выше  методы 

для  построения  парной  линейной  регрессии  по  данным 

статистики Крыма [1] за 2009-2013 гг. 

Результативный  признак 



y

 

–  "Объем  розничного 



товарооборота  предприятий,  осуществляющих  деятельность  в 

сфере  розничной  торговли  и  по  обеспечению  блюдами  и 

напитками, млн. грн."; фактор 

x

 – "Среднемесячная номинальная 

заработная  плата  в  расчете  на  одного  штатного  работника  (в 


 

413 


среднем  с  начала  года),  грн.".  Динамика  показателей 

представлена на рисунке 2. 

 

 

Рис. 6. Динамика показателей x и y. 



 

Данные за январь в статистике не представлены, кроме того, 

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

каждом году 8 значений, 



n

=40. 


Летом  за  счет  гостей  Крыма  показатель 

y

  существенно 

увеличивается:  в  июне  примерно  на  треть,  в  июле  в  1,9  раза,  в 

августе  –  более  чем  в  2  раза.  Это  можно  использовать  для 

косвенной оценки "удачности" курортного сезона [2]. 

Зависимость 



y

  от 


x

линейная  (Рис.7).  Получена  модель 

регрессии: 

 

15



,

610


743

,

0





x



y



(15) 

Значения  параметров:  коэффициент 



a

=0,743,  свободный 

член 

b

=  -610,15.  Оценка  коэффициента  является  статистически 

значимой: 

)

(a



t

=  9,98;  критическое  значение 



t

  на  уровне 

значимости 0,05 и числе степеней свободы 38 равно 2,02. 

 

 



Рис.7. Корреляционное поле и модель регрессии 

0

500



1 000

1 500


2 000

2 500


3 000

0

8



16

24

32



40

x, грн.


y, млн.грн.

y = 0,743x - 610,15

0

500


1 000

1 500


2 000

1400


1600

1800


2000

2200


2400

2600


2800

y, млн.грн.

Линейный (y, млн.грн.)


 

414 


Однако  хорошей  построенную  модель  назвать  нельзя: 

средняя  относительная  ошибка 



V

  =  20%  –  на  пределе 

допустимого, и остатки гетероскедастичны. На рисунке 4 видим, 

что  большим  значениям 



x

  соответствует  большее  рассеивание 

остатков.  

 

Рис. 8. Остатки (отклонения модели от фактических 



значений). 

 

Применив  формулы  (9)-(10)  к  исходным  данным, 



рассчитаем  границы  интервального  прогноза  с  доверительной 

вероятностью 0,95. Как видим на рисунке 5, для малых значений 



x

  ширина  интервала  слишком  завышена,  а  экстраполяция  в 

сторону увеличения 

x

 приведет к тому, что интервал не накроет 

истинных значений прогнозируемого показателя. 

 

Рис. 9. Доверительные границы интервального прогноза. 



 

Гетероскедастичность  остатков  подтверждает  тест  Уайта 

(квадраты  остатков  представлены  на  рисунке  6).  Значение 

2

R



=10,7; критическое значение 

2



 равно 6,0. 



Для применения взвешенного МНК используем полученную 

модель  квадратов  остатков: 

26167

587


,

12

0189



,

0

2



2





x

x

e

.  Как 



-500

-250


0

250


500

1400


1600

1800


2000

2200


2400

2600


2800

3000


0

500


1 000

1 500


2 000

1400


1600

1800


2000

2200


2400

2600


2800

 

415 


видим на рисунке 4, для малых значений 

x

 выполнить расчет по 

формуле 

2

i



i

e

 



не  удается,  так  как  модельные  значения 

отрицательны.  Из  рекомендуемых  моделей  от  этого  недостатка 

застрахованы,  пожалуй, только степенная  или экспоненциальная 

форма  зависимости.  Кстати,  тест  Парка  основан  на  степенной 

форме: 


v

a

e

x

b



2



Пропустив  четыре  значения,  для  остальных  наблюдений 

выполним  преобразование 

i

i

i

y

y



*



i

i

i

x

x



*

.  Связь  новых 



переменных  представлена  на  рисунке  7.  Видим,  что 

зависимости 

*

y

  от 


*

x

  нет,  наклон  прямой  получился  только 

из-за  неоднородности  данных.  Оценки  параметров  модели 

регрессии 



a

=0,267, 


b

 =2,303 использовать невозможно: они не 

подходят для моделирования исходного ряда. 

 

Рис. 10. Вспомогательная регрессия квадратов остатков на 



факторы и их квадраты. 

 

Рис. 11. Регрессия для новых переменных (взвешенный 



МНК). 

y = 0,0189x

2

 - 12,587x - 26167



R

2

 = 0,2667



-50 000

0

50 000



100 000

150 000


200 000

250 000


1400

1600


1800

2000


2200

2400


2600

2800


e2

Полиномиальный (e2)

3

4

5



6

7

8



9

10

8



10

12

14



16

18

20



22

24

26



 

416 


Предположения  о  пропорциональности 

2

i

  значениям 



x

  и 


значениям 

2

x

  привели  к  результату:  по  формуле  (13)  получены 

оценки: 


a

=0,733, 


b

 = -587,2, по формуле (14): 



a

=0,722, 


b

 = -564,9. 

Как  видим,  отличие  коэффициентов  от  модели,  полученной  с 

помощью обычного МНК, очень мало. 

 

 

Рис.12. Результаты моделирования по формулам (13) и (14). 



 

Результаты  расчета 



y

  с  применением  найденных  оценок 



представлены  на  рисунке  8.  Там  же  для  сравнения  –  линия, 

представленная на рисунке 3. 

Кроме 

того, 


устранения 

гетероскедастичности 

не 

произошло,  хотя  значения 



2

R

  несколько  уменьшились:  для 

модели,  полученной  по  формуле  (13),  получилось  9,2,  по 

формуле  (14)  –  7,5.  Но  оба  значения  превосходят  критическое 

значение 

2



, равное 6,0. На рисунке 9 представлена зависимость 

остатков 



v

,  полученных  по  формуле  (14),  от 



x

1

.  Как  видим, 



дисперсия не постоянна. 

 

 



Рис. 13. Зависимость остатков функции. 

y = 0,743x - 610,15

0

500


1 000

1 500


2 000

1400


1600

1800


2000

2200


2400

2600


2800

3000


y, млн.грн.

y1^


y2^

Линейный (y, млн.грн.)

-0,2

-0,1


0

0,1


0,2

0,0003


0,00035

0,0004


0,00045

0,0005


0,00055

0,0006


0,00065

0,0007


v

 

417 


 

Можно 


сделать 

вывод, 


что 

методы 


устранения 

гетероскедастичности,  доказанные  в  теории,  на  практике 

ожидаемого результата не приносят. 

Даже  если  удастся  получить  гомоскедастичные  остатки  с 

применением  ОМНК  или  взвешенного  МНК,  результатом  будут 

только  эффективные  оценки  параметров  модели  регрессии. 

Поскольку  в  этих  методах  минимизируется  сумма  квадратов 

остатков  не  исходного  ряда,  а  преобразованного,  применить 

результаты 

для 


получения 

интервального 

прогноза 

затруднительно.  По  крайней  мере,  описание  этих  методов  в 

литературе до интервальных оценок не доходит. 

Хорошую 


рекомендацию 

по 


устранению 

гетероскедастичности можно найти, например, в [3]: «Во многих 

случаях  дисперсия  отклонений  зависит  не  от  включенных  в 

уравнение  регрессии  объясняющих  переменных,  а  от  тех, 

которые не включены в модель, но играют существенную роль в 

исследуемой  зависимости.  В  этом  случае  они  должны  быть 

включены  в  модель".  И  далее:  "В  ряде  случаев  для  устранения 

гетероскедастичности  необходимо  изменить  спецификацию 

модели 

(например, 



линейную 

на 


лог-линейную, 

мультипликативную на аддитивную и т.п.)». 

Рекомендация  по  поводу  включения  дополнительных 

факторов  хорошая,  но  не  всегда  практически  реализуемая:  сбор 

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

вторичных  данных,  представленных  в  официальной  статистике, 

может просто не быть. 

Преобразование  данных  (например,  логарифмирование) 

приводит к изменению характеристик как исходных показателей, 

так и остатков, хотя, конечно, находит применение. 

Для  временных  рядов  переход  от  исходного  ряда  к  ряду 

цепных приростов или коэффициентов (темпов) роста во многих 

случаях  дает  хороший  результат  (проверено  нами  при 

моделировании реальных статистических показателей). 

На  последней  рекомендации  остановимся  более  подробно, 

хотя  выражение  перефразируем:  замена  аддитивной  модели  на 

мультипликативную. 

 


 

418 



Достарыңызбен бөлісу:
1   ...   18   19   20   21   22   23   24   25   26




©emirsaba.org 2024
әкімшілігінің қараңыз

    Басты бет