АКТУАЛЬНЫЕ НАПРАВЛЕНИЯ ИССЛЕДОВАНИЙ В ОЦЕНКЕ ТЕХНОГЕННОГО РИСКА ПРИ ЭКСПЛУАТАЦИИ ОБОРУДОВАНИЯ Текст научной статьи по специальности «Математика»

УДК 519.87

АКТУАЛЬНЫЕ НАПРАВЛЕНИЯ ИССЛЕДОВАНИЙ В ОЦЕНКЕ ТЕХНОГЕННОГО РИСКА

ПРИ ЭКСПЛУАТАЦИИ ОБОРУДОВАНИЯ

В. А. Острейковский

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

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

Ключевые слова: техногенный риск, исходное событие, ущерб от исходного события, модели оценки фактора времени.

CURRENT RESEARCH AREAS IN THE EVALUATION OF TECHNICAL RISK IN THE

EQUIPMENT OPERATION STUDIES

V. A. Ostreykovsky

Surgut State University

Generalization of results of the analysis of mathematical models of technogenic risk is executed. It is shown that in literature extremely a fractional information on analytical models of an assessment of a factor of time in the theory of technogenic risk. Mathematical methods and analytical models for taking note of a factor of time in problems of a quantitative assessment of technogenic risk of difficult dynamic systems are offered.

Keywords: technogenic risk, initial event, damage from an initial event, the model estimates the time

factor.

Введение

В довольно обширной литературе по теории риска мало внимания уделено моделированию фактора времени. Это частично объясняется, если вспомнить, что теория риска сложных динамических технических систем бурно развивалась с классической формулы Ф.Фармера R = QC, где и вероятность исходных событий аварий и катастроф Q и ущерб от них C понимались как множества, а не зависящие аналитически от времени. В настоящее время это выглядит архаично, но в шестидесятые годы прошлого века на учет времени в функции риска не было обращено должного внимания. Хотя уже и тогда четко оговаривалось что характеристики надежности и безотказности сложных систем суть функции времени. Поэтому в данной статье задача учета фактора времени при моделировании риска R рассматривается как одна из центральных в теории техногенного риска.

Постановка задачи

В работах [1-3,7] разработана классификация моделей техногенного риска сложных динамических систем, в которых детально сформулирован теоретико-множественный подход к анализу характеристик множеств R, Q и C при математическом моделировании техногенного риска. В общем виде графическая интерпретация соотношения множеств риска R, вероятностей исходных событий (отказов, аварий, катастроф) Q, ущерба от них C во времени T приведено на рис.1

В соответствии с агрегативным описанием функционирования сложных систем набор операторов H имеет вид:

R = H{Q*C*T} , (1)

где: H : H\ * Щ* H<$, реализуют следующие отображения:

H\ : Q * C R\ при независимых между собой Q и C;

H2 : Q * C/Q —> R2 при коррелированных Q и C; H3: Q * C * T R

Рис. 1. Множества Я,С и Т в моделях техногенного риска при количественной оценке безопасности сложных систем

Рис. 2. Ущерб как случайный процесс

При статистическом исследовании составляющие Я(1), Ц (1) и С(I) рассматриваются как случайные функции времени, в общем случае векторные. Если представить риск как возможный ущерб (рис. 2), то модель Я = /(Т) времени должна обладать следующими свойствами:

1) процесс изменения во времени Я = / (Т) рассматривается в общем случае как нестационарный стохастический процесс;

2) длительность наблюдения за системой в течение Т конечна;

3) энергетический спектр случайного процесса Я(Т) сплошной и отличен от нуля на всей оси частот

^оо < ш < оо;

4) интервал корреляции то случайной функции Я = / (Т) ограничен, причем то макс << Т. Под интервалом корреляции принимается промежуток времени, за который полностью затухают корреляционные связи между отдельными составляющими Я(1), Ц (^ и С(¿). При нестационарном характере изменения случайной функции Я = / (Т) интервал корреляции в общем случае зависит от рассматриваемого момента времени, т.е. то = /(¿).

Нестационарные случайные процессы изменения Я = / (Т) целесообразно представлять в виде суммы нескольких процессов Я = / (Т):

Я Ц) = А Ц) + В У) + еЦ),

где А (^ - нестационарный случайный процесс, характеризующий необратимые изменения в системе в результате старения, изнашивания, регулирования; В (¡) - стационарный случайный процесс, харак-

теризующий обратимые изменения Я (г) из-за колебаний внешних условий при эксплуатации системы; е (г) - стационарный случайный процесс ошибок измерений Ц (г) и С(г).

Рассмотрим методы оценки фактора времени, приведенные в [4], при планировании многофакторных испытаний на надежность сложных критичных систем (параметра Ц).

Основная часть

1. Метод фактора времени при моделировании риска в условиях ортогонального дрейфа.

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

Я = М [С/Ц] = У(г) (2)

может быть представлена полиномом невысокой степени. Используя априорную информацию о характере дрейфа, можно исключить его влияние, планируя наблюдение за системой ортогонально к дрейфу. В этом случае задача моделирования сводится к построению плана, обеспечивающего получение наилучших оценок эффектов управляемых факторов, которые были бы ортогональны эффектам дрейфа. Многофакторные планы в условиях дрейфа можно строить на базе любой подходящей ортогональной системы функций. Математическая модель риска Я = / (г) представляется в виде разложения по выбранной системе ортогональных функций, часть их используется для описания дрейфа, а часть — для варьирования управляемыми факторами, т.е. образования плана. Предложенный Боксом [4] метод построения планов для оценки поверхности отклика в условиях дрейфа основан на использовании полиномов Чебышева, Сущность метода заключается в следующем. Функция дрейфа Яг = У(г) на интервале Т описывается полиномом порядка й:

Яг = аоЬо (Т) + а\Ь\(Т) + ■■■ + айЬй (Т)

(3)

где Ь^ — ортогональные полиномы Чебышева, определенные для N = Т/Д г равностоящих абсцисс наблюдений интервала Т.

Чтобы выделить влияние управляемых переменных Ц и С независимо от дрейфа, вектор-столбцы планирования формируют в виде линейных комбинаций всех оставшихся полиномов порядка £ > й.В этом случае модель имеет вид

Я (Ц,С ,Т) = Я (X ,Т) = аоЬо (Т) + а\Ь\ (Т ) + ■■■ + айЬй (Т) +Ь\х\ + Ь2х2 + ■■■ + Ькхк,

4-V-'

Ьо(г)

(4)

где X С ЦхС; Ьо(г) определяет смещение поверхности отклика во времени, а коэффициенты Ь — ее наклон в исследуемой области Я(Ц,С). Анализ уравнения (4) сводится к оценке значимости эффектов дрейфа и управляемых факторов по F-критерию. Каталог планов для наиболее употребительных значений N d, k приведен в [4].

Применение ортогональных полиномов Чебышева особенно эффективно, когда отсчеты времени делаются через равные интервалы Д г, а для каждого значения г ставится одно и то же число параллельных опытов. Шкала времени при этом преобразуется в последовательность натуральных чисел l,N. В этом случае ортогональные полиномы равны [4]:

Фо = \,

Ф\ = {г-7)

Ф2 = № Фз = № Фа = Ш

(г-1) - (ы2^ \)/\2

{г-г) ^(зы2 г)/20

(т^т)4 - (зы2^ \з)(г-т)2/\а + з ы2-\) ш2^9)/560

фъ = иъ{ {г-г)5 ^5(7ы2 -7)( т-Т) 2/\8+ - 230ы2 + 407^ {г-г)/\008

Коэффициенты //1 зависят от числа наблюдений ЖПолиномы (5) табулированы [4].

Недостатки метода Бокса:

1) планы, основанные на применении ортогональных полиномов Чебышева, могут быть использованы только в том случае, когда отсутствуют корреляционные связи типа «время—фактор»;

2) построенные на его основе планы позволяют оценивать независимо от дрейфа только линейные эффекты варьируемых факторов Ц и С. Следовательно там, где взаимодействиями факторов и эффектами высших порядков пренебречь нельзя, применение метода Бокса не будет давать полной информации о системе;

3) необходимо, чтобы на стадии предварительных исследований был выявлен порядок полинома дрейфа d, что само по себе эквивалентно проведению громоздкого эксперимента.

2. Метод оценки фактора времени путем преобразования параметров регрессии из функций времени в числовые коэффициенты.

Функцию отклика Я(Ц,С,Т) системы в условиях динамической задачи можно представить в

виде:

Я (Ц,С ,Т) = Яо (Ц,С ,Т) + £ Ь (!) х1 + £ Ьи (!) хх + £ Ьц (!) х2 + ... (6)

¡=1 К] 1=1

Целесообразно преобразовать выражение (6) в уравнение с числовыми коэффициентами регрессии. Принципиально это можно сделать несколькими путями.

Рассмотрим один из методов для случая линейной модели:

к

Я (Ц,С,Т) = Яо (Ц,С,Т) + £ ь (!) X. (7)

1=1

Представим Яо (Ц,С,Т) и Ь (!) в виде функциональных рядов

N

Яо(Ц,С,Т)=^ ОоаФа (!) , (8)

и=1

N

bi (t) = Е ^иФи (t) ,

и=1

где 0и (t), и = 1,N - коэффициенты функции; Cou, Cou, u = 1,N, i = l,k - постоянные коэффициенты. Тогда

N k N

R (Q,C ,T ) = J2 соифи (t) + E E СшФи (t) Xi. (9)

0u u

u=l i=l и=1

Проинтегрируем выражение (9) справа и слева в пределах (0;7):

T N t k N t

R (Q,C,T) dt = V cou Фи (t) dt+YY] Ciu I фи (t) dtXi.

0 u=i Jo ti Jo

Введем обозначения:

Г R (Q,C,T) dt = Z, c0u = cou Î Фи (t) dt, c*iu = f фи (t) dt.

Jo Jo Jo

Тогда получим выражение:

N k N

z = E c0u + EE ciuXi. (io)

0u

u=l i=l u=l

В этом уравнении коэффициентами регрессии являются величины c0u и c*u, i = l,k и = l,N.

Оценки этих коэффициентов определяются соотношениями:

Сш —

Т

/о Фа У) йг

(11)

Функции фа (г) и значение Т выбираются такими, чтобы

Т

Фа (г) йг — о,

о

а — \,М.

3. Модели оценки фактора времени в эксперименте, когда одна из контролируемых переменных - время

Представим выражение для функции отклика в виде [4]:

я (Я,С,т) — [в (Я,С,В) + V (Я,С)Г/ (г) + е (г)

(12)

где / (г) ив* (Я,С ,В) — ||0\ (Я,С ,В),.. . ,вт (Я,С ,В)|| - известные функции; В - вектор неизвестных параметров; V (Я,С) - функция отражающая различие в системах, причем

М [V (Я,С)]—0, М [е(г)\— 0

' Dv (Я,С), Я,С — (Я,С)

М^ (Я,С) v(Q,C)'] —

0,

Я,С — (Я,С)

М [е(и)еШ— Кя г Ь),

(13)

где V (Я,С) и е (г) - между собой независимы; Кя (г\,г%) - корреляционная функция риска.

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

В качестве меры точности могут служить любые критерии оптимальности, например определитель дисперсионной матрицы оценок параметров \D(B)\ или сумма дисперсий этих оценок SpD(B). Случай 1.

Кд (г\ ь)</ (г) Dv (Я,С) Г (г), г е [0,Т] (14)

или

я (Я,С,Т) — [о (Я,С,В) + V (Я,С)]*/ (г).

(15)

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

Имея реализацию Я [(Я,С)^, г) можно определить значения в [(Я,С)^,Щ + VI при заданных (Я,Спо формуле:

в [(Я,С),В] + VI —

минимизируя интеграл

/ (г) Г (г) йг

я [(Я,С)гд]! (г) йг,

(16)

{\е [(Я,Су,В] + vi] 7 (г)^ я [(Я,Су, г) }2йг

При конечном числе замеров реализации я [(Я,С)i, г) вместо выражения (16) имеем:

в [(Я,С)i,В] + Vi —

Количество замеров / и положение точек должны быть таковы, чтобы матрица

£/(М / (Ь)

¡—\

^я [®,Су^Цг,).

1—\

А — Е / (М )Г (^)

!—\

С

Т

О

О

Т

О

была неособенной. Для уменьшения погрешностей вычислений последовательности , ] = 1, m следует располагать так, чтобы матрица А* была плохо обусловленной. Таким образом, при каждом I] можно получить случайные значения в [(,С)1,Щ + Vi и задача становится эквивалентной следующей регрессионной задаче:

7 = в(,С,В) + v (Ц,С).

Матрица Ь (Ц,С) предполагается известной. На основании измерений случайной величины X получаются оценки поверхностей 01 (Ц,С,В),.. . ,0т (Ц,С,В) или оценки параметров В. Если зависимость от неизвестных параметров линейна

да (Ц,С ,В) = В фа (,С) ,

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

В = М1Щ, (18)

где

М = £ В [Ю,С)] О-1 [(Ц,С)] Г [(,С)] ,

i=1

г = Е В [((,С1 [(,С\]х,, i=1

В [(,С= \\ф1 [(,С,...,фт [((,С)г]||.

Случай 2. Неравенство (14) не выполняется, но корреляционная функция кя ^1^2) известна. В этом случае моделирование для регрессионной задачи (12) так же удается свести к моделированию в случае одновременного измерения нескольких величин [4].

При известной функции кя (11^2) наилучшие линейные оценки для величин = в [(,С)i,В] + VI равны:

X = А\~1Щ,

(19)

где

т

А1=Т, нь)Кя 1 ],

!,!' = 1 т

Щ1 = у^ ¡{1,)Кя 1 [1]]Я [(,С)],

¡,г=1

Причем дисперсионная матрица оценок Х^ при заданном VI обратна матрице А1

(20)

Из выражения (20) следует, что

Ь(Х1/VI )= АТ1.

X =в[(Ц,С ) ,В] + Vi +11 i.

(21) (22)

При этом Ь (цI) = А1 1, VI и ¡лI независимы, и

м[г]=е [(,с)В , ь [X] = о [(,с)] + А\~х.

Таким образом, задача сведена к случаю 1 лишь после произведения замены Ьц [(,С)] на Ьц [(,С)] +

А1 1 .

Случай 3. Неравенство (14) не выполняется, корреляционная функция Кя ^1^2) неизвестна. В этом случае непосредственное применение формулы (19) невозможно. Однако во многих случаях имеется информация о скорости убывания |Кя ^1^2^. Если измерения проводятся через промежутки времени, за которые |Кя (становится близким к нулю, а амплитуда внешних помех не зависит

от I, т.е. Кд (^¿2) = сгд (это требование обычно выполняется при не слишком больших интервалах [0,7]), то выражение для А1 и Щпринимают вид:

А*«^Е Ш/* ^), (23)

1=1

т

Щ^^Е^СЪ1 , (24)

1=1

Если значение ст2неизвестно, то в выражении (23), (24) следует подставлять его оценку

2 1 Д №,с))/(1 )]2 (25)

5 = —--, (25)

т — п

где т — п - число степеней свободы.

В тех случаях, когда неизвестна матрица 0О ,но известно, что она одна и та же во всей области (Ц,С), можно воспользоваться оценкой для

г

(п - 1) й {%) = Е ¿г) (¿г,- %) *,

1=1

где ¿г, = г 4=1 ¿п; гI - число реализаций при заданном (Ц,С.

Точку ^,С),стремятся выбрать так, чтобы она принадлежала оптимальному плану.

Изложенным методом можно воспользоваться и в тех случаях, когда ,Т) зависит от

9(Ц,С,В) нелинейно. Это характерно для многих задач теории надежности и безопасности, в которых функция 9(Ц,С,В) является параметром того или иного закона распределения. В этих случаях оценки могут отыскиваться любым методом, подходящим в данной экспериментальной ситуации. Необходимо лишь знание дисперсионной матрицы й^^о,), соответствующей использованному методу [4].

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

4. Метод оценки фактора времени с применением спектрально-корреляционной теории случайных процессов

Продолжительность наблюдений можно рассчитать по виду реализации случайного процесса риска Я(0, принимая определенные предположения. Обычно принимается, что стационарный случайный процесс имеет нормальное распределение координат, корреляционная функция его может быть аппроксимирована суммой экспонент. Корреляционная функция стационарного процесса обычно имеет монотонный характер. При подходе корреляционной функции к оси абцисс колебания носят случайный характер. Флуктуации корреляционной функции обусловлены в основном конечностью длины реализации Поэтому случайные колебания, наблюдающиеся на «хвосте» корреляционной функции, могут быть исключены из рассмотрения.

Рассмотрим случай с линейной дискретной моделью риска Я(0, который задается Фурье-разложением по конечному числу частот:

т

Д (() = Е (а1 + ЬsmuJjt ), (26)

1=1

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

= М

а2 = М

которая является функцией частоты .

Корреляционная функция случайного процесса У равна

т

М [Я (!) Я(г + т)] = Кя (т) = о^соьш/т. (27)

1=1

Согласно теореме Хинчина—Винера [8] спектр случайного процесса и корреляционная функция являются сопряженными Фурье-функциями. Поэтому

1 г ^1

Кя (о)+ 2^ Кя (!)Кя (! + т)

=1

где сг2 и Кя - выборочные оценки. Для выборки объема N

2

°' = 2п

(28)

1 N-7

Кя (т) = ^^ Е Кя (!) Кя (! + т), (29)

Т Т=1

где принимает лишь дискретные значения.

Особенностью спектрального анализа является его непараметричность: метод не позволяет оценить параметры кривой спектральной плотности, ее приходится оценивать по каждой ординате отдельно. При сглаживании и обеспечении состоятельности оценки спектра (или спектральной плотности) оценки становятся смещенными. Обычно сглаживание производят, когда выполняют Фурье-преобразование корреляционной функции. В качестве сглаживающей функции (фильтра) хорошо известно усеченное сглаживание с помощью спектрального «окна» Бартлета

1

а ¡ = 2^

кя + а кя (т)

где = 1 /т, 1 т.

Остановимся на способах выбора длительности наблюдений Т и интервала съема данных Д В литературе по теории случайных процессов имеется ряд соображений по оценке необходимой длины реализации при записи случайного процесса. Оценка величины Т производится непосредственно по виду реализации Я(V) при принятии ряда предположений об исследуемом случайном процессе.

Для нормального случайного стационарного процесса средняя квадратичная погрешность вычислений корреляционной функции, возникающая из-за конечной длины реализации, при условии М [Я(01 = тг (!) = о равна [9]

2 рТ-т

°КВ = ^-72 (Т-т-в) Ц(в) + кя (в + т)+ кя (6-т) й (в), (30)

Я (Т ) о

где кя — корреляционная функция ВПО Я(!); Т — длина реализации (время наблюдения); — временной сдвиг; — приращение временного сдвига.

Обычно практически наиболее важен участок корреляционной функции при небольших значениях т, поэтому для дальнейшего расчета будем пользоваться окШ(У) = ®Ь —средней квадратичной погрешностью вычисления дисперсии случайного процесса

4 ГТ

аЬ = ^ (Т- в) КЯ (в) М, (31)

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

Так, для случайных процессов риска, характеризуемых корреляционной функцией

кя (т) = Ье~хм (32)

в [4] для значений Т 100 получена следующая приближенная аналитическая зависимость

ак< 2лЬ, (33)

где Ь - заданная погрешность вычисления дисперсии корреляционной функции, возникающей из-за конечной величины реализации Кд (т).

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

" ^ (34)

Кк(т) = Ь (ае-с^ - Ье-^

В случае монотонного характера кривой оценки корреляционной функции вычисления часто прекращаются при ттах, определяемом из условия

К(ттах) < 0,05 К(0)1 .

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

Кд (т) = Ь (4^^ - 1 е-^ (35)

приближенная зависимость (34) имеет вид:

Г)2

оЬ* 3—, (36)

где с - параметр формулы аппроксимации, характеризующий наклон функции.

Введя понятие относительной средней квадратичной погрешности определения дисперсии

ГЮ = вй/й и подставив его в формулу (36), получим

3

Т =-. (37)

Щв

Таким образом, задаваясь значением г]в и вычисляя по виду реализации с можно по (37) предварительно оценить значение Т времени наблюдения для оценки значений функций Q (Т), С (Т) и д(Т).

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

Т = ^, (38)

где 7 = тс; р = ^е-7 = — е-7)2.Требуемое значение времени наблюдения Т находится, если по выбранным значениям Дt и вероятности р определить параметры 7 и с. Заключение

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

1) Установить равную для всех режимов функционирования системы продолжительность наблюдения;

2) Через равные промежутки времени рассчитывать значения вероятности исходных событий ущерба от них для каждого режима эксплуатации системы;

3) По результатам расчета определять константы скоростей изменения вероятности исходных событий и ущерба во времени;

4) По полученным значениям констант изменения риска строить модель, связывающую константы с действующими факторами.

Получение наиболее полной информации о закономерностях изменения значений риска во времени при комплексном воздействии на систему факторов внешней среды возможно после интерпретации полученной таким образом зависимости константы скорости от действующих факторов. Работа выполнена при поддержке РФФИ (проект 14-01-00230)

ЛИТЕРАТУРА

1. Королев В.Ю., Бенинг В.Е., Шоргин С.Я. Математические основы теории риска: Учеб. Пособ. М. : ФИЗМАТЛИТ, 2007. 544 с.

2. Острейковский В.А. Теория техногенного риска: математические методы и модели: монография, Сургут. гос. ун-т ХМАО - Югры. Сургут : ИЦ СурГУ, 2013. 320 с.

3. Острейковский В.А., Шевченко Е.Н., Микшина В.С. Количественная оценка риска в техногенной безопасности сложных динамический систем: монография // Итоги науки. Том 1. Избранные труды международного симпозиума по фундаментальным и прикладным проблемам науки. М.: РАН, 2013. С. 12-31.

4. Острейковский В.А. Многофакторные испытания на надежность. М.: Энергия, 1978. 152 с., ил.

5. Острейковский В.А. Теория надежности: Учебник для вузов 2-е изд., испр. М. : Высш. шк., 2008. 463 с.: ил.

6. Жакот А.Д., Острейковский В.А., Челнокова Е.В. Математические модели регрессионного анализа и теории катастроф синдрома дыхательных путей / Надежность и качество. Тр. междунар. симпозиума / под ред. Н.К.Юркова. Пенза : Изд-во ПГУ, 2005. с.467-472.

7. Острейковский В.А. О некоторых классах моделей риска в теории техногенной безопасности // Надежность и качество. Тр. междунар. симпозиума: в 2 т. / под ред. Н.К.Юркова. Пенза : Изд-во ПГУ, 2013. т. 1. с. 46-49.

8. Гихман И.И., Скороход А.В. Введение в теорию случайных процессов. Изд. 2-е, Главная ред. физ.-мат. лит. изд-ва «Наука», М., 1977. 568 с.

9. Пугачев, В.С. Теория случайных функций и ее применение к задачам автоматического уравнения. Изд. 3-е, испр. М. : Гос. изд-во физ.-мат. лит., 1962. 883 с.