МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ
Доказательство. Достаточно заметить, что каждая координатная функция порождаемого отображения является самодвойственной.
Пример. Построим подстановку по квазиадамаровой матрице
(1110 ^
A - 1 -10 -1 10 -1 1 ■
10 -1 1 1J
Непосредственной проверкой можно убедиться в том, что координатные функции порождаемой подстановки имеют вид f1(x1, x2, x3, x4) = 1 о x1 + x2 + x3 > 3/2,
f2('X\, ^ X3, X4) = 1 0 X1 + (1 -X2) + (1 - X4) - 3/2, f3(x1, x2, x3, X4) = 1 о xj + (1 -X3) + X4 - 3/2,
f4(x1, x2, X3, X4) = 1 о (1 -X2) + X3 + X4 - 3/2. Сама подстановка в двустрочной записи имеет вид
_ /0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15\ 8л 3 5 1 0 2 8 9 6 7 13 15 14 10 12 11/
где двоичные вектора пространства V4 естественным образом представлены как целые числа. Легко убедиться, что подстановка gA
раскладывается на 8 независимых циклов
длины 2 и является обратной к самой себе.
Работа выполнена при поддержке гранта Президента РФ (НШ-4.2008.10)
Библиографический список
1. Логачев, О.А. Булевы функции в теории кодирования и криптологии / О.А. Логачев, А.А. Сальников, В.В. Ященко. - М.: МНЦМО, 2004.
2. Мак-Вильямс, Ф.Дж. Теория кодов, исправляющих ошибки / Ф.Дж. Мак-Вильямс, Н.Дж.А. Сло-эн. - М.: Связь, 1979.
3. Никонов, В.Г. Методы компактной реализации биективных отображений, заданных регулярными системами однотипных булевых функций / В.Г. Никонов, А.В. Саранцев // Вестник РУДН, серия Прикладная и компьютерная математика. - 2003.- Т. 2. - № 1. - С. 94-105.
4. Носов, В.А. Специальные главы дискретной математики / В.А. Носов. - М., 1990.
5. P. Delsarte, J.-M. Goethals, J.J. Seidel. Orthogonal matrices with zero diagonal II // Canadian Journal of Mathematics, 23 (1971) 816-832.
6. J.-M. Goethals, J.J. Seidel. Orthogonal matrices with zero diagonal // Canadian Journal of Mathematics, 19 (1969) 1001-1010.
НОВАЯ ДИСКРЕТНАЯ БАЙЕСОВСКАЯ МОДЕЛЬ НАДЕЖНОСТИ
программного обеспечения с использованием интервальных вероятностей
С.И. ЗАТЕНКО, ст. преподаватель каф. высшей математики СПБГЛТА,
Л.В. УТКИН, проф. каф. информатики и информационных систем СПБГЛТА, д-р техн. наук
В настоящее время в условиях развития и распространения сложных информационных систем, основной частью которых является программное обеспечение (ПО), качество и надежность становятся одними из важнейших характеристик программных средств. Анализ надежности ПО - один из наиболее важных этапов разработки программ и программных систем, так как именно от наличия ошибок в программе зависит эффективность ее дальнейшего использования. Сокращение времени разработки программ, включая их тестирование, существенно усложняет моделирование надежности ПО. Поэтому разработка моделей надежности, позволяющих учитывать недостаток статистической информации, является актуальной задачей [1].
s_lana2004@mail.ru, lev.utkin@mail.ru
Существует два основных подхода для моделирования надежности ПО: подход, основанный на нахождении максимума функции правдоподобия, и байесовский подход. Однако ни один из них нельзя считать универсальным. При использовании первого необходимо иметь большое количество наблюдений, кроме того, он плохо работает в случае большого количества параметров. Проблема, которая возникает при использовании второго, заключается в том, что произвольный выбор априорного распределения при малом количестве статистической информации в достаточно сильной степени может влиять на результаты прогнозирования.
В данной работе рассматривается новый метод моделирования надежности, кото-
158
ЛЕСНОЙ ВЕСТНИК 2/2009
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ
рый сочетает байесовский подход и принцип максимума правдоподобия. В результате этого сочетания решается проблема произвольности выбора априорного распределения, кроме того, получающаяся задача оптимизации функции правдоподобия становится однопараметрической. Используя этот метод, мы создаем новую интервальную модель надежности ПО дискретного времени.
Рассмотрим программное обеспечение на этапе его тестирования. Любой процесс функционирования ПО может быть разделен на последовательность запусков. Поэтому часто удобно считать, что времена между отказами являются дискретными случайными величинами, т.е. один запуск программы и ее выполнение является минимальной единицей дискретного времени безотказной работы. В результате тестирования и отладки получим поток случайных чисел k,...,k который мы будем рассматривать как набор значений случайных величин X Хп, где X - дискретная случайная величина, принимающая значение k, если программа отказала в течение k-го выполнения после (k - 1) успешных запусков и выполнений. После обнаружения отказа в ПО возможная ошибка, вызвавшая отказ, исправляется, что изменяет характеристики программы и соответственно случайных величин X Хп- Конечной задачей моделирования надежности ПО является нахождение вероятностных характеристик случайной величины X где Xn+1 - время безотказной работы программы после окончания отладки, т.е. в процессе эксплуатации.
В работе авторов [2], был представлен целый класс моделей, который можно использовать как инструмент для создания самых различных частных интервальных моделей. А также в рамках созданного инструментария было доказано утверждение
max max L(K 19) =
e r^.xr
= maxП№(k 19)-F,((k-1)19)}.
i=1
Здесь F_i(k | 9), Fi(k 19) - нижняя и верхняя функции распределения вероятностей случайной величины Xi, ограничивающие множество возможных распределений R.
Трудность, которая возникает при использовании полученных результатов, заключается в определении множества распределений или граничных функций распределения случайных величин Xi, а также в наличии большого числа параметров, по которым осуществляется максимизация функции правдоподобия.
Основная идея моделирования состоит в том, что все множество параметров модели мы разделим на два непересекающихся подмножества, D1 и D Подмножество D состоит из параметров, определяющих множество всевозможных функций распределения F(k), а второе подмножество - это параметр, характеризующий рост надежности ПО в процессе исправления ошибок. Применяя байесовский подход для нахождения параметров распределения принадлежащих D1, получаем последовательность функций распределения F(k\b,d), где b,d - векторы параметров, b е D1 и d е D2. Далее, фиксируя параметр d, находим нижнюю и верхнюю функции распределения F (k I d) и F (k | d), такие что F (k | d) < F (k | b, d) < F (k | d). После этого вектор d вычисляется с использованием принципа максимума правдоподобия.
Следующей основной задачей является определение конкретного вида F (k | d), F (k | d). Для этого используем обобщенный байесовский подход [3].
Построение интервальной модели надежности ПО на основе бета-геометрической байесовской модели
Пусть X дискретные случайные величины, имеющие геометрическое распределение P{Xt < k} = F(k\p) =
= Е Рг (1 - Рг ) 1 -1 = 1 - (1 - Рг f ,
j=1
где р .- параметр распределения или вероятность неуспешного запуска программы между (г - 1)-м и i-м отказами. Таким образом,
р(Х = к) = р(к) = Рг(1 - p)kA.
Байесовский подход к моделированию заключается в следующем. Предположим, что параметр р сам является случайной величиной с некоторой функцией плотности вероятности n(pi\(a.,b)), здесь a.,b. - параметры
ЛЕСНОЙ ВЕСТНИК 2/2009
159
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ
распределения. Для согласования априорного и апостериорного распределений, функция п должна быть бета-распределением, обозначим его Beta(a.,b.)
п(р) = Beta(a.,b.) =
= \ , РУ - (1 - Рг )Ь‘^ ,0 ^ Рг < 1.
B(ai, Ьг )
Тогда модифицированная функция
распределения X будет иметь вид:
1 1
F (k) = j (1 - (1 - Рг У )Beta(a., b, )dPi =
= 1 -
Г(аг + Ьг ) Г(Ьг + k)
Г(Ьг ) Г(аг + Ьг + k)
Здесь Г(х) - стандартная гамма-функ-
ция.
Если появилась информация о времени безотказной работы на i-м этапе тестирования, то, пересчитывая параметры на основе известной формулы Байеса, получим апостериорную плотность распределения
п(рг+1) = n(p.|k.) = Beta(a + 1, b + k. - 1), где (k . - 1) - фактическое число успешных запусков между (. - 1)-м и !-м отказами.
Пусть k ,k - реализация X Хп-До проведения эксперимента распределению случайной величины X1 будет соответствовать параметр р0 с плотностью распределения вероятности п(р0) = Beta(a;P). Значения параметров, бета-распределения (a; в) неизвестны. Поэтому предлагается вместо какого-либо конкретного распределения рассмотреть весь класс возможных распределений на множестве параметров [3,4].
Будем считать, что параметры (a; в) изменяются вдоль отрезка, граничные точки которого находятся в точках с координатами (s;0),(0;s), s > 0. Выбор отрезка обусловлен тем, что математическое ожидание вероятности отказа для бета-распределения равно Мр0 = a / (a + в). А значит, изменяя параметры (a; в) данным образом, получим значения Мр0 е (0;1). Здесь s - параметр осторожности, определяющий размер треугольника и соответственно границы множества распределений .
Пересчитывая параметры распределения Р по мере поступления статистической информации, получим на -м этапе тестирования
п(р) = Beta(a +г; в + (/ - Щ + 2 (kj -1)),
j=1
где в1 - параметр, характеризующий рост надежности ПО в процессе исправления ошибок.
А значит, модифицированная функция распределения случайной величины X будет иметь вид
F (k) = } (1 - (1 - р, )k) х
0
i
xBeta(a +1; в + (/ - 1)в1 + 2 (kj - 1))dP =
f
j =1
Г
= 1 --
a+.+в+(i' - 1)в1+2 (kj-1)
j=1
Л
f
Г
в+(/- 1)в1 +2(kj -1)
Л
j=1
f
Г
в+(/- 1)в1 +2(kj -1)+k
j=1
л
Г
a +1 + в + (. - 1)в1 + 2 (kj -1) + k
X
j=1
Обозначая
K=2 (*, - !)■
j=1
где K - суммарное число успешных запусков до .-го отказа, и D . = K. + (г - 1)в1 получим
F (k) = 1 Di) Г(в + Di + k)
Г(в + D) r(a +. + в + D + k) '
Фиксируя параметр роста в1, найдем минимальное F. (k | в1) и максимальное F (k | в1) значения функции F.(k).
Докажем, что минимальное значение функции распределения F (k) достигается в том случае, когда параметры (a; в) = (0;s), а максимальное при (a^) = (s;0).
Если a; в принадлежат отрезку, концы которого имеют координаты (0;s) и (s;0), тогда a = ss, а в = s - ss, где s изменяется от [0;1]. Рассмотрим
F (k;s 1 в1) =1 -
r(ss +. + s - ss + D )
r(s - ss + Di)
Г(s - ss + D + k)
x----------------------
r(ss +. + s - ss + D + k)
и докажем, что F.(k; s^) возрастает с ростом s.
Заметим, что F (k; s|в1) может быть записана как
х
х
160
ЛЕСНОЙ ВЕСТНИК 2/2009
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ
F,(k; s1 Pi) =1 -
Beta( s + i + Di, k)
k-i
Beta(s - ss + Dt, k)
А значит, F.(k; s|P,) возрастает на [0,1], так как Beta(xy), убывает с ростом х, когда х > 0.
Таким образом F,(k | Pi) = if F (k, s | Pi) =
= f (k ,0|Pl) = 1 - Beta(s+i+D- k >
__ Beta(s + Dt, k)
F (k \ Pi) = sup F (k, s | Pi) =
0<s<1
= F,(k,1\в) = 1 -Beta(s +' + Dk>.
г 1 Beta(D,, k)
Определим теперь параметр p находя максимум функции правдоподобия.
max L(K \ Pi) =
Pi
= m?xOIF(k, |Pi)-Fi(k, - 1|в)} =
e,
i = i
= max П
Pi 11
i=i
Beta(s + i + Di, k, -1)
Beta(s + Dt, kt -1) Beta(s + i + Dt, kt) ^
V
Beta( Di, kt)
Так как k, целые числа, то для упрощения данного выражения можно использовать следующее равенство
о . ч (n -1)!
Beta (a, n) = —-----.
* П (a + i)
А значит,
’ i=0
in (s + Dt + j) n(D, + j)
Ц (K\b i) = j=0 j=0 =
П (s + i + Di + j) П (s + i + Di + j)
j=0 j=0
(s + i + Di + ki - 1)П (s + Di + j) - П (Di + j)
j =0
j=0
П (s + i + Di + j)
j=0
Логарифмируя функцию правдоподобия, получим равенство
ln Ls) (K \ p i) = £ (-£ ln(s + i + Dt + j) + ln Gt),
i=1 j=0 где
Gi = (s + i + Di + ki - 1)П (s + Di + j) - П (Di + j) .
j=0 j=0
Заметим, что в случае, когда параметр s =1, функция правдоподобия значительно упрощается
(i + k )П (D + j)
L (k\P, ) =
ki -i
П (1 + i + Di + j )
j=0
Логарифмируя теперь функцию правдоподобия при s = 1, получим
ln L"( K i Pi) = £ (ln(i+k,) +
i=
+£ ln( D, + j)-£ ln(1 + i + D, + j)).
j=1 i=0
Тогда
maxln L(1)(K \ P i) =
f
Pi
k, -i
= mRax £ ln(i+К)+£ln(K, +(i - 1)e i + j) -
i = V
k-
j=i
-£ln(1+i+K,+(i - 1)e i + j)
i=0
Для того чтобы найти оптимальное значение Р, параметра Р i рассмотрим уравнение.
5lnL(1)(K|P ,) / dp , = 0.
Подставляя в него выражение для SlnL(1)(K|P ,) получаем
( k-1 ,
£- 1
£
K, + (i - 1)в i + j
-£ — 1
i=0 1 + i
j= k, - 1
(i -1) = 0.
V ,=0 1 +i + K, + (i - 1)e, + j’y Используя численные методы, получим оптимальное значение Р,, на основании которого сможем найти функции распределения Fn+, (k), Fn+, (k) , а значит и вероятность безотказной работы ПО после n-го отказа, которая определяется как интервал
R(k) = 1 - Fn+, (k), R(k) = 1 - Fn+,(k), где _
Beta(s + n +1 + Kn+, + nP,, k)
F' ■ )n+ (k) = 1 -
Fn“(k) = 1 -
Beta(s + Kn+, + nP,, k) Beta(s + n +1 + Kn+, + nP,, k)
Beta( Kn+, + nP,, k )
Сравнительный анализ, построенной обобщенной байесовской модели надежности По, с дискретным аналогом модели джелински-Моранда
Для сравнения моделей используем алгоритм усечения, который заключается в следующем. Рассмотрим некоторый набор статис-
i=1
ЛЕСНОЙ ВЕСТНИК 2/2009
161
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ
тических данных об отказах ПО. Из множества этих данных извлечем первые k - статистик, используя которые спрогнозируем время безотказной работы на (k + 1) шаге с помощью модели Джелински-Моранда и полученной обобщенной байесовской модели. Однако прогнозируя значения с использованием моделей надежности, мы не можем определить их точно, так как имеем дело со случайными величинами, поэтому любая модель дает только вероятностные характеристики прогнозируемого значения. Наиболее информативной с точки зрения оценки достоверности модели является ожидаемое значение (математическое ожидание) числа запусков до отказа EX Так как рассматриваемая новая модель надежности является обобщенной, то для сравнения будем использовать математическое ожидание числа запусков до отказа, которое будем вычислять, используя формулы
E'+X +1 =1 (1 - F k+i (m))
m=1
__( ) го
и E "+,Xt „=£ (1 - F “(m)).
m=
Заметим что в отличие от стандартной модели, для которой искомое математическое ожидание - это число, для байесовских интервальных моделей это интервал, определяемый нижним и верхним значениями математических ожиданий. Однако сравнивать интервал с точным значением (числом) нельзя, поэтому мы заменяем интервал точным аналогом, который вычисляем, используя формулу
-(')
E+! хм-У e(+) x, + + (i -y)E+i хм
Здесь у - коэффициент пессимизма в принятии решений о точном значении интервала. Если у = 0, то пессимизм минимальный и в качестве точной оценки имеем верхнюю оценку математического ожидания, если Y = 1, имеем максимальную степень пессимизма и нижнюю оценку.
Таким образом, выбрав у, дальше по k имеющимся значениям прогнозируем (k + 1), отметим, что эта точка на самом деле уже имеется в таблице данных. Следовательно, теперь возможно определить отклонение прогнозов для разных моделей от реальных точек. Проделав эту процедуру для (k + 1)-й, (k + 2)-й и так далее до (n - 1)-й точки, получаем статис-
тические данные отклонений прогнозируемых параметров данных моделей от реальных данных. Кроме того появляется возможность проследить тенденцию изменения числа запусков до отказа, а также сравнить усредненные значения отклонений для двух моделей. Средние значения отклонений R), R2, R3, R4 будем считать, используя формулы
Z E+) x,+i — к
R -
=k
П - k
z
=k
, R2 - -
E (+) X и - k
+ +
n — k
n—1
n-i
Z e“ x,+,—k.
Zl E,+i XM— k+
R3 =
i-k
n — k
R4 =
-k
n — k
Здесь R , R2 - соответственно нижние и верхние средние значения отклонений для обобщенной байесовской модели, R3, R4 - отклонения прогнозируемых средних значений времен до отказа от реальных данных для обобщенной байесовской модели и модели Джелински-Моранда.
Таблица
Данные о надежности программного обеспечения
Номер Число запусков до отказа
1 9
2 12
3 11
4 4
5 7
6 2
7 5
8 8
9 5
10 7
11 1
12 6
1 3 1
1 4 9
1 5 4
1 6 1
17 3
1 8 3
19 6
20 1
21 11
22 33
23 7
24 91
25 2
26 1
162
ЛЕСНОЙ ВЕСТНИК 2/2009
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ
Рис. 1. Изменение надежности ПО в процессе отладки для различных моделей
Номер этапа отладки i
Рис. 2. Изменение надежности ПО в процессе отладки для различных моделей (до 21-го этапа отладки)
Для сравнительного анализа моделей будут использованы данные о надежности программного обеспечения, которые содержатся в работе Джелински-Моранда [4], они размещены в таблице. Это информация о 26 отказах, происшедших на стадии разработки программного обеспечения для морского флота США.
Рассмотрим параметр осторожности ^ = 1, коэффициент пессимизма у = 0,5 и первые три статистики, используя которые начнем процесс прогнозирования.
Представим полученную информацию в виде графиков. По оси абсцисс откладываем номер испытания, а по оси ординат - среднее число запусков до отказа (£(1)Х. - для обобщенной байесовской модели и EX. для модели Джелински-Моранда).
График, изображенный тонкой линией с треугольными маркерами на рис. 1 и рис. 2 - это график изменения надежности ПО в процессе отладки, построенный по мо-
дели Джелински-Моранда, толстой линией - по обобщенной байесовской модели, линия средней толщины соответствует реальным статистическим данным.
Сравнивая графики, видим, что модель Джелински-Моранда ведет себя менее точно, так как она дает более завышенные значения средних значений числа запусков до отказа. Это происходит вследствие существенной «памяти» модели Джелински-Моранда. Посмотрев на первые три значения числа запусков ПО до отказа, видим, что они демонстрируют высокий уровень показателей надежности, на начальном этапе отладки модель Джелински-Моранда «запоминает» их, что, в свою очередь, приводит к завышению прогнозируемых значений средних по сравнению с реальными данными. В то время как рассматриваемая обобщенная байесовская модель является более гибкой, лучше отслеживает изменения показателей надежности и позволяет давать адекватные прогнозы, особенно в том случае, когда объем статистической информации мал.
Вычисляя средние значения отклонений R1, R2, R3, R4, получаем R1 = 8,448, R2 = 8,605, R3 = 8,485, R4 = 10,272.
Найдем средние величины относительных отклонений, используя формулы
R/ = R1 / к , R' = R2 / к ,
1 1 ср 2 2 ср’
R’ = R3 / к , R.' = R4 / к ,
3 3 ср 4 4 ср
где
к
26
Z к
i=1
26
Имеем R/ = 0,878, R/ = 0,895, R3’ = 0,882, R/ = 1,068.
Таким образом, средние значения отклонений и средние величины относительных отклонений меньше для обобщенной байесовской модели надежности по сравнению с моделью Джелински-Моранда. Эти факты указывают на явное преимущество обобщенной байесовской модели.
Вычислим средние значения отклонений по первым 13 прогнозам, т.е. для случая, когда объем статистической информации мал.
Получаем R3 = 3,280, R4 = 3,498 и R3' = 0,341, R; = 0,363.
ЛЕСНОЙ ВЕСТНИК 2/2009
163