CC BY
35
УДК 630*3
А.М. Меньшиков, А.М. Копейкин
Меньшиков Александр Михайлович родился в 1958 г., окончил в 1985 г. Архангельский лесотехнический институт, заведующий лабораторией транспорта леса ОАО «Научдревпром - ЦНИИМОД». Имеет более 10 печатных работ в области лесовозного транспорта.
L II 1
Копейкин Адольф Михайлович родился в 1936 г., окончил в 1959 г. Архангельский лесотехнический институт, доктор технических наук, профессор кафедры лесопильно-строгальных производств Архангельского государственного технического университета, заслуженный работник лесной промышленности РФ. Имеет более 100 печатных трудов в области прогнозирования, технологии лесопиления и деревообработки.
ПРИМЕНЕНИЕ СПЕКТРАЛЬНЫХ МЕТОДОВ В ИССЛЕДОВАНИЯХ ТЕХНОЛОГИЧЕСКИХ ПРОЦЕССОВ ЛЕСОЗАГОТОВИТЕЛЬНОГО ПРОИЗВОДСТВА
С позиций теории случайных функций рассмотрена внутренняя структура и динамика тренд-сезонных колебаний нестационарных технологических процессов вывозки древесины и производства круглых лесоматериалов предприятиями лесопромышленного комплекса. Предложена методика количественного определения силы, продолжительности действия и направленности взаимных связей параллельных технологических процессов лесозаготовок на основе кросс-спектрального анализа их динамических рядов.
Ключевые слова: лесозаготовительное производство, динамические временные ряды, спектральный анализ процессов.
Лесозаготовительное производство в целом и отдельные его фазы подвержены влиянию общей тенденции развития лесопромышленного комплекса, значительным сезонным колебаниям и воздействию случайных флуктуаций. Годовые ритмы идущих параллельно процессов лесосечных работ, вывозки древесины, производства круглых лесоматериалов, поставки автотранспортом сортиментов во двор потребителя взаимно обусловливают друг друга, а сезонные колебания объемов производства в них происходят почти синхронно. Пример годичных колебаний объемов производства на некоторых фазах лесозаготовок представлен на рис. 1.
Как видно, наиболее вероятные графики объемов производства у всех процессов идентичны. Они имеют вид '-образной осциллирующей кривой с двумя примерно равными по значению минимумами в периоды весенней и осенней распутицы. Периодичность колебаний явно выражена. Очевидно, что при таком характере протекания процессов их динамические
Рис. 1. Динамика объемов производства по фазам лесозаготовок в 2000 г.: 1 - вывозка древесины; 2 - производство круглых лесоматериалов; 3 - поставка сортиментов автотранспортом в Архангельский промузел
ряды не являются в статистическом отношении стационарными, содержат сильную автокорреляцию и потенциальную мультиколлинеарность уровней. Это исключает возможность прямого применения к ним распространенных методов многомерного статистического анализа, существенно затрудняет изолированный, тем более комплексный технологический анализ процессов.
Методологический подход к исследованию и анализу такого рода процессов в общем случае включает следующие этапы. На первой стадии формируют динамические временные ряды показателей. Далее предпринимают индуктивное разложение исходных временных рядов на детерминированные и стохастические компоненты с последующей аппроксимацией детерминированных составляющих функциональными моделями. Затем производят либо фильтрацию исходных рядов от случайных компонентов, либо элиминирование из исходных рядов систематических составляющих. При элиминировании ряды остатков должны соответствовать основным статистическим гипотезам, на которых базируются классические критерии согласия и методы многомерного статистического анализа, т. е. в максимальной степени отвечать свойствам белого шума и иметь наилучшее качество аппроксимирующих функций. Очевидно, что в результате изменяется лишь масштаб колебаний исходных данных, а остаточные ряды сохраняют всю информацию о вариации процесса и его внутренней структуре. В дальнейшем, в зависимости от цели исследования, используют сами функциональные модели, или остаточные ряды процессов [1-5].
Основываясь на данных методологических положениях, авторы ранее получили остаточный ряд обобщенного транспортно-технологического процесса вывозки древесины предприятиями лесопромышленного комплекса Архангельской области в 1998-2002 гг., обладающий свойствами белого шума [6]. На примере параллельных технологических процессов вывозки древесины и производства круглых лесоматериалов для одной и той же группы предприятий предпринята попытка количественно определить силу и направленность их взаимных связей. В основании анализируемых процессов лежат массивы месячных статистических данных из 3180 единиц каж-
Рис. 2. Развертки динамических временных рядов процессов: А - исходные ряды; Б - ряды остатков; 1 - вывозка древесины; 2 - производство круглых
лесоматериалов
дый, собранных и агрегированных в динамические ряды из 60 уровней по разработанной нами методике [6]. Отсюда же взят остаточный ряд процесса вывозки древесины с соответствующими статистическими характеристиками.
Развертки временных рядов представлены на рис. 2.
Приняты следующие обозначения: отрезок времени, равный одному месяцу, обозначен индексом 7 (7 = 1, 2, ..., I); одному году - индексом 7 (/ = 1, 2, ..., ^, текущие значения рассматриваемого периода времени -индексом 7 (7 = 1, 2,..., Т), причем Т = и. Уровни исходного ряда вывозки обозначены {Х(}, производства круглых лесоматериалов - Тогда, по аналогии с процессом вывозки, каждому показателю процесса, характеризующему объем производства круглых лесоматериалов в 7-м месяце 7-го года, ставится в соответствие уровень динамического ряда:
zt = уг + st +
(1)
где у, - средние значения соответственно тренда и сезонной компоненты на отрезке времени 7; ег - остаточные колебания.
Первые два слагаемых в правой части (1) представляют собой детерминированную систематическую составляющую процесса, которую необходимо определить какой-либо аналитической функцией времени с наилучшими аппроксимационными свойствами. Третье слагаемое в уравнении (1) является стохастической переменной, включающей остатки от аппроксимации функций тренда и сезонной составляющей, а также всевозможные случайные ошибки обмера, учета и первичной обработки исходных данных.
На рис. 3 показана динамика объемов вывозки древесины и производства круглых лесоматериалов в 1998-2002 гг.
Анализ динамики объемов производства ЛПК Архангельской области в 1998-2002 гг. показал, что для аппроксимации трендов в этом периоде времени наиболее подходящими являются квадратичные функции, а для се-
1998 1999 2000 2001 2001 1998 1999 2000 2001 2002
а 6
Рис. 3. Динамика годовых объемов вывозки древесины (а) и производства круглых лесоматериалов (б) в 1998 - 2002 гг.
зонных компонент - полигармонические функции Фурье [6]. Вследствие этого тренд ряда производства круглых лесоматериалов аппроксимировали параболой
y(j) = 7,840 +1,38 j - 0,138 j 2, (2)
коэффициент детерминации составил R2 = 0,977. Как видно на рис. 3, градиенты трендовых функций процессов почти аналогичны.
Произведя интерполяцию тренда внутри интервалов анализируемого периода, перейдем в формуле (2) от годичной переменной j = 1, 2,..., 12 к текущей переменной t. Выражение для математического ожидания среднего значения тренда имеет вид
1 Т
Е[y(t)]=- J(7,840 + 0,116t - 9,583-10-4 t2) dt. (3)
T 0
После интегрирования и преобразований формула (3) принимает окончательный вид
Е[ y(t)] = 7,84 + 5,79 -10"2 T - 3,19 -10"4T2. (4)
При подстановке Т = 60 в формулу (4) математическое ожидание среднего значения тренда ряда производства круглых лесоматериалов составило E[y(t)]=10,161 тыс. м3/мес, что на 1,267 тыс. м3/мес меньше, чем значение тренда ряда вывозки древесины, равное 11,428 тыс. м3/мес [6]. Таким образом, среднестатистический выход круглых лесоматериалов в рассматриваемом периоде составил 88,9 % от объема вывезенной древесины.
Для аппроксимации сезонной компоненты использовали ряд Фурье an N
s(t) = -0 + 2 К cos n(th + j -1) + bn sin n(th + j -1)] (5)
2 n=1
с коэффициентами
1 n 1 n 1 n
a0 = — Js(t) dt; a = _ J s(t)cos ntdt ; bn = — J s(t)sinntdt, (6)
n -n n -n n -n
где n - число гармоник ряда Фурье, n = 1, 2, ..., N.
Коэффициенты an, bn определяли для опорного годичного интервала
с использованием значений исходного ряда z(th), h = 1, 2, ..., 12, представ-
ляющих собой средние арифметические объемов производства, вычисленные для пяти одноименных месяцев каждого года рассматриваемого периода. Расчет ап, Ьп выполняли по итеративному алгоритму с использованием последовательно двух, трех и четырех первых гармоник ряда Фурье.
Получены следующие значения коэффициентов полинома (5): а0/2 = = 10,161; а1=3,621; Ь = 2,362; а2 = 1,373; Ь =1,769; а3 = - 0,286; Ьъ = - 0,728; а4 = - 0,122; Ь4 = - 0,946.
Вычисленные значения оценок ) в соответствующих частотах опорного интервала представлены в таблице.
Частоты, Месяц опорного интервала th
оценки I II III IV V VI VII VIII IX X XI XII
ю(4) 0 Ж 6 Ж 3 Ж 2 2 — Ж 3 5 — Ж 6 Ж 7 — Ж 6 4 — Ж 3 3 — Ж 2 5 — Ж 3 11 -Ж 6
§(4) 15,072 15,536 16,353 12,079 7,455 7,833 8,396 8,352 8,066 5,900 7,560 13,203
Далее, центрируя амплитуды полинома (5) в соответствующих частотах по математическому ожиданию функции тренда (2), получим аналитическое уравнение для тренд-сезонных изменений объемов производства круглых лесоматериалов в рассматриваемом периоде:
(y + S) = 7,840 + 0,1167 - 9,583-10-4 t2 +
n=4 (7)
+ Z [ a„ cos n (a(th ) + j -1) + bn sin n (a(th ) + j -1)].
n=1
Значения частот для месяцев th опорного интервала приведены в таблице.
Уровни остаточного ряда {et} определяли с помощью формулы (1), вычитая из zt соответствующую детерминированную составляющую (y + st), определенную по уравнению (7). Проверка нормальности распределения уровней с помощью показателей асимметрии и эксцесса, а также по R-критерию Романовского показала, что полученный остаточный ряд реализует гауссовский процесс. Вид остаточного ряда ez показан на нижнем графике рис. 2.
Одновременно с вычислением значений остаточного ряда на каждой итерации по известным формулам [5] определяли первый циклический коэффициент автокорреляции ^(z) и производили диагностику уравнения (7) с помощью теста Дарбина - Уотсона на автокорреляцию уровней ряда ez(t), а также проверку с помощью кумулятивной периодограммы Qk на наличие в ряде существенных периодических составляющих. Граничные значения критериев принимали из работ [1, 4]. Расчеты показывают, что положительные результаты тестов получаются, как правило, при включении в представление Фурье в уравнении (7) четырех гармоник, т. е. при n = 4.
Рис. 4. Эмпирические коррелограммы одномерных остаточных рядов: 1 - вывозка древесины; 2 - производство круглых лесоматериалов
Эмпирические коррелограммы остаточных рядов процессов {Х^ и ^} представлены на рис. 4.
Малая величина первого коэффициента автокорреляции г1(г) = 0,093
и оценка <Л = 1,661, превышающая верхнее табличное значение критерия Дарбина - Уотсона <и = 1,608 при расчетных степенях свободы свидетельствуют о том, что в рассматриваемом ряде остатков ег(1) нет автокорреляции уровней. Однако обращают на себя внимание коэффициенты автокорреляции ряда производства круглых лесоматериалов г23(г) = 0,494 и г25(г) = - 0,327, почти совпадающие с соответствующими коэффициентами ряда остатков вывозки г23(х) = 0,519 и г25(х) = - 0,372. Совпадают при этом и коэффициенты взаимной корреляции процессов, вычисленные по следующей формуле [3]:
ге (х, z) = Се (х,г)/ .¡С^х^С^г . (8)
Здесь 0 - временной сдвиг уровней, 0 = - п, - п + 1, ..., п - 1, п;
п - размер используемого далее спектрального окна, принятого с учетом рекомендаций в [5], равный одному годичному циклу, т. е. п =12;
С0(х, г) - значения взаимных ковариаций при сдвиге 0:
1 т-е _ _
се(хг) =-!(х "х)(^+е"
Т t=l (9)
Се (х,г) =С _е (г, х), С0 (х, г) = С0 (г, х).
Графики взаимных корреляционных функций представлены на рис. 5. Из них видно, что коэффициенты взаимной корреляции принимают наибольшие абсолютные значения при временном сдвиге, равном 23 . 25 мес.: г23(х, г) = 0,584 и г23(г, х) = 0,378, г25(х, г)= - 0,230 и г25(г, х) = = - 0,448. Это свидетельствует либо о присутствии в ряде элемента периодичности, проявляющегося с лагом примерно в два года, либо о так называемом «эффекте Слуцкого», когда в сравнительно гладких коррелограммах имеют место обособленные всплески, объясняемые лишь спецификой алгоритма суммирования случайных величин [3]. Второе предположение пред-
1 2 г\ /' ч
д Л 1 / \ /у \ / 1 / > Г П ч /
'V 7 V о / л V \/19\ г V л Г/ 2> А у'вж
\ 1 \ \ V ^ 1
Рис. 5. Графики взаимных корреляционных функций: 1 - вывозка - раскряжевка; 2 - раскряжевка - вывозка
ставляется в данном случае более вероятным, поскольку всплески обнаружились при временных сдвигах, превышающих предельное значение 0 < Т/3 ... Т/6, рекомендуемое при использовании эмпирических коррело-грамм в практических расчетах [3, 5]. Вопреки этим рекомендациям расчеты коррелограмм размером 0 = Т/2 = 30 выполнены лишь для того, чтобы проследить, как приспосабливается в рассматриваемом периоде процесс к изменениям {X}; в расчетах же спектральных характеристик использовали корреляционное окно, равное Т/5.
Случайность всплеска в коррелограмме остаточного ряда ег подтвердил и более мощный в отношении гармонических составляющих критерий кумулятивной периодограммы Qk, применение которого графически показано на рис. 6. Как видно, все значения эмпирической кумулятивной периодограммы Qk укладываются в 5 %-е доверительные интервалы, при этом они
О 0,2 Ор 0,6 0,8 1М 0 0,2 0,4 0,6 0,8 /,о Частота полу периода Частота полупериода
а 5
Рис. 6. Кумулятивные периодограммы остаточных рядов: а - вывозка древесины; б - производство круглых лесоматериалов
достаточно тесно группируются вблизи биссектрисы первого координатного угла квадрата, которая представляет собой теоретическую кумулятивную периодограмму абсолютно стационарного процесса. Этот критерий однозначно указывает на то, что ряд остатков г2 не содержит существенных периодических составляющих, т. е. на 5 %-м уровне значимости обладает свойствами белого шума и может далее использоваться для применения к нему спектральных методов (с этого момента термин «гармонические составляющие» по отношению к остаточным рядам будет применяться в смысле их спектрального представления Фурье).
Спектральную плотность остаточного ряда производства круглых лесоматериалов оценивали с использованием п - первых автоковариаций по той же формуле, что и для остаточного ряда вывозки [6]:
п
Л К) = +1 £ ^е С (г)0С8 9 + 2 § ЯеСе (г)ос8 ^ 9, (10)
2%
% 0=1
% 0= 2
где - весовая функция Парзена [3]:
—9 =<
1 _ ^ -0| при 0 <8< n 211
n
при — < 0 < n
(11)
Спектральное окно проф. Парзена было принято нами потому, что его весовая функция отличается от других известных весовых функций очень незначительными боковыми выбросами (не более 0,2 % от основного максимума), а также дает для рассматриваемого спектрального окна всегда положительные значения —q [2, 3].
Оценки ко-спектров функций для частоты юк с учетом наложенных на взаимные ковариации условий (9) определяли для обоих остаточных рядов по формуле [3, 5]:
А/ 1 п
c(&k) = —[C0(х,z) + C0(z,х)] +--Ё—eCe(х,z) + Ce(z, х)] cosш9, (12)
4%
2% 9
а оценки квадратурных спектров на той же частоте по формуле:
q(&k ) = — Z —е [C (х, z) - Ce (z, х)] sin ш9 .
2% 9=1
(13)
Оценки ко-спектров и квадратурных спектров представляют собой характеристики двухмерного случайного процесса, реализацией которого является ряд пар {вх, вг}.
1
Частота полу период а
Рис. 7. График когерентности процессов вывозки древесины и производства круглых лесоматериалов
Показатель когерентности, характеризующий тесноту связи между гармоническими составляющими двухмерного процесса, определяли по формуле
С(„.) = , (14)
/х (ю* )Л (ю* )
где /х (ю* ), /2 (ю* ) - оценки спектральных плотностей остаточных рядов соответственно вывозки и производства круглых лесоматериалов.
По аналогии с коэффициентом множественной корреляции, значение показателя когерентности находится в интервале 0 < С(юк ) < 1 (рис. 7).
Можно констатировать, что на графике когерентности в рассматриваемом периоде времени в целом наблюдается достаточно сильная взаимная связь процессов на 5 %-м уровне значимости. Эмпирический показатель когерентности принимает значения от 0,461 до 0,755, причем наиболее сильная и устойчивая по амплитуде когерентность процессов обнаруживается в частотном интервале, соответствующем периодам продолжительностью примерно от 3 до 6 мес. Когерентность средней силы с меняющейся амплитудой имеет место в двух интервалах частот: с периодом более полугода (наименьшее значение взаимосвязи соответствует 1 году) и менее одного квартала, когда когерентность снижается до минимума.
Сдвиг фаз гармонических составляющих двухмерного процесса, являющийся аналогом сдвига уровней временного ряда одномерного процесса, определяли из выражений:
ф(ю *) = -агс^ ф(ю*) =
ГЧ К
|0 при с (ю*) > 0|
+ ол, о = -
с(ю*)) [1 при с(ю*) < 0]'
Г л /2 при с (ю) = 0 и д (ю) < 0| [- л/2 при с (ю) = 0 и д (ю) > 0]
Графически изменение угла сдвига фаз показано на рис. 8.
(15)
Рис. 8. График опережающего сдвига фаз процесса вывозки по отношению к производству круглых лесоматериалов
Из графика видно, что в интервале низких и средних частот, соответствующих периодам более одного квартала, технологический процесс вывозки опережает по фазе производство круглых лесоматериалов, при этом значение угла опережения сильно колеблется и в целом имеет тенденцию к снижению. Размах колебаний угла опережения уменьшается от максимума, проявляющегося с периодичностью 15 мес и равного 21° (что соответствует временному лагу примерно 4 мес), до нуля в области 3-месячной периодичности. Это показывает, что уровень загрузки мощностей по производству круглых лесоматериалов и темпы его роста в рассмотренном периоде по отношению к вывозке значительно выше. В частотах, соответствующих 3-месячному периоду, график угла фазового сдвига гармонических составляющих процессов осциллирует вблизи нулевого значения. Это говорит о недостаточности межоперационных и межсезонных запасов древесины, в результате чего основные технологические линии нижних складов вынуждены приспосабливаться к режиму работы «с колес».
В целом на уровне лесопромышленного комплекса Архангельского региона такой характер взаимозависимости говорит о стагнации, а в определенные периоды - даже о регрессе производственных мощностей по вывозке древесины. Можно предположить, что у конкретных лесозаготовительных предприятий, использующих более или менее совершенные транспорт-но-технологические схемы, ритмичность и надежность лесообеспечения мощностей по производству круглых лесоматериалов будут отличаться от тех, что здесь рассмотрены. Предлагаемая методика позволяет вычислить и сравнить эти показатели и для конкретных предприятий, определить наиболее рациональные схемы. Однако в связи с большим объемом расчетов это выходит за рамки данной статьи и является предметом дальнейших исследований.
Таким образом, исследование технологических процессов лесозаготовок с позиций теории случайных процессов позволяет установить степень их разупорядоченности, причем даже при одинаковых объемах производства, оценить количественно силу их взаимосвязи с другими фазами лесозаготовок с учетом нестационарности, а также объективно определить направ-
ление и продолжительность отставания уровней взаимозависимых динамических рядов в различные периоды времени, что невозможно сделать, анализируя показатели порождающих процессов обычными способами.
Современное состояние математического и аппаратного обеспечения расчетов позволяет исследовать не только двухмерные процессы, но и по аналогии с множественной регрессией устанавливать методом кросс-спектрального анализа зависимость изучаемого ряда от двух и более рядов одновременно. Продуктивность спектрального анализа при исследовании технологических систем настолько велика, что иногда говорят о новом способе мышления, а не просто о техническом достижении в статистических методах. Очевидно, что применение спектральных методов позволяет производить технологический анализ и прогнозировать развитие процессов лесозаготовок на новом, более высоком уровне.
СПИСОК ЛИТЕРАТУРЫ
1. Бокс Дж. Анализ временных рядов. Прогноз и управление / Дж.Бокс, Г. Дженкинс. - М.: Мир, 1974. - Вып. 2. - 197 с.
2. Грибанов Ю.И. Выборочные оценки спектральных характеристик стационарных случайных процессов / Ю.И. Грибанов, В.Л. Мальков. - М.: Энергия, 1978. -149 с.
3. Гренджер К. Спектральный анализ временных рядов в экономике / К. Гренджер, М. Хатанака. - М.: Статистика, 1972. - 312 с.
4. Кендэлл М. Временные ряды / М. Кендэлл; пер. с англ. Ю.П. Лукашина. -М.: Финансы и статистика, 1981. - 199 с.
5. Маленво Э. Статистические методы эконометрии / Э. Маленво; пер. с франц. - М.: Статистика, 1976. - Вып. 2, 4. - 325 с.
6. Меньшиков А.М. Анализ динамических рядов транспортно-технологи-ческих процессов вывозки древесины (на примере Архангельской области) / А.М. Меньшиков, А.М. Копейкин; ОАО «Научдревпром - ЦНИИМОД». Лесн. и дерево-обраб. пром-сть. - Архангельск, 2003. - 40 с. - Деп. в ВИНИТИ 15.12.2003, № 2177-В2003.
ОАО «Научдревпром - ЦНИИМОД» Архангельский государственный технический университет
Поступила 22.01.04
A.M. Menshikov, A.M. Kopeikin
Use of Spectral Methods in Engineering Process Studies of Forest-harvesting Production
The internal structure and dynamics of trend-seasonal fluctuations of time-dependent engineering processes of wood-harvesting and round-wood production by the forest-industrial enterprises are analyzed based on the random functions' theory. Technique of quantitative determination of force, operation duration and interconnections directivity for parallel technological processes of forest harvesting is offered based on cross-spectral analysis of their dynamic rows.
