close

Вход

Забыли?

вход по аккаунту

?

Исследование максиминно эффективных планов для модели Михаэлиса—Ментена.

код для вставкиСкачать
УДК 519.24
Вестник СПбГУ. Сер. 1, 2006, вып. 2
В. Б. Мелас, Ю. М. Старосельский
ИССЛЕДОВАНИЕ МАКСИМИННО ЭФФЕКТИВНЫХ
ПЛАНОВ ДЛЯ МОДЕЛИ МИХАЭЛИСА—МЕНТЕНА∗
1. Введение. Регрессионная модель Михаэлиса—Ментена
E[Y |t] =
at
,
t+b
t ∈ [0, t0 ],
(1)
где a и b — оцениваемые параметры, широко используется для описания физических и
биологических феноменов (см. например, [1] и [2]). Задача построения планов эксперимента, локально оптимальных в смысле различных критериев, изучалась в работах
[3], [4] и многих других. В настоящей работе рассматриваются стандартизированные
максминно эффективные E-оптимальные планы (далее для краткости такие планы
называются СММЭ-планами). Такие планы были численно построены в работе [5] на
основе аналитического решения промежуточной задачи — построения стандартизированных локально E-оптимальных планов. Преимущество СММЭ-планов по сравнению
с локально оптимальными планами заключается в использовании не точечных, а интервальных оценок для нелинейно входящего параметра b.
В настоящей работе мы исследуем СММЭ-планы на основе функционального подхода, введенного в работах [6], [7] и развитого в монографии [8]. Этот подход позволяет
аппроксимировать опорные точки и весовые коэффициенты СММЭ-планов отрезками
рядов Тейлора.
В следующем разделе мы сформулируем задачу и представим предварительные результаты. В разделе 3 формулируются основные теоретические результаты, а в разделе
4 приводятся результаты численного исследования. Доказательство основных теоретических результатов дано в приложении.
2. Постановка задачи. Модель (1) представляет собой специальный случай общей
нелинейной регресионной модели
E(Y |t) = η(t, θ),
где θ = (θ1 , . . . , θm )T — вектор неизвестных параметров. Примем стандартное предположение о том, что ошибки измерений — независимые одинаково распределенные случайные величины. Под планом эксперимента будем понимать дискретную вероятностную
меру
x1 . . . xn
,
ξ=
w1 . . . wn
где ti = tj (i = j) — точки, в которых проводятся наблюдения,а wi — весовые коэффициn
енты: wi > 0,
wi = 1. Пусть N — общее количество экспериментов. Будем говорить,
i=1
что эксперимент проводится в соответствии с планом ξ, если в точках ti проведено ri
∗ Работа
c
выполнена при частичной финансовой поддержке РФФИ (грант № 06-01-00284).
В. Б. Мелас, Ю. М. Старосельский, 2006
41
экспериментов, где ri = wi N или ri = wi N + 1 так, что
n
ri = N. Пусть оценивание
i=1
параметров производится по методу наименьших квадратов,
θ̂(N ) = arg min
ri
N θ
(Yij − η(ti , θ))2 ,
i=1 j=1
где Yij — результат j–го измерения в точке ti , i = 1, . . . , n; j = 1, . . . , ri .
Обозначим через θ∗ истинное значение параметра θ. Хорошо известно (см., например, [9]), что при некоторых условиях
√ регулярности (эти условия нетрудно проверить
для модели (1)) при N → ∞ вектор N (θ̂N −θ∗ ) имеет асимптотически нормальное
распределение N (0, σ 2 D), где σ 2 — дисперсия ошибки измерения, а D = M −1 (ξ, θ) θ=θ∗ ,
где
m
∂η(t, θ) ∂η(t, θ)
ξ(dt)
M (ξ, θ) =
∂θi
∂θj
i,j=1
— информационная матрица Фишера.
Основная трудность планирования эксперимента для нелинейных по параметрам
моделей состоит в зависимости информационной матрицы от значений неизвестных
параметров.
Рассмотрим регрессионную модель (1). Информационная матрица Фишера для этой
модели имеет вид
M (ξ, a, b) = Da f (t, b)f (t, b)T ξ(dt) Da ,
где
(
Da = diag{1, a},
f (t, b) =
t
t
t+b , − (t+b)2
)T
.
Следуя работе [10], введем стандартизированный критерий E-оптимальности. Через A−
будем обозначать матрицу, обобщенно обратную для матрицы A. Положим dT A− d = ∞,
если вектор d не принадлежит подпространству, натянутому на столбцы матрицы A.
Обозначим eT1 = (1, 0), eT2 = (0, 1). Через ξj∗ (a, b), j = 1, 2, обозначим локально
ej -оптимальный план эксперимента, т. е. план, минимизирующий величину
eTj M − (ξ, a, b)ej
при фиксированных значениях a и b, являющихся начальными приближениями для
истинных значений a∗ , b∗ . Положим
+
*
1
1
Ka,b = diag (eT1 M − (ξ1∗ , a, b)e1 )− 2 , (eT2 M − (ξ2∗ , a, b)e2 )− 2 .
Заметим, что матрица
T
C(ξ, b) = (Ka,b
M − (ξ, a, b)Ka,b )−1
не зависит от a. План, максимизирующий величину λmin (C(ξ, b)), будем называть стандартизированным локально E-оптимальным планом.
Для преодоления зависимости критерия от начального приближения для параметра
b, следуя работе [5], введем следующее определение.
42
Определение 1. План, максимизирующий величину
min
b∈[b1 ,b2 ]
λmin (C(ξ, b))
,
max λmin (C(η, b))
η
где [b1 , b2 ] — интервал возможных значений параметра b, будем называть стандартизированным максиминно эффективным E-оптимальным планом (или, для краткости, СММЭ-планом).
В работе [5] стандартизированные локально E-оптимальные планы для модели (1)
найдены в явном виде. Кроме того, доказано, что при достаточно малой величине b2 −b1
СММЭ-план имеет вид
t̂
t0
(2)
ξ̂ =
ŵ 1 − ŵ
и максимизирует величину
min αλmin (C(ξ, b1 )) + (1 − α)λmin (C(ξ, b2 )).
α∈[0,1]
В следующем разделе мы исследуем СММЭ-планы с помощью функционального подхода, описанного в [8].
3. Аналитические свойства СММЭ-планов. Фиксируем одну из границ b1 и
будем изучать зависимость величин t̂ и ŵ, входящих в формулу (2), от границы z = b2 .
Аналогичные результаты получаются, если зафиксировать границу b2 и z = b1 . Пусть
Φ(α, ξ, z) = αλmin (C(ξ, b1 )) + (1 − α)λmin (C(ξ, z)).
t
t0
,
Обозначим u = (u1 , u2 , u3 ) = (t, w, α), ξu =
w 1−w
Ψ(u, z) = Φ(α, ξu , z).
Рассмотрим систему уравнений
∂
Ψ(u, z) = 0, i = 1, 2, 3.
∂ui
(3)
Пусть ũ(z) = (t̃(z), w̃(z), α̃(z)) — некоторое решение этой системы. Через J(z) обозначим
матрицу Якоби системы (3), оцененную для этого решения:
2
3
∂ Ψ(u, z)
J(z) =
.
(4)
∂ui ∂uj i,j=1 u=ũ(z)
Справедлива следующая теорема
Теорема 1. Если 0 < z − b1 < δ, где δ — достаточно малая положительная
величина, то система (3) имеет единственное решение такое, что 0 ≤ u1 < t0 ,
0 < ui < 1, i = 1, 2, причем матрица Якоби (4) обратима, а функции t̃, w̃, и α̃, сопоставляющие z единственное решение системы (3), являются вещественными аналитическими функциями. Кроме того, план
43
t̃
t0
ξ = ξũ =
,
w̃ 1 − w̃
где ũ = ũ(z) = (t̃(z), w̃(z), α̃(z)), есть СММЭ-план для модели (1) на отрезке [b1 , b2 ] =
[b1 , z].
Доказательство этой теоремы будет дано в разделе 5.
Хотя теорема утверждает только факты для малых δ, численные расчеты показывают, что при b1 = 1 и δ ≤ 7 матрица J(z) обратима. Отсюда вытекает, что все утверждения теоремы при b1 = 1 справедливы для δ ≤ 7. Далее для упрощения обозначений
знак «волна» у функции ũ(z) будем опускать.
Пусть u(z) = ũ(z), u0 = ũ(z0 ), и для любой вещественно аналитической функции
ϕ(z)
1 ds
ϕ(z) , s ≥ 1,
ϕ(s) =
s! dz s
z=z0
где z0 > b1 некоторая фиксированная точка.
Разложение функции u(z) в ряд Тейлора в окрестности этой точки в наших обозначениях принимает вид
∞
u(z) =
u(i) (z − z0 )i .
i=0
Коэффициент u(0) может быть найден численно, путем решения системы (3) при z = z0 .
Остальные коэффициенты могут быть найдены численно с помощью общих рекуррентных формул, введенных в работе [11] (см. также [8, гл. 2]).
Для системы (3) эти формулы принимают вид
u(s) = J0−1 (z0 )(g(u<s−1> (z), z))(s) ,
(5)
где
g(u, z) =
∂Ψ(u, z) ∂Ψ(u, z) ∂Ψ(u, z)
,
,
∂u1
∂u2
∂u3
u<k−1> (z) =
k−1
T
,
u(t) (z − z0 )t .
t=0
Формула (5) и теорема 1 позволяют предложить следующий алгоритм построения и
исследования СММЭ-планов.
10 Численно решаем систему уравнений (3) для некоторого значения z = z0 > b1 .
Полагаем k0 = 0.
20 С помощью рекуррентных формул (5) находим u(k) при k = k0 + 1.
30 Проверяем для плана ξ = ξu , u = u<k> (z) выполнение с заданной точностью условий теоремы 3.3 из работы [5] при интересующих нас значениях z. Если условия
не выполняются, то полагаем k0 := k0 + 1 и возвращаемся к шагу 20 .
В следующем разделе мы продемонстрируем работу этого алгоритма.
4. Численные результаты. Рассмотрим случай b1 = 1, z = b2 , t0 = 10, z0 = 5
(для других случаев мы получили аналогичные результаты, которые опускаются ради
краткости). С помощью пакета Matlab находим
u(0) = (1.1757, 0.4550, 0.5575).
44
После трех итераций алгоритма получаем коэффициенты u(j) , j = 1, 2, 3. Значения этих
коэффициентов приведены в таблице 1. Планы, построенные с помощью этих коэффициентов для z = 1, 2, 3, 4, 6, 7, 8, приведены в таблице 2. Все эти планы удовлетворяют
необходимым и достаточным условиям, накладываемым на СММЭ-планы с точностью
до 10−5 . Отметим, что они также совпадают с планами, численно найденными в работе
[5]. Заметим, что при z ≥ 8 СММЭ-планы уже не имеют вид (2), а число их опорных
точек равно 3. Однако планы, построенные с помощью описанного алгоритма, и при
z ≥ 8 имеют минимальную эффективность
λmin (C(ξ, b))
,
b∈[b1 ,b2 ] max λmin (C(η, b))
eff = min
η
близкую к минимальной эффективности СММЭ-планов, найденных численно в работе
[5], см. табл. 3.
Таблица 1. Коэффициены u(0) , u(1) , u(2) , u(3)
∞
P
разложения u∗ (z) =
u(s) (z − z0 )s при z0 = 5
s=0
u(0)
1.1757
0.4550
0.5575
u(1)
0.0816
−0.0064
0.0015
u(2)
−0.0084
0.0008
−0.0005
u(3)
0.0009
0.0001
0.0001
Таблица 2. Опорные точки и веса СММЭ-планов
для различных z
z
x1
x2
w1
w2
1
0.6045
10
0.4844
0.5156
2
0.8169
10
0.5120
0.4880
3
0.9693
10
0.5274
0.4726
4
1.0847
10
0.5375
0.4625
6
1.2497
10
0.5506
0.4494
7
1.3111
10
0.552
0.448
8
1.363
10
0.5587
0.4413
Заметим, что на практике обычно используются планы в равноотстоящих точках.
Численное исследование показывает, что минимум эффективности таких планов достигается на левой границе при b = b1 . Эта минимальная эффективность растет с
увеличением числа опорных точек, но довольно медленно. Для плана в 15-и равноотстоящих точках она равна 0.3685. Таким образом, применение СММЭ-планов позволяет
уменьшить число экспериментов в 2 и более раза в зависимости от величины b2 .
Таблица 3. Минимальные эффективности СММЭ-планов
найденных численно (ξ ∗ ) и планов полученных
с помощью рядов Тейлора (ξ2 )
z
eff(ξ ∗ )
eff(ξ2 )
2
0.9544
0.9544
4
0.8451
0.8451
6
0.7733
0.7733
8
0.7270
0.7253
10
0.7034
0.6911
12
0.6906
0.6657
14
0.6829
0.6460
16
0.6779
0.6302
Отметим, что используя коэффициенты, указанные в табл. 2, читатель может построить вручную план для интервала вида [b1 , b2 ] = [1, z] при t0 = 10. Первая опорная
точка и ее вес вычисляются по формулам
x1 = t̂
ω1 = ω̂
≈ 1.1757 + 0.0816(z − 5) − 0.0084(z − 5)2 + 0.0009(z − 5)3 ,
≈ 0.4550 − 0.0064(z − 5) + 0.0008(z − 5)2 + 0.0001(z − 5)3 .
Вторая опорная точка совпадает с t0 = 10, а ее вес равен 1 − ω1 .
45
5. Приложение: доказательство теоремы 1. Решающий пункт доказательства
заключается в проверке невырожденности матрицы Якоби.
Пусть сначала z — произвольное
число, удовлетворяющее
неравенству z > b1 .
t
t0
,
Обозначим τ = (t, w), ξτ =
w 1−w
H(τ, b) = λmin (C(ξτ , b)).
Непосредственное вычисление показывает, что матрица J = J(z) имеет следующую
структуру:
A d
J=
,
dT 0
где
A=
T
d =
∂ 2 Ψ(u, z)
∂ui ∂uj
2
,
i,j=1
∂
∂
[H(τ, b1 ) − H(τ, z)],
[H(τ, b1 ) − H(τ, z)] ,
∂t
∂w
u = u(z) = (u1 , u2 , u3 ) — некоторое решение системы (3), τ = (u1 , u2 ).
В силу известных формул Фробениуса (см., например, [9]), если det A = 0, d =
(0, 0)T , то
det J = −(det A)dT A−1 d = 0.
Покажем, что d = (0, 0)T . Пусть, напротив, d = (0, 0)T . Тогда, как показывает прямое
вычисление, производные по t и w функций
u3 H(τ, b1 ) + (1 − u3 )H(τ, z) и H(τ, b1 ) − H(τ, z)
равны нулю при τ = (u1 , u2 ). Отсюда вытекает, что
Ht (τ, b) = 0, Hw (τ, b) = 0
(6)
при b = b1 и b = z. В работе [5] доказано, что единственным решением уравнений (6)
по τ является вектор τ (b) = (t(b), w(b)) такой, что план
t(b)
w(b)
t0
1 − w(b)
есть стандартизированный локально E-оптимальный план. Из явной формулы [5, формула (3.4)] вытекает, что t(b1 ) = t(z). Таким образом, выполнение равенства (6) для
одного и того же вектора τ невозможно. Полученное противоречие доказывает, что
b = (0, 0)T .
Пусть теперь 0 < z −b1 < δ, где δ достаточно мало. В силу аргумента непрерывности
при δ → 0 матрица A стремится к матрице Якоби для системы уравнений
∂
∂
H(τ, b1 ) = 0,
H(τ, b1 ) = 0,
∂t
∂w
46
которые однозначно определяют единственный стандартизированный локально E-оп¯ 1 ):
тимальный план. Обозначим эту матрицу J(b
2
∂2
¯
J(b) =
H(τ, b)
,
∂τi ∂τj
i,j=1
t(b)
t0
w(b) 1 − w(b)
Обозначим
τ =(t(b),w(b))
— стандартизированный локально E-оптимальный план.
ξτ =
Q(τ, p2 , b) =
pT C(ξτ (b) , b)p
,
pT p
где p = (1, p2 )T , 0 < b < ∞.
Пусть τ ∗ (b) = (t∗ (b), w∗ (b)) таково, что ξτ ∗ (b) — стандартизированный локально
E-оптимальный план, (1, p∗2 (b))T — собственный вектор матрицы C(ξτ ∗ (b) , b), соответствующий ее минимальному собственному числу. В силу результатов [5] τ ∗ (b) и p∗2 (b)
определены однозначно.
В силу свойств минимального собственного числа имеем
Q(τ ∗ (b), p∗2 (b), b) = λmin (C(ξτ ∗ (b) , b)),
∂
Q(τ ∗ (b), p2 , b) = 0 при p2 = p∗2 (b),
∂p2
∂
Q(τ, p∗2 (b), b) = 0 при τ = τ ∗ (b), i = 1, 2.
∂τi
Обозначим
(7)
2 ,
∂2
∗
Q(τ (b), p2 , b) = T eT2 Ce2 − λ ,
d(b) =
∂p22
p
p
p2 =p∗ (b)
2
где
p = (1, p∗2 (b))T ,
C = C(ξτ ∗ (b) , b),
λmin = λmin (C),
e2 = (0, 1)T .
В силу единственности минимального собственного числа матрицы C (см. [5]), получаем, что d(b) > 0.
Пусть
( 2
)2
G(b) = ∂τ∂i ∂τj Q(τ, p2 , b)
, τ = τ ∗ (b), p2 = p∗2 (b),
)i,j=1
( 2
2
, τ = τ ∗ (b), p2 = p∗2 (b),
q(b) = ∂τ∂i ∂p2 Q(τ, p2 , b)
i=1
q2 (b) = (q(b))2 .
Нам потребуется следующая
Лемма. При любом b, 0 < b < ∞
T
¯ = G(b) − q(b)q (b) .
J(b)
d(b)
47
Доказательство этой леммы аналогично доказательству теоремы 2.4.2 в работе [8].
Непосредственное вычисление показывает, что
2
G11 (b) 0
¯ = −G11 (b) [q2 (b)] ,
G(b) =
, det J(b)
0
0
d(b)
где
−1
(
) ∂2
f (t) pT K1,b
−1
T
∗
G11 (b) = 2w (b) p K1,b f (t )
,
∗
∂t2
pT p
t=t
−1
T
∗ ∗
p
K
f
(t
)
t
t
1,b
0
−1
q2 (b) = −2(K1,b )22
+
,
pT p
(t∗ + b)2
(t0 + b)2
∗
t∗ = t∗ (b),
p = (1, p∗2 (b))T .
Из (7) следует, что
−1
[pT K1,b
f (t∗ )]2
pT p
= λmin (C(ξτ ∗ (b) , b))
и, значит, q2 (b) > 0.
Докажем, что G11 (b) = 0. Из соотношения (7) следует
g (t∗ (b)) = 0,
(8)
−1
где g(t) = pT K1,b
f (t), p = (1, p∗2 (b))T .
−1
−1
)11 , c2 = (K1,b
)22 p∗ (b).
Пусть c1 = (K1,b
Тогда имеем
c1 tb + c1 b2 − 2c2
c1 t2 + c1 bt + c2 ,
g
(t)
=
,
(t + b)2
(t + b)3
2
c1 b(t + b) − 3(c1 tb + c1 b − 2c2 )
.
g (t) =
(t + b)4
g(t) =
Из (8) получаем соотношение
c1 t∗ (b)b + c1 b2 − 2c2 = 0
и, значит,
0 < t∗ (b) + b =
2c2
,
c1 b
g (t∗ (b)) =
c1 b
∗
(t (b) +
b)3
= 0.
¯ = 0.
Отсюда вытекает, что G11 (b) = 0 и det J(b)
Следовательно, матрица J(z) обратима, а система (3) имеет единственное решение
для достаточно малых δ. Аналитичность функций t̃(z), w̃(z), α̃(z) при достаточно малых δ вытекает из теоремы о неявной функции (см., например, [12]).
48
Summary
V. B. Melas, Yu. M. Staroselsky. The maximin efficient designs for the Michaelis—Menten model.
The standardized maximin efficient E-optimal experimental designs are studied by the functional
approach. This approach allows one to approximate the support points and weight factors of designs
by the Taylor series.
Литература
1. Beverton R. J. H., Holt S. J. On the Dynamics of Exploited Fish Populations. London, 1957.
2. Johansen S. Functional relations, Random Coefficients and Nonlinear Regression, with application to Kinetic Data. New York: Springer, 1984.
3. Duggleby R. G. Experimental designs for estimating the kinetic parameters for enzymecatalysed reactions // J. Theoret. Biology, 1979. Vol. 81. P. 671–684.
4. Dunn G. Optimal designs for drug, neurotransmitter and hormone receptor assays // Statist.
Medicine, 1988. Vol. 7. P. 805–815.
5. Dette H., Melas V., Pepelyshev A. Standardized E-optimal designs for the Michaelis—Menten
model // Statistica Sinica, 2003. Vol. 13. P. 1147–1163.
6. Melas V. B. Optimal designs for exponential regression // Math. Operationsforsh. Statist.
1978. Vol. 9. P. 45–59.
7. Мелас В. Б. Оптимальное планирование эксперимента для экспоненциальной регрессии
// Математические методы планирования эксперимента / Под ред. В. В. Пененко. Новосибирск: Наука, 1981. C. 174–198.
8. Melas V. B. Functional Approach to Optimal Experimental Design. Heidelberg: Springer,
2006. 336 p.
9. Ермаков С. М., Жиглявский А. А. Математическая теория оптимального эксперимента.
М.: Наука, 1987. 320 c.
10. Dette H. E-optimal designs for regression models with quantitive factors — a reasonable
choice // Canad. J. Statist. 1997. Vol. 25. P. 531–543.
11. Мелас В. Б., Пепелышев А. Н. Степенные разложения неявных функций и локально
оптимальные планы эксперимента // Статистические модели с приложениями в эконометрике
и смежных областях / Под ред. С. М. Ермакова и Ю. Н. Каштанова. СПб., 1999. С. 108–117.
12. Ганнинг Р., Росси Х. Аналитические функции многих комплексных переменных. М.:
Наука, 1969.
Статья поступила в редакцию 12 февраля 2006 г.
49
Документ
Категория
Без категории
Просмотров
7
Размер файла
268 Кб
Теги
планово, михаэлиса, ментен, максиминное, эффективные, исследование, модель
1/--страниц
Пожаловаться на содержимое документа