Загребаев Методы обработки статистической информации в задачах контроля 2008

ции, т.е. как среднее арифметическое функций потерь по N измерениям:

~		~		1		N	~	~
J (c ) = M{F(ε(i, c ))}						∑F(ε(i, c ))		= J N (c ) .	(3.2.16)
					N
						i=1
							ˆ
Оптимальной	выборочной			оценкой				является	решение
							c(N )
уравнения:		~
	∂J
		(c )

=0 . (3.2.17)

∂~=ˆ

c c (N )

Вдальнейшем, если этого не требует решение задачи, будем опускать индекс N при написании выборочной оценки.N~c

3.2.3.Инвариантность оптимального решения относительно четных функций потерь (для линейного объекта)

Естественно принять допущение, что оптимальным решением уравнения (3.2.15а) является истинное значение параметра c . Тогда условие минимума (3.2.15а) средних потерь (3.2.13), опреде-

				ˆ
ляющее оптимальное решение c = c , можно записать в виде равен-
ства нулю градиента средних потерь:
~				~
J (c )		~		= M{ F[ε(i, c )]}				~	=
		c =c						c =c
		′						= 0	.	(3.2.18)
= M{F [ε(i,c )] ε(i,c )}
						c	=c

и положительной определенности матрицы Гессе – матрицы вторых производных:

	2		2			> 0 .	(3.2.19)
		J (c ) = M{		F[ε(i,c )]}
					c	=c

~

невязка ε(i, c)

между выходом объек-

Учитывая, что при c = c

та y(i) и выходом оптимальной настраиваемой модели

~

y(i) в i-й

момент времени равна возмущению η(i) , можно записать

= M

′

}= 0 .

(3.2.20)

J (c )

{F [η(i)] ε(i,c )

c =c

что η(i)

~

c

=c

статистически неза-

Обратим внимание,

~

и ε(i, c )

c

=

c

висимы.

181

Тогда

J (c )			= M{F′ [η(i)] ε(i,c )} =
		c =c
= M{F′ [η(i)]}M{ ε(i,c )} = 0 .				(3.2.21)
Тривиальное решение уравнения (3.2.21) будет при условии
			M{F′ [η(i)]} = 0 ,
или
		∞	F′ [η]p(η)dη = 0 .
		∫		(3.2.22)
		−∞

Так как p(η) , как правило, четная функция η, то для выполне-

ния тождества (3.2.22) необходимо и достаточно, чтобы		′	бы-
	F (η)
ла нечетной функцией η или F(η) – четной функцией η.
Таким образом, условие оптимальности (3.2.21) при	~	= c	вы-
	c

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

Однако на практике, как уже отмечалось, вместо математических ожиданий используются эмпирические средние, вычисляемые по наблюдаемым данным. В результате можно находить лишь те или

иные оценки cˆ(N) истинного параметра c . Причем даже если по-

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

изложенного возникает задача нахождения такой функции потерь
~	ˆ	обладает максимально воз-
F * (ε(c )) , при которой оценка	c(N )

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

настоящего пособия.

Рассмотрим наиболее часто используемые в настоящее время в теории оценивания функции потерь и штрафа:

1) квадратичные:

182

т

(3.2.23)

ϕ(c − c ) = (c − c )

R(c − c ) ;

2

;

(3.2.24)

F(ε(i,c )) = r(i)( y(i) − y

(i,c ))

где R – положительно-определенная симметричная матрица весо-

вых коэффициентов;

2) модульные:

ϕ(c − c ) =

c − c

;

(3.2.25)

F(ε(i,c )) = r(i) y(i) − y(i,c ) ;

3) простые:

i = 0, 1, 2, … ;

< δ / 2;

0,

c − c

ϕ(с − с) =

c − c

≥ δ / 2;

1 / δ,

y(i) − y(i)

< δ / 2;

0,

F(ε(i,c ))

=

y(i) − y(i)

≥ δ / 2.

1 / δ,

(3.2.26)

(3.2.27)

(3.2.28)

Можно предложить и другие функции потерь и штрафа.

3.3. Статистические методы обработки данных при использовании полного объема информации

3.3.1. Метод наименьших квадратов

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

Функция потерь, критерий качества оценки

Как уже отмечалось, методу наименьших квадратов соответствует квадратичная функция потерь. Для объекта с одним выходом эта функция имеет вид:

			2		.	(3.3.1)
F(ε(i,c )) = r(i)( y(i) − ψ(i,c ))
Критерий качества имеет вид:
		2		}.
J (c ) = M{r(i)( y(i) − ψ(i,c ))

183

Для эргодических процессов в качестве критерия качества используются средние потери:

J (c ) = 1 ∑N r(i)( y(i) − ψ(i,c ))2 .

N i=1

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

	количество используемых для оценки значений выходов объекта и
	модели и независящее от оценки			~				~	можно пе-
				c , то функцию				J (c )
	реписать в виде:			2
			N
			J (c ) = ∑r(i)( y(i) − ψ(i,c ))				,		(3.3.2)
			i=1
	или в векторной форме:
	Оценка	ˆ	J (c ) = ( y − ψ(c))т R( y −		ψ	(c)) .			(3.3.2а)
			обеспечивает минимальное				значение		критерия
		cis

(3.3.2), т.е.

ˆ	N				2
cls = arg min ∑( y(i) − ψ(i,c ))						r(i)
c	i=1
и является решением уравнения:
	∂J (c )			= 0 .
	∂c
				ˆ
			c	=cls

(3.3.3)

(3.3.4)

В общем случае, когда правая часть модели представляет собой

нелинейную функцию ψ ~ , задача минимизации критерия

(u(i), c )

(3.3.2) как правило, не может быть решена в явном виде и требует использования методов нелинейного программирования. В случае линейной модели вида (3.1.18) получены явные формулы оценки.

3.3.2. Метод наименьших квадратов для линейных регрессионных объектов

Как уже отмечалось в разд. 3.1.1, регрессионный объект и соответствующая ему модель имеют вид:

y(i) = b1u1(i) + b2u2 (i) +... + bmum (i) + η(i) ;							(3.3.5)
y(i) = b u (i) + b		u	2	(i) +... + b u	m	(i) .	(3.3.6)
1 1	2			m

184

Запишем уравнения объекта и модели в векторно-матричной форме:

				+ η ;				(3.3.7)
	y =Ub
	~				~
		=Ub					,	(3.3.8)
	y

~

где y , y – N-мерные вектора выходов объекта и модели, соответственно, U – матрица входов размерности N × m , η – N-мерный вектор шума измерений.

Представим критерий качества (3.3.2) также в векторной форме:
	~	~	~	~	~		(3.3.9)
J (b ) = ( y − y(b ))т R( y − y(b )) ,
R – диагональная матрица весовых коэффициентов.

Подставляя уравнение модели (3.3.8) в выражение							~	(3.3.9), полу-
чим явную зависимость критерия качества от оценки							b	:
	~		~		~
J (b ) = ( y −Ub )т R( y −Ub ) .									(3.3.10)

ˆ		является корнем уравнения:
Как уже отмечалось, bis
	∂J (			b		)			= 0 .
		∂	b
									ˆ
								b	=bls

Преобразуем выражение (3.3.10) к виду, удобному для диффе-

ренцирования. Для этого произведем указанные действия
	~		~			~			~
J (b ) = y тRy + 2b					тU тRy + b			тU тRUb .			(3.3.11)

Дифференцируя (3.3.11) и приравнивая полученное выражение нулю, запишем:

−U

т

Ry +U

т

ˆ

= 0

RUbLS

или

ˆ

U

т

т

Ry .

(3.3.12)

RUbls =U

ˆ

Разрешая (3.3.12) относительно bls , получаем искомую оценку

наименьших квадратов

ˆ

т

RU )

−1 т

Ry .

(3.3.13)

bls = (U

U

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

185

рованная форма удобна и в том случае, если в уравнениях (3.3.5) и (3.3.6) присутствует свободный член, т.е.

y(i) = b0 +b1u1 (i) +... +bmum (i) + η(i) , i =				;	(3.3.14)
			1, N
~	~	~			(3.3.15)
y(i) = b0		+ b1u1 (i) + ... + bmum (i) .

Получим центрированную форму на примере использования объекта (3.3.5). Для объектов типа (3.3.14) решение задачи приведено в конце данного раздела. Просуммируем поэлементно урав-

нения (3.3.14) по всем i =1, N :

N	N	N	N
∑y(i) = b1	∑u1 (i) +... +bm ∑um (i) + ∑η(i) .
i=1	i=1	i=1	i=1

Разделив последнее выражение на N, получим уравнение объекта, записанное относительно средних значений входов и выхода:

yср = b uср +... +b uср + ηср .		(3.3.16)
1 1	m m

Вычитая из (3.3.5) (3.3.16), получим уравнение объекта относительно центрированных входов и выходов:

yо(i) = b uо(i) +... +b			u	о	(i) + ηо(i) ;
1	1	m	m
yо(i) = y(i) − yср(i) , uоj	(i) = u j (i) − u		срj	,	j =		, i =		. (3.3.17)
						1, m		1, N

Модель, соответствующая (3.3.17), будет иметь вид:

~о

~

о

(i) +

~

о

(3.3.18)

y

(i) = b1u1

... + bmum (i) .

ˆ

можно записать как

Тогда оценку bls

ˆ

от

RU

о

)

−1 от

Ry ,

(3.3.19)

где

bls = (U

U

U о =U −U ср ;

y о = y − y ср.

Очевидно, уравнения (3.3.13) и (3.3.19) совершенно эквивалентны. Условием единственности решения этих уравнений является

условие невырожденности матрицы (U тRU ) , или, что эквивалент-

но, (U отRU о) , т.е.:

rank(U тRU ) = rank(U отRU о) = m ,

где m – число оцениваемых параметров.

186

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

измерений N меньше m, то матрица (U тRU ) оказывается вырож-

денной. В том случае, если N = m – оценивание по достаточному числу измерений. В этом случае уравнение (3.3.13) (или (3.3.19)) может быть преобразовано к тривиальной форме:

ˆ

= U

−1

y

(3.3.20а)

или

bls

ˆ

о−1

y .

(3.3.20б)

bls =U

Задача. Найти оценку наименьших квадратов для линейного

объекта

y(i) = b0 +b1u1 (i) +... +bmum (i) + η(i) .

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

ет вид:

yо (i) = b~uо (i) +... + b~ uо

(i) .

1 1

m m

~

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

. Оценку параметров

~

~

, … ,

~

найдем по формуле

b0

b1 ,

b2

bm

(3.3.19):

ˆ

(U

от

RU

о

)

−1

от

Ry

о

,

bls′ =

U

ˆ	ˆ
где bls′	– вектор оценок без b0 .

Для нахождения оценки ˆ b0

ˆ	воспользу-
и ошибки оценки b0 − b0

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

yср = b0 + b1u1ср +... +bmumср + ηср ;

~ср	~	~ ср	~ ср
y	= b0	+ b1u1	+... + bmum .

Тогда оценка и ошибка оценки параметра b0 могут быть определены по формулам:

ˆ	= y	ср	ˆ ср	ˆ ср	,
b0			− b1u1	−... −bmum

187

ˆ	ˆ	ср	ˆ	ср	+ η	ср	.
b0 −b0	= (b1 −b1 )u1		+... + (bm −bm )um
Как видно из последней формулы, ошибка оценки b0 включает

в себя ошибки оценок других параметров, а также среднестатистическое значение шума измерений. Как правило, оценка параметра

b0 оказывается наиболее неточной.

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

Ковариационная матрица ошибки оценки

Найдем ковариационную матрицу ошибок оценки, при этом бу-

ˆ

ˆ

– несмещенная, т.е.

= M{b}.

дем считать, что оценка bLS

M{bLS }

ˆ

ˆ

ˆ

т

}.

var{bLS − b} = M {(bLS − b )(bLS − b )

Воспользуемся формулой для оценки (3.3.13), тогда

ˆ

т

+Qη −b )(QUb +Qη −b )

},

(3.3.21)

var{bLS −b} = M{(QUb

где Q = (U тRU )−1U тR .

Для конкретных N реализаций матрица входа U является детер-

минированной величиной,

кроме того, учитывая, что QU = I , по-

лучим

ˆ

т

т

− b ) = QM{ηη

}Q

.

(3.3.22)

var(bLS

Но M{ηηт} – ковариационная матрица помехи η .

Таким образом, имеем

ˆ

т

RU )

−1

U

т

RDηRU

т

(U

т

RU )

−1

. (3.3.23)

var(bLS − b ) = (U

Выражение (3.3.23) определяет ковариационную матрицу ошибок несмещенных оценок наименьших квадратов при линейном уравнении объекта.

Если измерения некоррелированы и равноточны (Dη = ση2 I ) , то

в качестве матрицы весов R можно принять единичную матрицу,

т.е. R = I .

188

В этом случае выражение (3.3.23) существенно упрощается:

ˆ	− b) = (U	т	U )	−1	2	(3.3.24)
var(bLS					ση .

Нетрудно заметить, что если в формулу для ковариационной матрицы (3.3.21) вместо выражения (3.3.13) подставить эквивалентное ему выражение (3.3.20), то получим выражение для ковариационной матрицы ошибки оценки, записанное относительно центрированных входов и выходов:

ˆ

от

RU

о −1

U

от

RDηRU

о

от

RU

о

var{bLS − b} = U

U

. (3.3.25)

Или для некоррелированных равноточных измерений и R = I :

ˆ

от

U

о −1

2

(3.3.26)

var{bLS

− b} = U

ση .

Очевидно, выражения (3.3.23) и (3.3.25), а также (3.3.24) и (3.3.26) совершенно эквивалентны. В дальнейшем, в зависимости от решаемой задачи, мы будем использовать то или иное выражение для ковариационной матрицы ошибки оценки.

Свойства оценки наименьших квадратов

вслучае линейных регрессионных объектов

1.Линейность. Очевидно, оценка (3.3.13) (или эквивалентная ей оценка 1.3.19) линейны.

2.Несмещенность. Найдем математическое ожидание оценки. Для этого воспользуемся формулой (3.3.13):

ˆ

т

RU )

−1

U

т

Ry}.

M{bLS } = M{(U

Подставим вместо y его выражение (3.3.7):

ˆ

т

−1

т

RU )

R(Ub + η)} =

M{bLS } = M{(U

U

= M{

}+ (U тRU )−1U тRM{η} .

(3.3.27)

b

Очевидно, что оценка

ˆ

– несмещенная, если

M{η} = 0 . В

bLS

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

= (U тRU )−1U тRM{η} .

(3.3.28)

δ

189

3. Состоятельность. Покажем, что оценка наименьших квадратов состоятельна, Воспользуемся формулой (3.3.26) для ковариационной матрицы ошибки оценки при некоррелированных равноточных измерениях. Констатируем факт, что оценка получена по N измерениям:

	ˆ			от	U	о −1	2
var{bLS (N) − b} = U							ση .

Возьмем предел при N → ∞ от правой и левой частей последнего выражения:

					от			−1
	ˆ					U	о			2
lim var{bLS (N ) − b} =				lim U					ση .
N →∞				N →∞

т

U

о

Представим матрицу U о

в поэлементной форме:

т

=

U о

U о

N

N

)(u2

(i) −μu )

...

∑(u1(i) −μu )2

∑(u1(i) −μu

i=1

1

i=1

1

2

=

...

...

...

.

N

N

∑(um(i) −μu

)(u1(i)

−μu )

...

∑(um (i)

−μu

)2

m

1

i=1

m

i=1

Нетрудно заметить, что эта матрица представляет собой оценку ковариационной матрицы входов, умноженную на число измерений N, т.е.:

(U отU о ) Du N .

(3.3.29)

Учитывая, что элементы матрицы Du ограничены, получим:

lim	(U отU о )−1 ση2 = lim	1	Duση2 = 0 .
N →∞	N →∞ N

Последнее соотношение эквивалентно условию состоятельности оценки. Используя приближенное равенство (3.3.29), можно записать оценочную формулу для ковариационной матрицы ошибки оценки:

190