advanced.tex

\chapter{Приложение}

\small

\section{Корреляции}
\label{sec:correlation}

Напомним, если величины $x$ и $y$ независимы, то среднее значение (математическое ожидание) произведения отклонений $\Delta x = x-\limaverage{x}$ и $\Delta y = y - \limaverage{y}$
равно нулю:
\begin{equation*}
\limaverage{\Delta x\cdot\Delta y}=\limaverage{\Delta x}\cdot\limaverage{\Delta
y}=0.
\end{equation*}
Если же $x$ и $y$ не являются полностью независимыми, среднее значение произведения
их отклонений может быть использовано как количественная мера их зависимости.
Наиболее употребительной мерой зависимости двух случайных величин
является \emph{коэффициент линейной корреляции}:
\begin{equation}
r_{xy}=\frac{\average{\Delta x\cdot\Delta
y}}{\sigma_{x}\cdot\sigma_{y}}.\label{eq:pearson}
\end{equation}
Нетрудно проверить (с помощью неравенства Коши\textendash Буняковского),
что $-1\le r\le1$. В частности, для полностью независимых величин
коэффициент корреляции равен нулю, $r=0$, а для линейно зависимых
$y=kx+b$ нетрудно получить $r=1$ при $k>0$ и $r=-1$ при $k<0$.
Примеры промежуточных случаев представлены на рис. TODO.

\note{Угловой коэффициент прямой в задаче линейной регрессии \eqref{eq:MNK} выражается
через коэффициент корреляции как
\begin{equation*}
k = r_{xy}\frac{\sigma_y}{\sigma_x}.
\end{equation*}}

Если коэффициент $r_{xy}$ близок к единице, говорят, что величины
\emph{коррелируют} между собой (от \emph{англ.} correlate ---
находиться в связи).

\todo[inline,author=ppv]{А не объяснить ли заодно, что чудный $R^2$, который выдаёт эксель,
это для прямой и есть $r^2$?}

\paragraph{Отсутствие корреляции $\protect\not\Rightarrow$ независимость.}

Отметим, что (\ref{eq:indep}) --- необходимое,
но не достаточное условие независимости величин. На рис. TODO приведён
пример очевидно зависимых $x$ и $y$, для которых $r\approx0$.

\paragraph{Корреляция $\protect\not\Rightarrow$ причинность.}

Ещё одна типичная ошибка --- исходя из большого
коэффициента корреляции ($r\to1$) между двумя величинами сделать
вывод о функциональной (причинной) связи между $x$ и $y$. Рассмотрим
конкретный пример. Между током и напряжением на некотором резисторе
имеет место линейная зависимость $U=IR$, и коэффициент корреляции
$r_{UI}$ действительно равен единице. Однако \emph{обратное
в общем случае неверно}. Например, ток в резисторе коррелирует
с его температурой $T$, $r_{IT}\to1$ (больше ток --- больше
тепловыделение по закону Джоуля\textendash Ленца), однако ясно, что
нагрев резистора извне не приведёт к повышению тока в нём (скорее
наоборот, так как сопротивление металлов с температурой растёт). Ошибка
отождествления корреляции и причинности особенно характерна при исследовании
сложных многофакторных систем, например, в медицине, социологии и
т.п.

\section{Свойства точечных оценок} \label{sec:point}
Если измеряется одна физическая величина $x$, то можно поставить задачу
по конечному набору данных $\mathbf{x}=\{x_i\}$ ($i=1\ldots n$) оценить параметры
случайного распределения, которому подчиняется $x$. В частности,
найти среднее значение (математическое ожидание)~$\limaverage{x}$ и
дисперсию~$\sigma^2$.

Если результатом оценки параметра является просто число~---
без указания интервала, в котором может лежать истинное значение,~---
такую оценку называют \emph{точечной}.
Пример точечных оценок дают формулы для выборочного среднего
\eqref{eq:average}:
\begin{equation}
\limaverage{x} \approx \average{x} = \frac{1}{n}\sum_i x_i\qquad
\label{eq:estimate_x_avg}
\end{equation}
и выборочной дисперсии \eqref{eq:sigma}:
\begin{equation}
\sigma^2 \approx s^2_n = \frac{1}{n}\sum_i (x_i - \average{x})^2.
\label{eq:estimate_x_sigma}
\end{equation}

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

Оценка параметров должна давать правильное значение
хотя бы в пределе большого числа измерений.
Если при $n\to \infty$ оценка стремится к истинному значению параметра,
\[
\lim_{n\to \infty} \hat{\theta}(\mathbf{x}) \to \limaverage{\theta},
\]
то её называют \emph{состоятельной}.
Можно показать (см. \cite{idie}), что если у распределения,
которому подчиняется случайная величина,
существуют конечные средние и дисперсия, то оценки
\eqref{eq:estimate_x_avg}, \eqref{eq:estimate_x_sigma} являются состоятельными.

\paragraph{Несмещенные оценки.}
Рассмотрим случай малого числа измерений ($n\gtrsim 1$).
Тогда даже если оценка состоятельна, она может давать довольно большую ошибку.
При фиксированном $n$ функцию оценки $\hat{\theta}(x_1,\ldots,x_n)$
можно рассматривать как случайную величину с некоторым распределением,
отличающимся от распределения измеряемой величины.
Естественно потребовать, чтобы среднее (математическое ожидание) этого
распределения совпадало с истинным значением искомого параметра:
\[
\limaverage{\strut\hat{\theta}(\mathbf{x})} = \limaverage{\theta}.
\]
В таком случае оценку называют \emph{несмещённой}.

Нетрудно показать, что выборочное среднее \eqref{eq:estimate_x_avg}
является несмещённой оценкой. А вот оценка $s_n^2$ из
\eqref{eq:estimate_x_sigma} таким свойством не обладает. Математическое
ожидание для величины $s_n^2$ при фиксированном $n$ оказывается равно
$\limaverage{s_n^2} = \frac{n-1}{n} \sigma^2$ (предлагаем в качестве упражнения
проверить данное утверждение самостоятельно). Именно поэтому при малых $n$
для оценки дисперсии рекомендуется использовать формулу
\eqref{eq:sigma_straight}:
\[
\sigma^2 \approx s^2_{n-1} = \frac{1}{n-1}\sum_i (x_i - \average{x})^2.
\]

% \example{Рассмотрим формулу \eqref{eq:estimate_x_sigma} для среднеквадратичного
%     отклонения при $n=2$. Выборочное среднее равно
%     $\average{x}=\frac{x_{1}+x_{2}}{2}$, выборочная дисперсия:
%     \[
%     s_{x}^{2}=\frac{1}{2}\left[\left(x_{1}-\frac{x_{1}+x_{2}}{2}\right)^{2}+\left(x_{2}-\frac{x_{1}+x_{2}}{2}\right)^{2}\right]=\frac{1}{4}x_{1}^{2}-\frac{1}{2}x_{1}x_{2}+\frac{1}{4}x_{2}^{2}.
%     \]
%     Рассмотрим выражение для $s_x^2$ как случайную функцию двух случайных переменных
%     $x_1$ и $x_2$, и найдём её математическое ожидание
%     (то есть усредним $s_x^2$ по большому числу опытов, в каждом из которых
%     проведено по по $n=2$ измерений).
%     Учитывая, что $x_1$ и $x_2$~--- независимы
%     ($\limaverage{x_1x_2}=\limaverage{x}_1\cdot \limaverage{x}_2$), получим
%     \[
%     \limaverage{s_{x}^{2}}=\frac12 \limaverage{x^2}~--- \frac12 \limaverage{x}^2.
%     \]
%     Сравнивая с известной формулой
%     $\sigma^2 = \limaverage{x^2}~--- \limaverage{x}^2$,
%     видим, что среднее значение оценки отличается от истинного в 2 раза.
% }


\paragraph{Эффективность оценки.}
Для сравнения разных методов оценки очень важным свойством является
их \emph{эффективность}. На качественном уровне эффективность~--- величина,
обратная разбросу значений $\hat{\theta}(\mathbf{x})$ при применении к разным
наборам данных $\mathbf{x}$. Как обсуждалось выше, оценка $\hat{\theta}(\mathbf{x})$
есть случайная величина, подчиняющаяся некоторому, в общем случае неизвестному,
распределению. Среднее $\overline{\hat{\theta}(\mathbf{x})}$ по этому распределению
определяет смещение оценки.
А его дисперсия $\sigma^2\left(\hat\theta\right)$~--- как раз мера ошибки
в определении параметра. Выбирая между различными методами (минимума хи-квадрат, максимального правдоподобия, наименьших квадратов и т. д.), мы, естественно,
хотим, чтобы ошибка была минимальной. Разные статистические методы обладают
разной эффективностью и в общем случае при конечном $n$ величина
$\sigma^2\left(\hat\theta\right)$ никогда не будет равна нулю.

Теорема, устанавливающая максимальное значение эффективности оценки,
рассмотрена в п.~\ref{sec:rao}.
% Разумеется, встает
% вопрос о том, можно ли построить оценку с максимальной возможной
% эффективностью.


\section{Максимальная эффективность оценки (граница Рао--Крамера)}\label{sec:rao}

Максимальная эффективность оценки ограничена теоремой Рао--Крамера.
\paragraph{Утверждение.}
Пусть оценка $\hat\theta$ параметра $\theta$ является несмещённой, тогда
всегда выполняется неравенство:
\begin{equation}
  \sigma^2(\hat\theta) \geq \frac{1}{I(\theta)},
\end{equation}
где
\begin{equation}
    I(\theta) =
    \limaverage{\left(\frac{\partial \ln L}{\partial \theta}\right)^2}.
\end{equation}
Здесь $L(\mathbf{y},\,\theta)$~--- введённая в п.~\ref{sec:chi2} функция правдоподобия
(вероятность получить набор результатов $\mathbf{y}$ при заданном параметре
$\theta$). Функцию $I(\theta)$ также называют \emph{информацией Фишера}.

\paragraph{Доказательство в одномерном случае.}
Обозначим
\[
U \equiv \frac{\partial \ln L} {\partial \theta} =
\frac{1}{L}\frac{\partial L}{\partial \theta}
\]
и найдём математическое ожидание этой функции:
\[
\limaverage{U} =
\int U\cdot L  d\mathbf{y} =
\int \frac{\partial L}{\partial \theta} d\mathbf{y} =
\frac{\partial}{\partial \theta} \int{L d\mathbf{y}} = 0.
\]
Теперь рассмотрим ковариацию параметра $\theta$ и функции
$U$:
\begin{equation}
\limaverage{\hat\theta\cdot U} 
= \frac{\partial}{\partial \theta} \int{\hat\theta L d\mathbf{y}} =
\frac{\partial \limaverage{\hat\theta}}{\partial \theta}.
\end{equation}
Для несмещенных оценок математическое ожидание оценки параметра равно
самому значению параметра: $\limaverage{\hat{\theta}} = \theta$,
поэтому последнее выражение есть просто единица.
Согласно неравенству Коши--Буняковского имеем
\[
\sigma^2 (\hat\theta) \cdot \sigma^2 (U) \geq
\left| \limaverage{\hat\theta \cdot U} \right| = 1,
\]
откуда и следует сделанное утверждение.

\paragraph{Следствие.}
Максимальная эффективность достигается в том случае, если величины
$\hat\theta$ и $U$ \emph{коррелируют} друг с другом.
Оценка, максимизирующая функцию $L(\mathbf{y},\theta)$
(метод максимального правдоподобия), является \emph{состоятельной},
\emph{несмещенной}, кроме того совпадает с оценкой вида
$U(\mathbf{y}, \hat\theta) = 0$, а значит является \emph{максимально эффективной}.

% \section{Принцип максимального правдоподобия и наилучшая оценка среднего}
%
% Пусть при измерениях одно и той же величины два студента
% независимым образом получили результаты
% \[
% x_{1}\pm\sigma_{1}\qquad\text{и}\qquad x_{2}\pm\sigma_{2}.
% \]
% Можно ли как-то объединить их ответы и таким образом улучшить оценку
% измеряемой величины?
%
% Первое, что может прийти на ум --- найти среднее
% арифметическое $\average{x}=\frac{x_{1}+x_{2}}{2}.$ Однако нетрудно
% понять, что если, скажем, измерение 2 сильно хуже, чем 1 ($\sigma_{2}\gg\sigma_{1}$),
% то разумнее было бы значение $x_{2}$ вообще отбросить и использовать
% $x_{1}$ как \textquote{наилучшую} оценку.
%
% Предположим, что измеряемая величина имеет нормальное распределение
% с некоторым средним $x_{0}$. Вероятность получить значение $x_{1}$
% при первом измерении согласно (\ref{eq:normal}) пропорциональна величине
% \[
% P_{1}\propto\frac{1}{\sigma_{1}}e^{-\left(x_{1}-x_{0}\right)/2\sigma_{1}^{2}},
% \]
% вероятность получить значение $x_{2}$ при втором измерении:
% \[
% P_{2}\propto\frac{1}{\sigma_{2}}e^{-\left(x_{2}-x_{0}\right)^{2}/2\sigma_{2}^{2}}.
% \]
% Вероятность получить пару значений $\left\{ x_{1},\,x_{2}\right\} $
% пропорциональна произведению $P_{1}P_{2}$:
% \[
% P\left(x_{1},x_{2}\right)=P_{1}P_{2}\propto\frac{1}{\sigma_{1}\sigma_{2}}e^{-\left(x_{1}-x_{0}\right)^{2}/2\sigma_{1}^{2}-\left(x_{2}-x_{0}\right)^{2}/2\sigma_{2}^{2}}.
% \]
%
% Рассмотрим выражение, оказавшееся в показателе экспоненты:
% \[
% \chi^{2}=\left(\frac{x_{1}-x_{0}}{\sigma_{1}}\right)^{2}+\left(\frac{x_{2}-x_{0}}{\sigma_{2}}\right)^{2}.
% \]
% Назовём \emph{наилучшей} такую оценку параметра
% $x_{0}$, при котором полученные в опытах результаты имеют \emph{максимальную
% вероятность} ($P\to\mathrm{max}$). Такой подход, называемый
% \emph{принципом максимального правдоподобия}, используется
% во многих задачах статистики.
%
% Из полученного выше видно, что $P\to\mathrm{max}$, если
% величина $\chi^{2}$ (\emph{хи-квадрат}) будет иметь
% \emph{минимум}. Дифференцируя по $x_{0}$ и приравнивая
% результат к нулю, запишем
% \[
% \frac{d\chi^{2}}{dx_{0}}=-2\frac{x_{1}-x_{0}}{\sigma_{1}^{2}}-2\frac{x_{2}-x_{0}}{\sigma_{2}^{2}}=0,
% \]
% откуда найдём
% \begin{equation}
% x_{0}=\frac{w_{1}x_{1}+w_{2}x_{2}}{w_{1}+w_{2}},\qquad\text{где }w_{1,2}=\frac{1}{\sigma_{1,2}^{2}}.
% \end{equation}
%
% Таким образом, для вычисления \emph{наилучшей} (максимально правдоподобной) оценки
% среднего нужно вычислить \emph{взвешенное среднее} с весами, \emph{обратно
% пропорциональными квадратам соответствующих погрешностей}.
%
% Результат непосредственно обобщается на произвольное число измерений:
% \begin{equation}
% x_{\text{наил}}=\frac{\sum\limits _{i}w_{i}x_{i}}{\sum\limits _{i}w_{i}},\qquad w_{i}=\sigma_{i}^{-2}.
% \end{equation}

% \section{Метод максимального правдоподобия для построения наилучшей
% прямой\label{subsec:MMP}}
%
% При описании метода наименьших квадратов мы не обосновали,
% почему и в каком смысле именно этот метод является \textquote{наилучшей}
% оценкой для коэффициентов линейной регрессии $y=kx+b$. Кроме того,
% мы получили формулы только для частного случая, когда погрешности
% всех экспериментальных точек одинаковы: $\sigma_{y}=\mathrm{const}$.
%
% Рассмотрим более общий случай. Пусть по-прежнему погрешности
% по оси абсцисс малы, $\sigma_{x}\to0$, а погрешности по оси $y$
% различны для каждой точки и равны $\sigma_{y_{i}}$. Пусть теория
% предсказывает \emph{линейную}\footnote{Отметим, что метод легко обобщается и на нелинейные зависимости общего вида $y=f\left(x;a,b,\ldots\right)$. Хотя формулы
% получаются существенно более громоздкими, при вычислении на компьютере
% оперировать с ними не сложнее, чем с линейной регрессией. В учебном
% практикуме мы рекомендуем всегда делать замену переменных, сводящую
% теоретическую зависимость к линейной, поскольку проведение прямой
% наиболее наглядно и может быть в грубом приближении проделано просто
% \textquote{по линейке}.}
% зависимость $y=kx+b$.
%
% Отклонение точки от теоретической зависимости обозначим как
% \[
% \Delta y_{i}=y_{i}-\left(kx_{i}+b\right).
% \]
%
% Воспользуемся \emph{принципом максимального правдоподобия}
% и построим такую прямую, чтобы вероятность обнаружить наблюдаемые
% в опыте отклонения $\left\{ \Delta y_{i}\right\} $ от неё была максимальна.
%
% Обозначим вероятность отклонения на величину $\Delta y_{i}$
% при известном $\sigma_{y_{i}}$ как $P\!\left(\Delta y_{i};\sigma_{y_{i}}\right)$.
% Предположим, что ошибки измерения для всех экспериментальных точек
% можно считать \emph{случайными} и \emph{независимыми}.
% В таком случае вероятность отклонения для всех $n$ точек равна произведению
% вероятностей, так что метод максимального правдоподобия сводится к
% поиску максимума выражения
% \begin{equation}
% \prod\limits _{i=1}^{n}P\!\left(\Delta y_{i};\sigma_{y_{i}}\right)\to\mathrm{max}.\label{eq:MMP_general}
% \end{equation}
% Максимизация производится по параметрам аппроксимирующей функции (на
% нашем случае это $k$ и $b$).
%
% Рассмотрим частный случай, когда погрешности имеют \emph{нормальное}
% (гауссово) распределение (\ref{eq:normal}) (напомним, что нормальное
% распределение применимо, если отклонения возникают из-за большого
% числа независимых факторов, что на практике реализуется довольно часто).
% Тогда, поскольку гауссова функция распределения пропорциональна величине
% $\propto e^{-\Delta y^{2}/2\sigma^{2}}$, выражение (\ref{eq:MMP_general})
% достигает максимума, если минимальна сумма
% \begin{equation}
% \boxed{\chi^{2}=\sum_{i=1}^{n}\frac{\Delta y_{i}^{2}}{\sigma_{y_{i}}^{2}}\to\mathrm{min}}.\label{eq:chi2}
% \end{equation}
% Здесь мы ввели стандартное обозначение для такой суммы ---
% $\chi^{2}$ (\emph{хи-квадрат)}.
%
% Таким образом, задача построения наилучшей прямой сводится
% к минимизации суммы квадратов отклонений, нормированных на соответствующие
% дисперсии $\sigma_{y_{i}}^{2}$. Если все погрешности одинаковы, $\sigma_{y_{i}}=\mathrm{const}$,
% мы приходим к методу наименьших квадратов.
%
% Получим выражения для наилучших коэффициентов $k$ и $b$.
% Заметим, что сумма (\ref{eq:chi2}) является \emph{взвешенной}
% суммой квадратов отклонений с весами
% \begin{equation}
% w_{i}=\frac{1}{\sigma_{y_{i}}^{2}}.
% \end{equation}
%
% Можно определить \emph{взвешенное среднее} от
% некоторого набора значений $\left\{ x_{i}\right\}$ как
% \[
% \left\langle x\right\rangle ^{\prime}=\frac{1}{W}\sum_{i}w_{i}x_{i},
% \]
% где $W=\sum\limits _{i}w_{i}$ --- нормировочная константа.
% Далее в этом разделе штрих будем для краткости опускать.
%
% Потребуем, согласно (\ref{eq:chi2}), чтобы была минимальна
% сумма
% \[
% \sum\limits _{i=1}^{n}w_{i}\Delta y_{i}^{2}\to\mathrm{min}.
% \]
% Повторяя процедуру, использованную при выводе (\ref{eq:MNK}), можно
% получить совершенно аналогичные формулы для оптимальных коэффициентов:
% \begin{equation}
% \boxed{k=\frac{\left\langle xy\right\rangle -\left\langle x\right\rangle \left\langle y\right\rangle }{\left\langle x^{2}\right\rangle -\left\langle x\right\rangle ^{2}},\qquad b=\left\langle y\right\rangle -k\left\langle x\right\rangle },\label{eq:MMP}
% \end{equation}
% с тем отличием, что под угловыми скобками $\left\langle \ldots\right\rangle $
% теперь надо понимать усреднение с весами $w_{i}=1/\sigma_{y_{i}}^{2}$.
%
% Найденные формулы позволяют вычислить коэффициенты линейной
% регрессии, \emph{если} известны величины $\sigma_{y_{i}}$.
% Значения $\sigma_{y_{i}}$ могут быть получены либо из некоторой теории,
% либо измерены непосредственно (многократным повторением измерений
% при каждом $x_{i}$), либо оценены из каких-то дополнительных соображений
% (например, как инструментальная погрешность).

\section{Погрешности коэффициентов построения прямой}\label{sec:sigma_kb}

Проведём подробный вывод для погрешностей коэффициентов наилучшей
прямой $\sigma_{k}$ и $\sigma_{b}$. Воспользуемся общей формулой
(\ref{eq:sigma_general}) для погрешности косвенных измерений. Считая,
что величины $x_{i}$ известны точно, запишем для погрешности углового
коэффициента
\[
\sigma_{k}^{2}=\sum\limits _{i}\left(\frac{\partial k}{\partial y_{i}}\right)^{2}\sigma_{y_{i}}^{2}.
\]
Продифференцируем (\ref{eq:MMP}) по $y_{i}$:
\[
\frac{\partial k}{\partial y_{i}}=\frac{1}{D_{xx}}\frac{\partial}{\partial y_{i}}\left(\frac{1}{W}\sum w_{i}x_{i}y_{i}-\left\langle x\right\rangle \frac{1}{W}\sum w_{i}y_{i}\right)=\frac{w_{i}\left(x_{i}-\left\langle x\right\rangle \right)}{WD_{xx}},
\]
где $D_{xx}=\average{x^{2}} -\average{x}^{2}$, $W=\sum_i \sigma_{y_i}^{-2}$.
Под угловыми скобками здесь понимается выборочное среднее с весами $w_i=1/\sigma_{y_i}^2$.
Тогда
\[
\sigma_{k}^{2}=\frac{1}{W^{2}D_{xx}^{2}}\sum\limits _{i}w_{i}^{2}\left(x_{i}-\left\langle x\right\rangle \right)^{2}\sigma_{y_{i}}^{2}.
\]
Учитывая, что $w_{i}\sigma_{y_{i}}^{2}=1$, получим
\begin{equation}
\sigma_{k}^{2}=\frac{1}{W D_{xx}}.\label{eq:MMP_sigma_k}
\end{equation}

Аналогично, для погрешности свободного члена имеем
\[
\sigma_{b}^{2}=\sum_{i}\left(\frac{\partial b}{\partial y_{i}}\right)^{2}\sigma_{y_{i}}^{2},
\]
где 
\[
\frac{\partial b}{\partial y_{i}}=\frac{w_{i}}{W}+\frac{\partial k}{\partial y_{i}}\left\langle x\right\rangle =\frac{w_{i}}{W}\left(1-\frac{x_{i}-\left\langle x\right\rangle }{\left\langle x^{2}\right\rangle -\left\langle x\right\rangle ^{2}}\left\langle x\right\rangle \right)=\frac{w_{i}}{W}\frac{\left\langle x^{2}\right\rangle -x_{i}\left\langle x\right\rangle }{D_{xx}}.
\]
Отсюда, пользуясь (\ref{eq:MMP_sigma_k}), приходим к формуле (\ref{eq:MNK_sigma_b}):
\begin{equation}
\sigma_{b}^{2}=\sigma_{k}^{2}
\frac{\left\langle \left(\left\langle x^{2}\right\rangle -x\left\langle x\right\rangle \right)^{2}\right\rangle }%
{D_{xx}}=
\sigma_{k}^{2}\left\langle x^{2}\right\rangle .
\end{equation}

\paragraph{Случай $\sigma_{y}=\mathrm{const}$.}
В частном случае метода наименьших квадратов (п. \ref{sec:MNK}),
формула (\ref{eq:MMP_sigma_k}) упрощается:
\begin{equation}
\sigma_{k}^{2}=\frac{\sigma_{y}^{2}}{nD_{xx}},\qquad\sigma_{b}^{2}=\sigma_{k}^{2}
\average{x^{2}}.\label{eq:MMP_sigma_k_simple}
\end{equation}
Здесь величина $\sigma_{y}$ может быть оценена непосредственно из
экспериментальных данных:
\begin{equation}
\sigma_{y}\approx\sqrt{\frac{1}{n-2}\sum_{i}\Delta y_{i}^{2}},\label{eq:MMP_sigma_y}
\end{equation}
где $n-2$ --- число \textquote{степеней свободы}
для приращений $\Delta y_{i} = y_i - (kx_i+b)$ ($n$ точек за вычетом двух связей
(\ref{eq:MMP})).

Формул (\ref{eq:MMP_sigma_k_simple}) и (\ref{eq:MMP_sigma_y}),
вообще говоря, достаточно для вычисления погрешности величины $k$
по известным экспериментальным точкам. Однако часто их объединяют
в одно упрощённое выражение. Для этого преобразуем (\ref{eq:MMP_sigma_y})
следующим образом: учитывая, что $\left\langle y\right\rangle =k\left\langle x\right\rangle +b$,
запишем
\[
\Delta y_{i}=y_{i}-kx_{i}-b=\left(y_{i}-\left\langle y\right\rangle \right)-k\left(x_{i}-\left\langle x\right\rangle \right).
\]
Возведём в квадрат, усредним и воспользуемся выражением для $k$ в
форме (\ref{eq:MNK_short}):
\[
\left\langle \Delta y^{2}\right\rangle =D_{yy}+k^{2}D_{xx}-2kD_{xy}=D_{yy}-k^{2}D_{xx}.
\]
Таким образом,
\[
\sigma_{y}=\sqrt{\frac{n}{n-2}\left(D_{yy}-k^{2}D_{xx}\right)},
\]
и с помощью (\ref{eq:MMP_sigma_k_simple}) получаем формулы (\ref{eq:MNK_sigma_k}),
(\ref{eq:MNK_sigma_b}):
\[
\boxed{\sigma_{k}=\sqrt{\frac{1}{n-2}\left(\frac{D_{yy}}{D_{xx}}-k^{2}\right)},\qquad\sigma_{b}=\sigma_{k}\sqrt{\left\langle x^{2}\right\rangle }}.
\]

\section{Многопараметрические оценки}
\label{sec:multiparam}

Однопараметрические оценки просты для понимания и реализации, но относительно
редко встречаются на практике. Даже при оценке параметров линейной зависимости
$y = k x + b$ требуется уже два параметра: наклон $k$ смещение $b$.
Все рассмотренные выше методы нахождения оптимальных параметров работают и
в многомерном случае, но поиск экстремума функций
(например, максимума функции правдоподобия или минимума суммы квадратов)
и интерпретация результатов требуют, как правило, использования численных методов.

\subsection{Двумерный случай}
Остановимся подробнее на построении прямой. Пусть некоторым методом получены
точечные оценки для наилучших значений $\hat{k}$ и $\hat{b}$.
Однако самих значений мало --- нас интересует область, в которой могут
оказаться параметры $k$, $b$ с некоторой доверительной вероятностью
(например, $P=0,68$) --- двумерная доверительная область.

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

Если бы $k$ и $b$ были независимы, достаточно было бы найти среднеквадратичные
отклонения $\sigma_k$ и $\sigma_b$, как это сделано в п. \ref{sec:sigma_kb}:
тогда искомая доверительная область значений параметров на плоскости $(k,b)$
представляла бы собой эллипс, оси которого параллельны координатным
(см. рис.~\ref{fig:kb}а).

Однако, если взглянуть, к примеру, на рис. \ref{fig:graph-method}б, иллюстрирующий
графический метод построения прямой, можно убедиться, что при варьировании
наклона $k$ обязательно меняется и смещение $b$. То есть параметры
$(k,\,b)$ вообще говоря \emph{скореллированы}. Количественно отклонения параметров
будут характеризоваться ковариационной матрицей:
\[
D = \left(\begin{matrix}
    D_{kk} & D_{kb} \\
    D_{bk} & D_{bb}
\end{matrix}\right),
\]
где $D_{kk}=\sigma_k^2$, $D_{bb}=\sigma_b^2$ --- дисперсии искомых параметров,
а 
\begin{equation*}
    D_{kb}=D_{bk} = \average{(k-\average{k})\cdot(b-\average{b})} = \rho_{kb}\sigma_k\sigma_b.  
\end{equation*}

Коэффициент $\rho_{kb}$ называют коэффициентом корреляции и он служит показателем \textquote{связанности} параметров. Для полностью независимых параметров он равен нулю, а в случае, если параметры нельзя отличить друг от друга --- единице.

По известной теореме линейной алгебры, симметричную матрицу можно
привести к диагональному виду поворотом координатных осей. Поэтому доверительная
область в таком случае будет представлять собой \emph{наклонный} эллипс
(см. рис.~\ref{fig:kb}б), а наклон его осей будет определяться
коэффициентом корреляции $r_{kb}$.
% Уравнение эллипса:
% \[
% \frac{k^2}{D_{kk}} - 2r_{kq}^2 k b + \frac{b^2}{D_{bb}}
% \]

\begin{figure}[h]
    \centering
    \input{images/kb.pdf_tex}
    \caption{Доверительная область значений коэффициентов прямой а)~$k$ и~$b$ независимы,
б)~$k$ и~$b$ скоррелированы.}
    \label{fig:kb}
\end{figure}


\subsection{Многомерный случай}

Принцип построения доверительной области в многомерном случае точно
такой же, как и для одномерных доверительных интервалов. Требуется найти
такую областью пространства параметров $\Omega$, для которой
вероятностное содержание для оценки параметра $\hat \theta$
% (или
% самого параметра $\theta$ в зависимости от того, какой философии вы
% придерживаетесь)
будет равно некоторой наперед заданной
величине $\alpha$:
\begin{equation}
    P(\theta \in \Omega) = \int\limits_\Omega{L(\mathbf{x} | \theta)}d\Omega = \alpha.
\end{equation}

Построение многомерной доверительной области на практике сталкивается с тремя
проблемами:
\begin{enumerate}
\item Взятие многомерного интеграла от произвольной функции~--- не тривиальная
  задача. Даже в случае двух параметров требуется владение
  методами вычислительной математики. Соответствующие методы реализованы
  в специализированных программных пакетах.
\item Определение центрального интервала для многомерной гиперобласти
является неоднозначной задачей.
\item Даже если удалось получить доверительную область, описать многомерный
объект в общем случае непросто, так что представление результатов
  составляет определенную сложность.
\end{enumerate}

Для решения этих проблем пользуются следующим приемом: согласно
центральной предельной теореме, усреднение большого количества одинаково
распределенных величин дает нормально распределенную величину. Это же
верно и в многомерном случае. В большинстве случаев, мы ожидаем, что
функция правдоподобия будет похожа на многомерное нормальное
распределение:
\begin{equation}
    L(\theta) = \frac{1}{(2 \pi)^{n/2}\left|D\right|^{1/2}} e^{-\frac{1}{2}
        (\mathbf{x} - \limaverage{\mathbf{x}})^T D^{-1} (\mathbf{x} - \limaverage{\mathbf{x}})},
\end{equation}
где $n$~--- размерность вектора параметров, $\limaverage{\mathbf{x}}$~--- вектор
наиболее вероятных значений, а $D$~--- ковариационная матрица распределения.

Для многомерного нормального распределения, линии постоянного уровня (то
есть поверхности, на которых значение плотности вероятности одинаковые)
имеют вид гиперэллипса, определяемого уравнением
$(\mathbf{x} - \limaverage{\mathbf{x}})^T D^{-1} (\mathbf{x} - \limaverage{\mathbf{x}}) = \mathrm{const}$.
Для любого вероятностного содержания $\alpha$ можно подобрать эллипс, который будет
удовлетворять условию на вероятностное содержание. Интерес, правда,
представляет не сам эллипс (в случае размерности больше двух, его просто
невозможно наглядно изобразить), а \emph{ковариацонная матрица}. Диагональные элементы
этой матрицы являются дисперсиями соответствующих параметров (с учетом
всех корреляций параметров).

\subsection{Использование пакета \texttt{scipy} для построение оценки}

Существует огромное количество программных пакетов для построения численной оценки параметров. Наиболее доступным и широко используемым является пакет \texttt{scipy} на языке Python. Приведем здесь только пример вызова процедуры оптимизации.

Пусть есть экспериментальные данные, представленные в виде трех колонок: $x$, $y$ и $err$. Требуется построить наилучшую прямую, описывающую эти данные.
Код для этого будет выглядеть следующим образом:
\begin{minted}[frame=lines, fontsize=\footnotesize]{python}
    from scipy.optimize import curve_fit

    function = lambda x, a, b: a*x + b

    popt, pcov = curve_fit(function, xdata = x, ydata = y, sigma = err)
\end{minted}

После выполнения этого кода, переменная \texttt{popt} содержит массив из двух значений, соответствующих оценке \texttt{a} и \texttt{b}, а переменная \texttt{pcov} содержит ковариационную матрицу для полученных параметров. 

Погрешности параметров можно получить как корни из диагональных элементов ковариационной матрицы:
\begin{minted}[frame=lines, fontsize=\footnotesize]{python}
    import numpy as np
    
    sigma_a = np.sqrt(pcov[0,0])
    sigma_b = np.sqrt(pcov[1,1])
\end{minted}

\note{Следует отметить, что существует огромное количество способов оценки оптимальных значений параметров и ковариационной матрицы. Поэтому при использовании того или иного инструмента, всегда следует сверяться с документацией и выяснять, что именно он делает. Также следует всегда проверять результаты обработки из качественных, \textquote{наивных} соображений.}



% \section{Проверка гипотез}
%
% Предположим, что теория предсказывает некоторую зависимость
% \[
% y=f\!\left(x;a,b,\ldots\right),
% \]
% а в эксперименте получен набор значений $\left\{ x_{i},\,y_{i}\right\} $.
% Метод максимального правдоподобия позволяет получить параметры
% $\left\{ a,b,\ldots\right\} $
% функции $f$, \textquote{наилучшим} образом приближающие
% экспериментальные значения. Причём, как бы плохо экспериментальные
% точки не ложились на теоретическую кривую, ответ будет получен в любом
% случае. Как проверить, действительно ли измеряемые величины можно
% считать связанными зависимостью $y=f\!\left(x\right)$?
%
% Задача может быть решена, если, как обычно, сделать ряд упрощающих
% предположений\footnote{Отметим, что в общем случае такая проверка не осуществима:
% в частности, если погрешности экспериментальных точек велики или не
% известны, то через них с равным \textquote{успехом}
% можно провести почти \emph{любую} функцию!}.
% Пусть опять все ошибки измерения \emph{независимы},
% распределены \emph{нормально} и нам известны их
% среднеквадратичные значения $\sigma_{y_{i}}$
% (или хотя бы их грубые оценки).
%
% Рассмотрим определённую выше сумму \emph{хи-квадрат} (\ref{eq:chi2})
% как функцию $n$ переменных $\left\{ y_{i}\right\}$:
% \begin{equation}
% \chi^{2}\!\left(y_{1},\,y_{2},\,\ldots,\,y_{n}\right)=\sum\limits _{i=1}^{n}\left(\frac{y_{i}-f\!\left(x_{i}\right)}{\sigma_{y_{i}}}\right)^{2}.\label{eq:chi2-1}
% \end{equation}
% Пусть функция $f\!\left(x\right)$ содержит $p$ \textquote{подгоночных}
% параметров (например, $p=2$ для линейной зависимости $f\!\left(x\right)=kx+b$).
% Найдём их наилучшие значения по методу максимального правдоподобия
% ($\chi^{2}\to\mathrm{min}$), и зафиксируем их. После этого $\chi^{2}$
% можно рассматривать как функцию
% \[
% m=n-p
% \]
% независимых переменных. Величину $m$ назовём \emph{числом степеней свободы} задачи.
%
% Попробуем сперва качественно ответить на вопрос, какое значение
% величины $\chi^{2}$ можно ожидать, если зависимость $y=f\!\left(x\right)$
% справедлива? Ясно, что если распределение ошибок нормальное, можно
% ожидать отклонений порядка среднеквадратичного: $\Delta y_{i}=y_{i}-f\!\left(x_{i}\right)\sim\sigma_{y_{i}}$.
% Поэтому значение суммы (\ref{eq:chi2-1}) должно оказаться порядка
% числа входящих в неё независимых слагаемых: $\chi_{m}^{2}\sim m$.
% В теории вероятностей доказывается (см., например, \cite{hudson}), что ожидаемое среднее значение (математическое ожидание) $\chi^{2}$ в точности равно числу степеней свободы: $\average{\chi_{m}^{2}}=m$.
%
% Теперь можно сформулировать качественный критерий проверки
% гипотезы о наличии некоторой функциональной зависимости (его называют
% \emph{критерий хи-квадрат}:
% \begin{itemize}
% \item  если $\chi^{2}\sim m$, согласие эксперимента с теорией \emph{удовлетворительное}
% (гипотеза не опровергнута);
% \item если $\chi^{2}\gg m$ --- \emph{согласия нет},
% то есть гипотеза о зависимости $y=f\!\left(x\right)$ скорее всего не верна.
% \end{itemize}
% Заметим, что если вдруг $\chi^{2}\ll m$, то совпадение \emph{слишком}
% хорошее, и скорее всего имеет место завышенная оценка для случайных
% погрешностей измерения $\sigma_{y_{i}}$.
%
% Для того, чтобы дать строгий \emph{количественный}
% критерий, с какой долей вероятности гипотезу $y=f\!\left(x\right)$
% можно считать подтверждённой или опровергнутой, нужно детально исследовать
% вероятностный закон, которому подчиняется функция $\chi^{2}$. В теории
% вероятностей он называется \emph{распределение хи-квадрат}
% (с $m$ степенями свободы). В элементарных функциях распределение
% хи-квадрат не выражается, но может быть легко найдено численно: функция
% встроена во все основные статистические пакеты, либо может быть вычислена
% по таблицам. Как правило, определяется вероятность
% $P\left(\chi^{2}>\chi_{0}^{2}\right)$
% того, что хи-квадрат имеет значение больше некоторого $\chi_{0}^{2}$,
% вычисленного из эксперимента. Если эта вероятность достаточно мала
% (например, $P<5\%$; конкретная величина доверительной вероятности
% всегда остаётся на усмотрение исследователя), соответствующую гипотезу
% следует признать несостоятельной.
%
% {\footnotesize
% \textbf{Пример.} Для данных на рис.~\ref{fig:correct}
% при $\sigma_{y}=0{,}2$ см величина хи-квадрат равна $\chi^{2}\approx4{,}7$.
% За вычетом двух параметров линейной аппроксимации имеем $m=6-2=4$
% степеней свободы. По графику TODO определяем, что вероятность того,
% что согласие окажется хуже, чем на TODO (т.е. разброс точек относительно
% \emph{той же} прямой будет больше),
% составляет $P\sim30\%$. Оснований для отказа от <<гипотезы
% о линейной зависимости>> нет.
%
% Если бы погрешность каждой точки была равна $\sigma_{y}=0{,}14$
% см, то мы получили бы $\chi^{2}=9{,}7$. Это соответствует $P<5\%$,
% так что считать зависимость линейной, по-видимому, было бы нельзя
% (либо не верна оценка для погрешности $\sigma_{y}$).\par
% }%\footnotesize
%
% Напоследок еще раз подчеркнём, что критерий хи-квадрат, во-первых,
% статистический и не может дать однозначного ответа --- только
% вероятностную оценку, а во-вторых, он работает корректно при условии,
% что ошибки разных точек независимы и каждая имеет нормальное распределение.
% Эти предположения, вообще говоря, выполняются далеко не всегда и,
% по-хорошему, требуют отдельной проверки.