Интерполяция, экстраполяция и сглаживание
или
Ложь, наглая ложь и статистика[1]

В.Ф. Очков, Е.П. Богомолова

Сначала дадим расшифровку некоторых терминов, которые помогут нам понять суть проблемы, поднимаемой в статье [1].

Функция (от лат. functio – исполнение, осуществление) – это закон f, по которому каждому элементу x из одного числового множества ставится в соответствие некоторый определенный элемент y из другого числового множества, так что y=f(x).

Корреляция (от лат. correlatio – соотношение) – вероятностная или статистическая зависимость, не имеющая строго функционального характера (см. выше) из-за невозможности точно учесть влияние множества одновременно меняющихся факторов.

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

Регрессия (от лат. regressio – обратное движение) – такая зависимость среднего значения какой-либо величины от некоторой другой величины или от нескольких величин, что одному и тому же значению независимой переменной x могут соответствовать в зависимости от случая различные значения величины y.

А теперь к делу.

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

Но есть и третий вариант: две величины действительно связаны функциональной зависимостью, и теоретически каждому значению аргумента x соответствует ровно одно значение y. Но при проведении эксперимента не удается получить значения с достаточной степенью точности, и экспериментальные результаты представляются в виде таблиц значений (x_i, y_i), y_i = f(x_i) ± ε_i, где ε_i – погрешности измерений. Если в этом случае просто соединить экспериментальные точки отрезками прямых, то получится ломаная, не имеющая ничего общего с той функциональной зависимостью, которая реально существует. Причем форма этой ломаной из-за ошибок измерений при повторном эксперименте не воспроизводится.

В первом из рассмотренных случаев для получения приближенной формулы зависимости двух величин применяют интерполяцию. Во втором и третьем – регрессивный анализ. Приведем примеры.

Студентам, решающим задачи по физике, химии и другим учебным дисциплинам, часто приходится заглядывать в справочники, где некая функциональная зависимость представлена в виде таблицы: теплопроводность какого-либо металла в зависимости от температуры, плотность водного раствора какой-либо соли при фиксированной концентрации зависит от той же температуры и т.д. и т.п.[2]. В этих таблицах пары чисел «температура-свойство» показаны лишь для некоторых значений температуры, например, 0, 10, 20, 30 и т.д. градусов по шкале Цельсия. А значение свойства (теплопроводности металла или плотности раствора) нужно, допустим, взять при 17 °С. Подразумевается, что величины, указанные в таблицах, являются точными (правильнее сказать, получены с достаточной степенью точности). Поскольку физические законы позволяют считать рассматриваемые величины функционально зависимыми, то применяют интерполяцию. Простейшую – линейную или более сложную – интерполяцию многочленами либо сплайнами. При линейной интерполяции мы мысленно, на бумаге или на дисплее компьютера фиксируем две точки, проводим через них прямую линию и по ней находим нужное промежуточное значение. При нелинейной интерполяции отыскивается обладающая определенной гладкостью функция, график которой проходит через указанные точки (интерполяционные узлы).

В случае отсутствия каких-либо определенных функциональных связей между двумя экспериментальными величинами интерполяция может привести к ложному результату.

Когда-то перед лекцией на тему «Регрессионный анализ» по курсу «Информационные технологии» первый автор статьи подбирал пример статистической выборки для такого анализа. Но когда он «взошел на кафедру» и взглянул на аудиторию (рис. 1), то он понял, что эта выборка находится прямо перед его глазами. Была проведена перекличка студентов. Юноши при вызове их по фамилии вставали и сообщали свой вес[3] и рост (студентки естественно, были исключены из этой процедуры). Данные заносились в два вектора с именами Вес и Рост в среде математической программы Mathcad (скачать Mathcad 15 файл с этими данными >>>). В векторах получилось по 50 элементов. Эти массивы чисел послужили хорошей затравкой для лекции.

Рис. 1. Студенты – "подопытные статистические кролики" (правая половина аудитории: часть выборки)

Линейная интерполяция, формально примененная к паре величин «Рост – Вес», показана на рис. 2. Были взяты данные двух студентов («маленький студент» и «большой студент») и через соответствующие точки (узлы интерполяции) была проведена прямая линия на плоскости. Были также определены параметры этой прямой – коэффициенты a и b уравнения y = a + b x.

Рис. 2. Линейная интерполяция по двум студентам, пардон, точкам

Как рассчитывались коэффициенты a и b, показано на рис. 3. Можно в справочниках поискать соответствующие формулы, а можно в среде Mathcad решить пару линейных уравнений (рис. 3) и получить ответ.

Рис. 3. Численное решение системы линейных уравнений прямой, проходящей через две точки на плоскости

На рисунке 3 в таблицу вводятся вес и рост "маленького и большого студентов", затем записывается система двух линейных алгебраических уравнений прохождения прямой через две точки, формируется квадратная матрица M коэффициентов при неизвестных a и b и вектор свободных членов v. Решение системы найдено с помощью встроенной в Mathcad функции lsolve. Нашу систему двух линейных уравнений можно решить и аналитически (рис. 4), получив формулу для определения коэффициентов a и b, входящих в уравнение прямой линии, проходящей через две заданные точки с координатами на плоскости x₁-y₁и x₂-y₂.

Рис. 4. Символьное решение системы линейных уравнений прямой, проходящей через две тоски на плоскости

Тут мы намерено смешали интерполяцию (проведение линии через точки) с аппроксимацией (проведение линии вблизи точек). Если выбирать пару студентов случайным образом, то наша прямая линия может быть различными способами ориентирована на плоскости «Рост – Вес», поэтому доверия к такому способу описания зависимости величин в этом случае нет. Линейная интерполяция тут не годится: она еще раз «льет воду на мельницу» тех, кто в шутку или всерьез утверждает, что есть «ложь, наглая ложь и… статистика». На график, показанный на рис. 2, нужно «высыпать» все 50 точек со статистическими данными студентов и постараться как-то провести прямую линию, отображающую взаимозависимость (корреляцию) роста и веса человека, не противоречащую основным тенденциям поведения измеренных величин – см. рис. 5.

Рис. 5. Линейное сглаживание (аппроксимация)

Когда-то давно в докалькуляторную и докомпьютерную эру такая работа делалась примерно так: на миллиметровой бумаге в определенном масштабе проставлялись точки исходных статистических данных, а затем с помощью прозрачной линейки выставлялась некая прямая линия, наилучшим (оптимальным) образом описывающая связь роста человека и его веса. При этом каждый проводящий такую линию мог руководствоваться какими-то своими субъективными понятиями об оптимальности положения этой прямой на плоскости, то есть о значениях коэффициентов а и b, задающих уравнение прямой линии у = а + b∙х (константы а и b, кстати, рассчитаны и прописаны на рис. 2 и 5). На компьютере эту работу можно автоматизировать, если выбрать какой-то объективный критерий оптимизации. Одним из таких критериев является критерий минимума суммы квадратов «вертикальных» отклонений точек от кривой (в нашем случае от прямой линии). Сам же метод такого расчета называется методом наименьших квадратов: МНК. Если возникнет спор о том, чья прямая линия, проведенная в массиве точек, более правильная, то можно циркулем и линейкой замерить отклонения по вертикали точек от прямой, возвести каждое такое отклонение в квадрат, просуммировать эти квадраты и сравнить полученные суммы. У кого эта сумма окажется меньше, у того и прямая линия правильнее иллюстрирует исследуемую статистическую зависимость. Если у нас есть только две точки (рис. 2), или точек больше двух, но все они лежат на одной прямой, то, очевидно, сумма квадратов отклонений будет равна нулю. Но если же не все точки лежат на одной прямой (рис. 5), то можно предположить, что есть такое положение сглаживающей прямой (такие значения коэффициентов а и b), при котором сумма квадратов отклонений будет минимальной. На рис. 6 показан Mathcad-документ, решающий эту задачу оптимизации.

Рис. 6. Метод наименьших квадратов как задача оптимизации

В Mathcad-документе, показанном на рис. 6, формируется функция пользователя с именем СКО (сумма квадратов отклонений масс студентов – целевая функция оптимизации), имеющая два аргумента – а и b. Далее задаются начальные значения этих переменных оптимизации с опорой на известную эмпирическую формулу, гласящую, что рост взрослого человека в сантиметрах – это его вес в килограммах плюс сто. Наша задача – проверить эту гипотезу на данной статистической выборке – на 50 студентах. Затем в расчете на рис. 6 вызывается встроенная в Mathcad функция Minimize, которая начинает менять заданные нами значения переменных а и b так, чтобы функция СКО приняла минимальное значение. Как это делается – это отдельный разговор. Мы же видим, что при начальных значениях коэффициентов (переменных) а и b функция СКО возвращала значение 3396 kg², а при оптимальных значениях – 1505 kg², и верим что меньшего значения в пределах заданной точности найти не удаться. На рис. 5 можно видеть эти оптимальные значения коэффициентов а и b и «оптимальную» прямую линию, «пронизывающую» массив точек. К оптимальности этого решения мы еще вернемся ниже.

Задачу, сформулированную на рис. 6, можно решить и аналитически, если вспомнить, что у функции двух аргументов в точке минимума частные производные по обеим переменным равны нулю. На рисунке 7 показано символьное решение системы двух уравнений равенства нулю частных производных функции СКО(a, b) по аргументам a и b.

Рис. 7. Символьное решение задачи оптимизации – поиска минимума функции двух аргументов

У функции СКО суммируются квадраты «вертикальных»[4] отклонений точки от прямой линии. И это понятно: если убрать квадрат, то суммирование приведет к некорректной постановке задачи: будут суммироваться отклонения с разными знаками. А можно ли возведение в квадрат заменить взятием абсолютного значения? Для этого достаточно внести соответствующее изменение в первый оператор расчета, показанного на рис. 6. Но в этом случае функция Minimize будет иметь плохую сходимость и выдавать совершенно разные ответы при разных начальных значениях переменных a и b. Дело в том, что функция абсолютного значения угловатая, если так можно выразится. А природа не любит острых углов. Кроме того, стоит взглянуть на рис. 8, где сделана неудачная попытка аналитического решения системы двух уравнений – равенство нулю частных производных не от функции СКО(a, b), а от функции СМО(a, b): Сумма Модулей Отклонений.

Рис. 8. Неудачное решение системы двух уравнений

На сайте http://www.math.wpi.edu/Course_Materials/SAS/lablets/7.3/7.3c/index.html можно видеть описание и сравнение методов наименьших квадратов и наименьших модулей.

Одна точка на рис. 5 обведена кружочком – это, можно сказать, некий типичный студент, отклонение веса которого от средней линии минимально.

В связи с этим возникает интересное предложение по подведению итогов различных конкурсов красоты. Сейчас в них слишком много субъективности, а значит – обид, слез и даже судебных тяжб. В финалы таких конкурсов обычно попадают «красавицы-раскрасавицы» (см., например, рис. 1), из которых довольно трудно выбрать самую оптимальную, пардон, самую красивую мисс или миссис. Так вот, можно у этих финалисток замерить вес и рост или другие размеры (пресловутые 90-60-90, например), провести через точки линию и назвать победительницу так, как это показано кружочком на рис. 5.

Но вернемся к нашим более серьезным делам. Встроенная в Mathcad функция Minimize, показанная на рис. 6, позволяет нам реализовать метод наименьших квадратов в его общей постановке: пользовательская (целевая) функция с именем СКО в принципе может иметь любое число аргументов, а ее правая часть – любую функциональную зависимость для сглаживания точек. Но для частных случаев в среде Mathcad есть более простые инструменты решения задачи статистической обработки данных методом наименьших квадратов. Так, на рисунке 9 показан вызов двух встроенных в Mathcad функций intercept и slope, возвращающих значение коэффициентов а и b сглаживающей прямой у = а + b∙х. В эти две функции заложены формулы, выведенные нами на рис. 7. Они также использованы в расчете на рис. 9. Кроме того, на рис. 9 показана работа еще одной функции, возвращающей коэффициенты линейной регрессии, – функции line. Английские слова an intercept и aslope можно перевести как «пересечение» и «наклон»: значение коэффициента а – это значение абсциссы при пересечении нашей прямой оси Y, а значение константы b – это тангенс угла наклона прямой к оси Х. Вернее так. Про тангенс и угол можно было бы упомянуть, если б наши исходные векторы были бы безразмерными. При размерных векторах (масса и длина, как в нашей задаче) это просто некий наклон – отношение приращения по вертикали (граммы) к приращению по горизонтали (сантиметры).

Рис. 9. Инструменты пакета Mathcadдля решения линейной регрессионной задачи методом наименьших квадратов

Через точки, показанные на рис. 5, можно провести только прямую линию. Но если исследователь видит в исходных разбросанных точках не прямую, а, например, параболу или полином более высокой степени с коэффициентами a, b, c, d и т.д., то он может вызвать Mathcad-функцию regress, у которой уже не два (см. рис. 9), а три аргумента: вектор ординат точек, вектор абсцисс точек и степень полинома: 1 (наш случай, показанный на рис. 5, 7, 8 и 9), 2 (парабола), 3 (полином третьей степени), 4 – четвертой (см. рис. 10 ниже) и т.д. до того момента, когда степень полинома станет равна числу точек минус единица. В этом случае аппроксимация (проведение кривой вблизи точек) перейдет в интерполяцию (проведение кривой через точки). Но тут не так все просто.

Авторская анимация на сайте http://communities.ptc.com/videos/1473 (рис. 10) показывает, как меняется форма аппроксимирующей кривой вблизи 15 точек и как добавляются новые коэффициенты k_i при изменении степени полинома от 1до 14.