В прикладной статистике метод наименьших полных квадратов (МНПК, TLS — англ. Total Least Squares ) — это вид , техника моделирования данных с помощью метода наименьших квадратов , в которой принимаются во внимание ошибки как в зависимых, так и в независимых переменных. Метод является обобщением регрессии Деминга и ортогональной регрессии и может быть применён как к линейным, так и нелинейным моделям.

Аппроксимация данных методом наименьших полных квадратов в общем виде эквивалентна лучшей по норме Фробениуса матрицы данных .

Линейная модель

Основы

В методе наименьших квадратов моделирования данных минимизируется функция потерь S ,

${\displaystyle S=\mathbf {r^{T}Wr} ,}$

где r — вектор , а W — весовая матрица. В линейном методе наименьших квадратов модель содержит уравнения, которые линейны от параметров в векторе ${\displaystyle {\boldsymbol {\beta }}}$ , так что отклонения вычисляются по формуле

${\displaystyle \mathbf {r=y-X{\boldsymbol {\beta }}} .}$

Имеется m наблюдений в векторе y и n параметров в β при m > n . X является m × n матрицей, элементы которой либо константы, либо функции от независимых переменных x . Весовая матрица W , в идеале, является обратной к дисперсионно-ковариационной матрице ${\displaystyle \mathbf {M} _{y}}$ наблюдений y . Предполагается, что независимые переменные ошибок не имеют. Параметры оценки находятся путём приравнивания градиента нулю, что приводит к уравнению

${\displaystyle \mathbf {X^{T}WX{\boldsymbol {\beta }}=X^{T}Wy} }$

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

Предположим теперь, что как x , так и y наблюдаются с ошибками с дисперсионно-ковариационными матрицами ${\displaystyle \mathbf {M} _{x}}$ и ${\displaystyle \mathbf {M} _{y}}$ соответственно. В этом случае функция потерь записывается как

${\displaystyle S=\mathbf {r_{x}^{T}M_{x}^{-1}r_{x}+r_{y}^{T}M_{y}^{-1}r_{y}} }$ ,

где ${\displaystyle \mathbf {r} _{x}}$ и ${\displaystyle \mathbf {r} _{y}}$ являются отклонениями для x и y соответственно. Ясно, что эти отклонения не могут быть независимыми и между ними должна быть какая-то связь. Если записать функцию как ${\displaystyle \mathbf {f(r_{x},r_{y},{\boldsymbol {\beta }})} }$ , ограничения выражаются m условиями .

${\displaystyle \mathbf {F=\Delta y-{\frac {\partial f}{\partial r_{x}}}r_{x}-{\frac {\partial f}{\partial r_{y}}}r_{y}-X\Delta {\boldsymbol {\beta }}=0} }$

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

${\displaystyle \mathbf {X^{T}M^{-1}X\Delta {\boldsymbol {\beta }}=X^{T}M^{-1}\Delta y} ,}$

или, альтернативно, ${\displaystyle \mathbf {X^{T}M^{-1}X{\boldsymbol {\beta }}=X^{T}M^{-1}y} }$

Здесь M — дисперсионно-ковариационная матрица, относящаяся как к независимым, так и зависимым переменным.

${\displaystyle \mathbf {M=K_{x}M_{x}K_{x}^{T}+K_{y}M_{y}K_{y}^{T};\ K_{x}=-{\frac {\partial f}{\partial r_{x}}},\ K_{y}=-{\frac {\partial f}{\partial r_{y}}}} }$

Пример

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

${\displaystyle f(x_{i},\beta )=\alpha +\beta x_{i}\!}$

И в этом случае

${\displaystyle M_{ii}=\sigma _{y,i}^{2}+\beta ^{2}\sigma _{x,i}^{2}}$

что показывает, как дисперсия в i -ой точке определяется дисперсией независимых и зависимых переменных, а также моделью, используемой для согласования данных. Выражение можно обобщить, если заметить, что параметр ${\displaystyle \beta }$ является наклоном прямой.

${\displaystyle M_{ii}=\sigma _{y,i}^{2}+\left({\frac {dy}{dx}}\right)_{i}^{2}\sigma _{x,i}^{2}}$

Выражение такого вида используется для аппроксимации , когда малые ошибки в x дают большие ошибки y в случае большого наклона.

С алгебраической точки зрения

Прежде всего следует заметить, что задача МНПК в общем случае решения не имеет, что было показано ещё в 1980 . Рассмотрим простой случай, где единственное решение существует без каких-либо предположений.

Вычисление МНПК с помощью сингулярного разложения описан в стандартных текстах . Мы можем решить уравнение

${\displaystyle XB\approx Y}$

относительно B , где X — матрица m -на- n , а Y — матрица m -на- k

То есть мы пытаемся найти матрицу B , минимизирующую матрицы ошибок R и F для X и Y соответственно. То есть

${\displaystyle \mathrm {argmin} _{R,F}\|[R\;F]\|_{F},\qquad (X+R)B=Y+F}$ ,

где ${\displaystyle [R\;F]}$ — с R и F рядом и ${\displaystyle \|\cdot \|_{F}}$ является нормой матрицы , квадратным корнем из суммы квадратов всех элементов матриц, что эквивалентно квадратному корню из суммы квадратов длин строк или столбцов матрицы.

Это можно переписать как

${\displaystyle [(X+R)\;(Y+F)]{\begin{bmatrix}B\\-E_{k}\end{bmatrix}}=0.}$

Где ${\displaystyle E_{k}}$ является ${\displaystyle k\times k}$ единичной матрицей. Целью является нахождение матрицы ${\displaystyle [R\;F]}$ , которая уменьшает ранг ${\displaystyle [X\;Y]}$ на k . Определим ${\displaystyle [U][\Sigma ][V]*}$ как сингулярное разложение расширенной матрицы ${\displaystyle [X\;Y]}$ .

${\displaystyle [X\;Y]=[U_{X}\;U_{Y}]{\begin{bmatrix}\Sigma _{X}&0\\0&\Sigma _{Y}\end{bmatrix}}{\begin{bmatrix}V_{XX}&V_{XY}\\V_{YX}&V_{YY}\end{bmatrix}}^{*}=[U_{X}\;U_{Y}]{\begin{bmatrix}\Sigma _{X}&0\\0&\Sigma _{Y}\end{bmatrix}}{\begin{bmatrix}V_{XX}^{*}&V_{YX}^{*}\\V_{XY}^{*}&V_{YY}^{*}\end{bmatrix}}}$ ,

где V разбита на блоки, соответствующие формам матриц X и Y .

Если использовать теорему Экарта-Янга , аппроксимация, минимизирующая норму ошибки, это такая аппроксимация, что матрицы ${\displaystyle U}$ и ${\displaystyle V}$ не меняются, в то время как ${\displaystyle k}$ наименьших сингулярных значений заменяются нулями. То есть мы хотим

${\displaystyle [(X+R)\;(Y+F)]=[U_{X}\;U_{Y}]{\begin{bmatrix}\Sigma _{X}&0\\0&0_{k\times k}\end{bmatrix}}{\begin{bmatrix}V_{XX}&V_{XY}\\V_{YX}&V_{YY}\end{bmatrix}}^{*}}$

так что, ввиду линейности,

${\displaystyle [R\;F]=-[U_{X}\;U_{Y}]{\begin{bmatrix}0_{n\times n}&0\\0&\Sigma _{Y}\end{bmatrix}}{\begin{bmatrix}V_{XX}&V_{XY}\\V_{YX}&V_{YY}\end{bmatrix}}^{*}.}$

Мы можем удалить блоки из матриц U и Σ, упростив выражение до

${\displaystyle [R\;F]=-U_{Y}\Sigma _{Y}{\begin{bmatrix}V_{XY}\\V_{YY}\end{bmatrix}}^{*}=-[X\;Y]{\begin{bmatrix}V_{XY}\\V_{YY}\end{bmatrix}}{\begin{bmatrix}V_{XY}\\V_{YY}\end{bmatrix}}^{*}.}$

Это даёт R и F , таки что

${\displaystyle [(X+R)\;(Y+F)]{\begin{bmatrix}V_{XY}\\V_{YY}\end{bmatrix}}=0.}$

Теперь, если ${\displaystyle V_{YY}}$ не вырождена, что не всегда верно (заметим, что поведение МНПК в случае вырожденности ${\displaystyle V_{YY}}$ не вполне понятно), мы можем умножить справа обе части на ${\displaystyle -V_{YY}^{-1}}$ , чтобы привести нижний блок правой матрицы к отрицательной единичной матрице, что даёт

${\displaystyle [(X+R)\;(Y+F)]{\begin{bmatrix}-V_{XY}V_{YY}^{-1}\\-V_{YY}V_{YY}^{-1}\end{bmatrix}}=[(X+R)\;(Y+F)]{\begin{bmatrix}B\\-E_{k}\end{bmatrix}}=0,}$

а тогда

${\displaystyle B=-V_{XY}V_{YY}^{-1}.}$

Имплементация в системе GNU Octave :

function B = tls(X,Y)

[m n]   = size(X);            % n является шириной матрицы X (X[m x n])
Z       = [X Y];              % Z является расширением X на Y.
[U S V] = svd(Z,0);           % находим [[Сингулярное разложение|SVD]] матрицы Z.
VXY     = V(1:n,1+n:end);     % Берём блок матрицы V, состоящий из первых n строк и n+1 последних столбцов
VYY     = V(1+n:end,1+n:end); % Берём нижний правы блок матрицы V.
B       = -VXY/VYY;

end

Метод решения задачи, описанный выше и требующий, чтобы матрица ${\displaystyle V_{YY}}$ не была вырожденной, может быть слегка расширен так называемым классическим МНПК алгоритмом .

Вычисление

Стандартная имплементация классического алгоритма МНПК доступна на , см. также статьи . Все современные имплементации, базирующиеся, например, на использовании обычного метода наименьших квадратов, аппроксимируют матрицу ${\displaystyle B}$ (которая в литературе обозначается как ${\displaystyle X}$ ), как это делают Ван Хуффель и Вандевалле. Стоит заметить, однако, что полученная матрица ${\displaystyle B}$ во многих случаях не является решением МНПК .

Нелинейная модель

Для похожие рассуждения показывают, что нормальное уравнение для итерационного цикла может быть переписано как

${\displaystyle \mathbf {J^{T}M^{-1}J\Delta {\boldsymbol {\beta }}=J^{T}M^{-1}\Delta y} .}$

Геометрическая интерпретация

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

Масштабно-инвариантные методы

Серьёзная трудность появляется, если переменные не измеряются в тех же самых единицах. Сначала рассмотрим измерение расстояния между точками данных и кривой — какова будет единица измерения для расстояния? Если мы будем измерять расстояние на основе теоремы Пифагора, ясно, что нам придётся складывать единицы, измеряемые в различных единицах, что приводит к бессмысленным результатам. Если мы заменим масштаб одной из переменных, например, будем измерять в граммах, а не килограммах, мы получим другие результаты (другую кривую). Чтобы избежать этой проблемы несоизмеримости, иногда предлагается переводить их в безразмерные величины — это можно назвать нормализацией или стандартизацией. Существуют, однако, различные пути сделать это, приводящие к неэквивалентным моделям. Один из подходов — нормализовать с помощью известной (или оценочной) точности измерения, минимизируя тем самым расстояние Махаланобиса до точек на линии и обеспечивая решение с максимальным правдоподобием . Неизвестные точности измерения могут быть найдены с помощью дисперсионного анализа .

Кратко, метод наименьших полных квадратов не имеет свойства инвариантности по единицам измерения, т.е. он не является масштабно инвариантным . Для полноценности модели мы требуем, чтобы это свойство выполнялось. Дальнейшее продвижение вперёд, это понимание, что отклонения (расстояния), измеряемые в других единицах, могут быть скомбинированы, если используется умножение, а не сложение. Рассмотрим аппроксимацию прямой, для каждой точки данных произведение горизонтального и вертикального отклонений равно удвоенной площади треугольника, образованного отрезками отклонений и аппроксимирующей прямой. Мы выбираем прямую, минимизирующую сумму этих площадей. Нобелевский лауреат Пол Самуэльсон доказал в 1942, что в двумерном случае эта прямая выражается исключительно в терминах отношений квадратических отклонений и корреляции коэффициентов, которые (1) удовлетворяют уравнению, если наблюдения находятся на прямой линии; (2) обнаруживают масштабную инвариантность, (3) обнаруживают инвариантность при обмене переменных . Эта прямая переоткрывалась в различных дисциплинах и известна как стандартизованная главная ось , приведённая главная ось, функциональная зависимость средних геометрических , регрессия наименьших квадратов, диагональная регрессия и прямая наименьших площадей. Тофаллис расширил этот подход для работы с несколькими переменными.

См. также

• Регрессия Деминга
• Линейная регрессия
• Метод наименьших квадратов

Примечания

Литература

Дополнительная литература

• Paige C. C., Strakoš Z.,. Core problems in linear algebraic systems // SIAM J. Matrix Anal. Appl. — 2006. — Т. 27 . — С. 861—875 .
• Jo S., Kim S. W. Consistent normalized least mean square filtering with noisy data matrix. — 2005. — Т. 53. — С. 2112—2123. — (IEEE Trans. Signal Processing).
• DeGroat R. D., Dowling E. M. The data least squares problem and channel equalization. — 1993. — Т. 41. — С. 407—411. — (IEEE Trans. Signal Processing).
• Abatzoglou T., Mendel J. Constrained total least squares. — 1987. — Т. 12. — С. 1485—1488. — (Proc. IEEE Int. Conf. Acoust., Speech, Signal Process. (ICASSP’87)).
• de Groen P. . — 1996. — С. 237—253. — (Nieuw Archief voor Wiskunde, Vierde serie, deel 14).
• at MathPages
• Amiri-Simkooei A. R., Jazaeri S. // Journal of Geodetic Science. — 2012. — Т. 2 (2) . — С. 113—124 .