CC BY
267
Электронный журнал Cloud of Science. 2015. T. 2. № 1
http:/ / cloudofscience.ru ISSN 2409-031X
Интерполяция, экстраполяция и сглаживание или «ложь, наглая ложь и статистика»1
В. Ф. Очков, Е. П. Богомолова
Национальный исследовательский университет «Московский энергетический институт», 111250, Москва, ул. Красноказарменная, 14
e-mail: [email protected]
Аннотация. В статье рассмотрены инструменты математической программы Mathcad, позволяющие решать задачи интерполяции и сглаживания. Обсуждаются вопросы правильного применения средств программы Mathcad для корректного решения задач с различными типами зависимости данных.
Ключевые слова: интерполяция, сглаживание, Mathcad, функциональная зависимость, статистическая зависимость.
Как справедливо было отмечено в [3, с. 241], «преувеличение роли информатики столь же, а может быть, и более опасно, чем преуменьшение этой роли». Подводя итоги последнего десятилетия, мы можем с уверенностью сказать, что современные студенты доверяют Википедии больше, чем своему лектору, а их вера в возможности вычислительных пакетов вообще не имеет границ. Не будучи способными разобраться в математических методах и алгоритмах, «спрятанных» в недрах математических вычислительных программ, молодые исследователи в своих работах допускают принципиальные и неустранимые ошибки. Выходом из такого положения является как можно более раннее (еще в школе) знакомство учащихся с базовыми основами численных методов на примерах, доступных пониманию школьников. Решение обыденных математических задач, переработанных для образования школьников с помощью информатики, и будет реальным шагом на пути фундаментализации обучения информатике [4].
Продемонстрируем реальную попытку объединения информатики и математики в рамках единой образовательной области, определим точки соприкосновения математики и информатики [5]. Приведем примеры результативного взаимодействия этих наук, используя понятные всем объекты и задачи.
Сначала дадим расшифровку некоторых терминов, которые помогут нам понять суть проблемы, поднимаемой в статье [6].
1 Эту фразу приписывают многим известным людям — запустите ее в Google и посмотрите ответ.
Функция (от лат. functio — исполнение, осуществление) — это закон f по которому каждому элементу x из одного числового множества ставится в соответствие некоторый определенный элемент у из другого числового множества, так что У = f (х).
Корреляция (от лат. correlatio — соотношение) — вероятностная или статистическая зависимость, не имеющая строго функционального характера (см. выше) из-за невозможности точно учесть влияние множества одновременно меняющихся факторов.
Интерполяция (от лат. interpolation — изменение, переделка) — приближенное или точное нахождение какой-либо величины по известным отдельным значениям этой величины, т. е. восстановление (точное или приближенное) функции (см. выше) по ее нескольким известным значениям.
Регрессия (от лат. regressio — обратное движение) — такая зависимость среднего значения какой-либо величины от некоторой другой величины или от нескольких величин, что одному и тому же значению независимой переменной x могут соответствовать в зависимости от случая различные значения величины у.
А теперь к делу.
При обработке экспериментальных данных мы обычно сталкиваемся с двумя принципиально разными случаями: или полученные данные связаны некоторой функциональной зависимостью (формула, которая нам неизвестна), или таковой функциональной зависимости не существует, хотя и наблюдается корреляция.
Но есть и третий вариант: две величины действительно связаны функциональной зависимостью, и теоретически каждому значению аргумента x соответствует ровно одно значение y. Но при проведении эксперимента не удается получить значения с достаточной степенью точности, и экспериментальные результаты представляются в виде таблиц значений (х, y ), y = f (х)± ef, где sf — погрешности измерений. Если в этом случае просто соединить экспериментальные точки отрезками прямых, то получится ломаная, не имеющая ничего общего с той функциональной зависимостью, которая реально существует. Причем форма этой ломаной из-за ошибок измерений при повторном эксперименте не воспроизводится.
В первом из рассмотренных случаев для получения приближенной формулы зависимости двух величин применяют интерполяцию. Во втором и третьем — регрессионный анализ. Приведем примеры.
Студентам, решающим задачи по физике, химии и другим учебным дисциплинам, часто приходится заглядывать в справочники, где некая функциональная зависимость представлена в виде таблицы: теплопроводность какого-либо металла в зависимости от температуры, плотность водного раствора какой-либо соли при
фиксированной концентрации зависит от той же температуры и т. д., и т. п.2 В этих таблицах пары чисел «температура-свойство» показаны лишь для некоторых значений температуры, например, 0, 10, 20, 30 и т. д. градусов по шкале Цельсия. А значение свойства (теплопроводности металла или плотности раствора) нужно, допустим, взять при 17 °С. Подразумевается, что величины, указанные в таблицах, являются точными (правильнее сказать, получены с достаточной степенью точности). Поскольку физические законы позволяют считать рассматриваемые величины функционально зависимыми, то применяют интерполяцию. Простейшую — линейную или более сложную — интерполяцию многочленами либо сплайнами. При линейной интерполяции мы мысленно, на бумаге или на дисплее компьютера фиксируем две точки, проводим через них прямую линию и на ней находим нужное промежуточное значение. При нелинейной интерполяции отыскивается обладающая определенной гладкостью функция, график которой проходит через указанные точки (интерполяционные узлы).
В случае отсутствия каких-либо определенных функциональных связей между двумя экспериментальными величинами интерполяция может привести к ложному результату.
Когда-то перед лекцией на тему «Регрессионный анализ» по курсу «Информационные технологии» первый автор статьи подбирал пример статистической выборки для такого анализа. Но когда он «взошел на кафедру» и взглянул на аудиторию (рис. 1), то понял, что эта выборка находится прямо перед его глазами. Была проведена перекличка студентов. Юноши при вызове их по фамилии вставали и сообщали свой вес3 и рост (студентки естественно, были исключены из этой процедуры). Данные заносились в два вектора с именами Вес и Рост в среде математической программы Mathcad. В векторах получилось по 50 элементов. Эти массивы чисел послужили хорошей затравкой для лекции.
Линейная интерполяция, формально примененная к паре величин «Рост — Вес», показана на рис. 2. Были взяты данные двух студентов («маленький студент» и «большой студент») и через соответствующие точки (узлы интерполяции) была проведена прямая линия на плоскости. Были также определены параметры этой прямой — коэффициенты a и b уравнения y=a+bx.
2Раньше приходилось определять по таблицам значение тригонометрических, показательных и других элементарных функций.
3
Имелась в виду, конечно, масса студента, а не его вес. В середине XX века, когда наша страна переходила на международную систему измерений (СИ), предписывающую вес измерять не в килограммах силы, а в ньютонах, школьники и студенты шутили примерно так — заходили в магазин и просили продавщицу отвесить 500 ньютонов колбасы.
Рисунок 1. Студенты — «подопытные статистические кролики» (правая половина аудитории — выборка из выборки)
Рисунок 2. Линейная интерполяция по двум точкам
Как рассчитывались коэффициенты a и b, показано на рис. 3. Можно в справочниках поискать соответствующие формулы, а можно в среде Mathcad решить пару линейных уравнений (см. рис. 3) и получить ответ.
На рис. 3 в таблицу вводятся вес и рост «маленького и большого студентов», затем записывается система двух линейных алгебраических уравнений прохождения прямой через две точки, формируется квадратная матрица M коэффициентов при неизвестных a и b и вектор свободных членов v. Решение системы найдено с помощью встроенной в Mathcad функции ¡solve. Нашу систему двух линейных уравнений можно решить и аналитически (рис. 4), получив формулу для определения коэффициентов a и b, входящих в уравнение прямой линии, проходящей через две заданные точки с координатами на плоскости (Xj, y ) и (х2 ,y2).
Рисунок 3. Численное решение системы линейных уравнений прямой, проходящей через две точки на плоскости
Рисунок 4. Символьное решение системы линейных уравнений прямой, проходящей через две тоски на плоскости
Тут мы намеренно смешали интерполяцию (проведение линии через точки) с аппроксимацией (проведение линии вблизи точек). Если выбирать пару студентов случайным образом, то наша прямая линия может быть различными способами ориентирована на плоскости «Рост — Вес», поэтому доверия к такому способу описания зависимости величин в этом случае нет. Линейная интерполяция тут не годится: она еще раз «льет воду на мельницу» тех, кто в шутку или всерьез утверждает, что есть «ложь, наглая ложь и... статистика»4. На график, показанный на рис. 2, нужно «высыпать» все 50 точек со статистическими данными студентов и
Это изречение приписывают премьер-министру Великобритании Дизраэли.
постараться как-то провести прямую линию, отображающую (коррелирующую) взаимозависимость роста и веса человека (см. рис. 5).
90 90 70 60 ••• Студенты __ • » • • • ^— Вес(Рос -,
ОО Типичный студент -55 "=50 •
ф ш • * • ' • • ¿а,--- ' * » » • т] = а + Ь Рост
• 1 -1 = -79.958 кд = 854.12 — ст
ЭО
170 180 190 Рост(ст) 2
Рисунок 5. Сглаживание прямой линией 50 точек
Когда-то давно в докалькуляторную и докомпьютерную эру такая работа делалась примерно так: на миллиметровой бумаге в определенном масштабе проставлялись точки исходных статистических данных, а затем с помощью прозрачной линейки выставлялась некая прямая линия, наилучшим (оптимальным) образом описывающая связь роста человека и его веса. При этом каждый проводящий такую линию мог руководствоваться какими-то своими субъективными понятиями об оптимальности положения такой прямой линии на плоскости, т. е. о значениях констант а и Ь, формирующих уравнение прямой линии у=а+Ъх (константы а и Ь, кстати, рассчитаны и прописаны на рис. 2 и 5). На компьютере эту работу можно автоматизировать, если выбрать какой-то объективный критерий оптимизации. Одним из таких критериев является критерий минимума суммы квадратов «вертикальных» отклонений точек от кривой (в нашем случае от прямой линии). Сам же метод такого расчета называется методом наименьших квадратов: МНК. Если у двух людей возникает спор на тему, чья прямая линия, проведенная вблизи точек, более правильная, то они могут циркулем и линейкой замерить отклонения по вертикали точек от прямой, возвести каждое такое отклонение в квадрат, просуммировать эти квадраты и сравнить полученные суммы. У кого эта сумма окажется меньше, у того и прямая линия правильнее иллюстрирует исследуемую статистическую зависимость. Если у нас есть только две точки (см. рис. 2), или точек больше двух, но все они лежат на одной прямой, то, очевидно, сумма квадратов отклонений будет равна нулю. Но если же не все точки лежат на одной прямой (см. рис. 5), то
можно предположить, что есть такое положение сглаживающей прямой (такие значения коэффициентов а и Ъ), при котором сумма квадратов отклонений будет минимальной. На рис. 6 показан Mathcad-документ, решающий эту задачу оптимизации.
п СКО (а, b) := ^ [Beq - (а + Ь-Рост,)]2 i = 1
а :=-100kg b := 1— cm CKO(a,b) = ЭЗЭ6 kg2
: ^ : := Minimize(CKO,a,b)
а = -79.Э53 kg b = 854 12 ^ cm CKO(a, b) = 1505 kg2
Рисунок 6. Метод наименьших квадратов как задача оптимизации
В Mathcad-документе, показанном на рис. 6, формируется функция пользователя с именем СКО (сумма квадратов отклонений — целевая функция оптимизации), имеющая два аргумента — а и b. Далее задаются начальные значения этих переменных оптимизации с опорой на известную эмпирическую формулу: рост человека в сантиметрах — это его вес в килограммах плюс сто. Наша задача — проверить эту гипотезу на данной статистической выборке — на 50 студентах. Затем вызывается встроенная в Mathcad функция Minimize, которая начинает менять заданные нами значения аргументов а и b так, чтобы функция СКО приняла минимальное значение. Как это делается — это отдельный разговор. Мы же видим, что при начальных значениях констант а и b функция СКО возвращала значение 3396 kg2, а при оптимальных значениях — 1505 kg2, и верим, что меньшего значения в пределах заданной точности найти не удаться. На рис. 5 можно видеть оптимальные значения а и b и «оптимальную» прямую линию, «пронизывающую» массив точек. К оптимальности этого решения мы еще вернемся ниже.
Задачу, сформулированную на рис. 6, можно решить и аналитически, если вспомнить, что у поверхности в точке минимума частные производные по двум переменным равны нулю. На рис. 7 показано символьное решение системы двух уравнений равенства нулю частных производных функции СКО(а, b) по аргументам а и b.
Рисунок 7. Символьное решение задачи оптимизации — поиска минимума функции двух аргументов
У функции СКО суммируются квадраты «вертикальных»5 отклонений точки от прямой линии. И это понятно: если убрать квадрат, то суммирование приведет к некорректной постановке задачи: будут суммироваться отклонения с разными знаками. А можно ли возведение в квадрат заменить взятием абсолютного значения? Для этого достаточно внести соответствующее изменение в первый оператор расчета, показанного на рис. 6. Но в этом случае функция Minimize будет иметь плохую сходимость и выдавать совершенно разные ответы при разных начальных значениях переменных а и b. Дело в том, что функция абсолютного значения угловатая, если так можно выразиться. А природа не любит острых углов. Кроме того, стоит взглянуть на рис. 8, где сделана неудачная попытка аналитического решения системы двух уравнений — равенство нулю частных производных не от функции СКО(а, b), а от функции СМО(а, b): Сумма Модулей Отклонений.
На сайте http://www.math.wpi.edu/Course_Materials/SAS/lablets/7.3/7.3c/ index.html можно видеть описание и сравнение методов наименьших квадратов и наименьших модулей.
Одна точка на рис. 5 обведена кружочком — это, можно сказать, некий типичный студент, отклонение веса которого от средней линии (стандарта) минимально.
5 Можно попробовать просуммировать квадраты «горизонтальных» отклонений — значений не масс студентов, а их роста. С «перпендикулярными» отклонениями поработать не удастся, так как наши два вектора с данными студентов имеют разные размерности — масса и длина. Но «вертикальные» отклонения тут более уместны — рост человека (аргумент) влияет на его вес (функция), а не наоборот.
п СМО{а,Ь):=
--1 СМо(а,Ь) d а СМО (а, Ь) db solve, " [Ь\ °
Решение не было найдено.
Рисунок 8. Неудачное решение системы двух уравнений
В связи с этим возникает интересное предложение по подведению итогов различных конкурсов красоты. Сейчас в них много субъективности, а значит — обид, слез и даже судебных тяжб. В финалы таких конкурсов обычно выходят «красавицы-раскрасавицы» (см., например, рис. 1), из которых довольно трудно выбрать самую оптимальную, пардон, самую красивую мисс или миссис. Так вот, можно у этих финалисток замерить вес и рост или другие размеры (пресловутые 90-60-90, например), провести через точки прямую и назвать победительницами тех, чьи данные попадут на построенную линию и могут быть обведены кружочком как на рис. 5.
Но, вернемся к нашим более серьезным делам. Встроенная в Mathcad функция Minimize, показанная на рис. 6, позволяет нам реализовать метод наименьших квадратов в его общей постановке: пользовательская функция СКО может иметь любое число аргументов, а ее правая часть — любую функциональную зависимость для сглаживания точек. Но для частных случаев в среде Mathcad есть более простые инструменты решения задачи статобработки данных методом наименьших квадратов. Так, на рис. 9 показан вызов двух встроенных в Mathcad функций intercept и slope, возвращающих значение констант а и b сглаживающей прямой y=a+bx. В эти функции заложены формулы, выведенные нами на рис. 7. Они также использованы в расчете на рис. 9. Кроме того, на рис. 9 показана работа еще одной функции, возвращающей коэффициенты линейной регрессии, — функции line. Английские слова an intercept и a slope можно перевести как «пересечение» и «наклон»: значение константы а — это значение абсциссы при пересечении нашей прямой оси 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 точек и как добавляются новые константы при изменении степени полинома от 1 до 14.
При приближении значения степени регрессионного полинома к значению числа точек (рис. 11) наблюдается так называемая осцилляция (лат. oscillatio качание, раскачивание): кривая регрессии приближается к точкам (при n = 14 аппроксимация переходит в интерполяцию), но при этом в промежутках между точками
кривая может «раскачиваться» поляцию и аппроксимацию.
далеко уходить от точек, сводя на нет всю интер-
Рисунок 10. Кадр-анимация работы функции regress: через 15 точек проводится линия полинома 4-й степени
Рисунок 11. Кадр-анимация работы функции regress: через 15 точек проводится линия (часть ее обрезана) полинома 12-й степени
Между двумя правыми точками на рис. 11 из-за осцилляции мы видим не «статистику» и не просто «ложь», а «наглую ложь» (см. второе название статьи).
Вернемся к терминам, описанными нами в самом начале статьи — к корреляции и регрессии. Взаимозависимость роста и веса человека — это типичный пример корреляции, а наш поиск значений констант а и b линейной функции — типичный пример регрессионного анализа.
В традиционной математике, как правило, мы имеем дело не с корреляцией, а с четкими функциональными зависимостями. Каждому значению плоского угла, например, соответствует строго определенное значение синуса, а каждому значению радиуса окружности — значение ее длины и т. д. Но в реальной жизни все намного сложнее. Как зависит теплопроводность металла, о которой мы упоминали выше, от температуры?! С одной стороны — это четкая функциональная зависимость, протабулированная6 во многих справочниках по свойствам веществ. Но с другой стороны, на это физическое свойство металла могут существенно влиять и некие неучтенные факторы — незначительные примеси, форма образца при замере теплопроводности, возраст металла и т. д., и т. п. На сайте http://twt.mpei.ac.ru/MCS/ Worksheets/Thermal/T-T-2-Tab-3 -04.xmcd можно провести аппроксимацию табличных данных по теплопроводности разных металлов с помощью полинома переменной степени и вычислить искомое значение. Плотность раствора соли также может неоднозначно зависеть от температуры и концентрации: соль может иметь примеси, вода, в которой ее растворили, может быть тоже не совсем чистой и т. д. На сайте http://twt.mpei.ac.ru/TTHB/1/Water-chem/C-NaCl.html расположены онлайн расчеты плотности водного раствора соли в зависимости от температуры и концентрации. Тут уже статистическая обработка ведется не на плоскости, а в объеме: у соответствующей функции не один, а два аргумента. Технологии получения подобных зависимостей по табличным и прочим данным посвящена первая глава «Свойства рабочих тел и теплоносителей для теплотехнических расчетов» книги [7].
Видя в книге или Интернете таблицу со значениями, не так просто определить, что это — функция или корреляция, и что применить к этим данным — интерполяцию (проведение линии строго через точки — см. рис. 1) или аппроксимацию (сглаживание, проведение линии вблизи точек — см. рис. 5 и 10).
Мы уже упоминали линейную интерполяцию — см. рис. 2. В среде Mathcad есть встроенная функция linterp(X, Y, х), которая формирует функцию пользователя с аргументом х для интерполяции отрезками прямой линии точек с координатами, хранящимися в вектора X и Y. При такой линейно-кусочной интерполяции
6 В справочниках, по традиции, идущей еще от докомпьютерных времен, помещаются не формулы для компьютерного расчета, а именно таблицы для ручных расчетов. В современных интернетовских справочниках, которые создают в том числе и авторы этой статьи, применяется компромиссный вариант: дается таблица, по которой интерполяцией находится нужное промежуточное значение.
необходимо чтобы эти векторы были одного размера, хранили вещественные числа и в векторе X элементы располагались в порядке возрастания. На авторском сайте http://communities.ptc.com/videos/1512 показана анимация работы функции ¡Метр — поведение прямой линии через пару очередных точек при линейной интерполяции данных.
Рисунок 12. Кадр-анимация линейной интерполяции
Линейная интерполяция очень проста в реализации. Но «иная простота хуже воровства»: при такой упрощенной интерполяции функция получается угловатой, а природа, повторяем, не терпит острых углов и старается их сгладить. Вспомним, как морские волны превращают острые обломки камней в гладкую гальку. Но главное то, что «ломаная» функция, состоящая из отрезков прямой линии, имеет в качестве первой производной набор констант. Вторая же производная такой функции вообще равна нулю (украдена!) во всех точках, кроме «угловых», в которых вторая производная просто не существует. Это часто затрудняет применение к такой функции некоторых численных методов, например, поиск нуля или минимума — см. пример на рис. 8, где была задействована «угловатая» функция.
Если через две точки можно провести прямую линию (см. рис. 2 и 12), то через четыре точки — кубическую параболу. На этом основана кубическая сплайн-интерполяция табличных зависимостей. Она более сложна в компьютерной реализации, чем линейная интерполяция, но лишена «угловатых» недостатков, отмеченных выше. На сайте http://communities.ptc.com/videos/1418 показана авторская анимация интерполяции кубическими сплайнами. Один кадр этой анимации отображен на рис. 13.
А •
Рисунок 13. Анимация сплайн-интерполяции
В среде Mathcad сплайн-интерполяцию можно создать с помощью функций ¡spline, pspline и cspline. Префиксы l, p и с отмечают, по какой зависимости будет вестись экстраполяция (нахождение значений воссозданной функции вне интервала дискретных значений вектора X): l — линейно с опорой на две крайние точки; р — по параболе с опорой на три крайние точки и с — по кубической зависимости с опорой на четыре крайние точки.
Сплайн-интерполяция широко применяется в компьютерной графике, когда, например, необходимо через фиксированные точки какого-либо объекта провести гладкую линию или поверхность. Есть такая технология создания компьютерных мультфильмов. На человеке крепят светящиеся точки, он двигается и эти точки компьютером фиксируются в пространстве. Затем компьютер проводит через эти точки криволинейную поверхность или кусочки поверхности и создает... виртуального движущегося человека.
Экстраполяцию по заданным точкам в среде Mathcad можно реализовать и с помощью встроенной функции Predict (предсказание). На рис. 14 показано, как эта функция предсказывает значения функции пользователя y(x): на первых десяти точках значение заданной функции довольно хорошо предсказывается, а потом происходит срыв. Пусть читатель попробует сам с помощью сайта, показанного на рис. 14, предсказать температуру воздуха или курс доллара на ближайшую неделю, опираясь на данные прошедшей недели, месяца, года и т. д. Затем по прошествии этой «будущей недели» сравнит прогноз с реальными наступившими событиями.
Рисунок 14. Предсказание в среде Mathcad
В среде Mathcad есть много и других функций, позволяющих обрабатывать дискретные данные — сглаживать их, вести по ним интерполяцию, аппроксимацию или экстраполяцию. На рис. 15 показан сайт Интернета, созданный с помощью пакета Mathcad, на котором проводится аппроксимация табличных данных с помощью произвольной функции трех аргументов. На рисунке показана работа так называемой логистической7 функции, с помощью которой можно описать многие явления живой или неживой природы: возникновение, развитие и умирание. Численность паровозов в мире сначала увеличивалось медленно, наблюдался бурный рост их производства, а затем их вытеснили тепловозы и электровозы...
Логистическое уравнение, также известное, как уравнение Ферхюльста (по имени впервые получившего его бельгийского математика), изначально появилось при рассмотрении модели роста численности населения.
Рисунок 15. Сглаживание табличной зависимости по произвольной формуле
В текстовое окно с меткой f(x, a, b, c), показанное на рис. 15, можно вставлять аппроксимирующие функции другого вида. В сам же пакет Mathcad встроены следующие функции с постфиксом fit (fitting), предназначенные для такой работы: expfit(vx, vy, [vg]8) — аппроксимация по формуле (зависимости) a ebb +c;
lgsfit(vx, vy, vg) — l+be~cx; ¡nfit(vx, vy) — a ln(x) + b; ¡ogfit(vx, vy, vg) — a ln(x+b)+c;
medfit(vx, vy) — a + bx (альтернативный алгоритм, реализующий не минимизацию суммы квадратов отклонений (см. рис. 6), а медиан-медианную линейную регрессию для расчета коэффициентов а и b (рис. 16);
pwrfit(vx, vy, vg) — a-xb + c; sinfit(vx, vy, vg) — a-sin(x + b) + c;
¡infit(vx, vy, F) — линейная (lin) комбинация элементарных функций; genfit(vx, vy, vg, F) — регрессионный анализ в общем (gen-general) виде. Из этого списка «фитинговых» функций несколько выпадает функция medfit, в имени которой сокращение med означает не вид аппроксимирующей функции, как у других функций в списке, а метод линейной аппроксимации. Функция medfit как
8 Аргумент vg — это вектор первых предположений при решении данной задачи регрессии. Векторы vx и vy хранят координаты узлов анализируемой выборки.
и функция line (см. рис. 9) предназначена для решения задачи линейной регрессии. Но в функцию medfit заложен не метод наименьших квадратов, а метод медиан-медианной регрессии с минимизацией суммы абсолютных значений медиан ошибок в узлах9. На рис. 16 показано сравнение регрессионного анализа параметров наших студентов, проведенного методом наименьших квадратов и медиан-медианным методом.
Вес
Рисунок 16. Сравнение результатов обработки данных медиан-медианной регрессией (использована функция medfit) и методом наименьших квадратов (line)
Какая прямая на рис. 16 ближе к истине? Или спросим несколько иначе — в какой прямой меньше лжи (см. второе название статьи)? Ответ на этот вопрос10, если, конечно, он вообще существует, зависит, в частности, от того, какое распределение заложено в анализируемые случайные величины. В среде Mathcad, кстати говоря, есть набор встроенных функций для генерации случайных величин по разным законам распределения.
Но нередки случаи, когда не наблюдается никакой корреляции данных, например, так, как это выявлено при исследовании статистической зависимости между баллами ЕГЭ по математике и результатами контрольной работы [1]. Диаграмма рассеяния (рис. 17) ясно дает понять, что остаточный уровень школьных математических знаний слабо зависит от полученных на ЕГЭ баллов, причем дальнейший регрессионный анализ данных не имеет смысла.
9Мы здесь не приводим описание этого метода и отсылает читателя к справочникам по статистике или к Интернету.
10 Наша идея объективного подведения итогов конкурса красоты (см. выше) тут может «споткнуться» о споры о том, какой метод регрессионного анализа тут применять. Можно предложить победительницу конкурса определять одним методом (по одной прямой), а вице-победительницу — другим.
Баллы ЕГЭ по математике ♦ * *
« ; * ! ♦ и ф ♦♦♦ ♦ ♦ ♦♦ ♦ i;:- « ^ » !
j * * ♦ ♦♦ ♦ ♦ * » * * * ♦
Баллы за контрольную работу
_L/U/I/IUI J и I I [^U/l UI I у ^/ЦЦЦ I у_
О 20 40 60 80 100
Рисунок 17. Диаграмма рассеяния между баллами за контрольную работу и баллами ЕГЭ по математике
С другой стороны, предварительный анализ данных может помочь выявить и слабые корреляции. Так, например, была обнаружена, хотя и незначительная, но отрицательная корреляционная связь между временем, проводимым студентами за компьютером в целях развлечений, и успешностью сдачи ЕГЭ по математике [2]. Для студентов из небольших городов она равна -0.04, а для студентов из крупных городов она равна -0.25.
Мы видим, что, приступая к статистической обработке случайных величин, желательно предварительно построить графики и диаграммы для визуального анализа данных. Так, встроенная в Mathcad функция histogram позволяет построить частотную характеристику (гистограмму) выборки. На рис. 18 графически показана частота появления того или иного значения роста у нашей полусотни студентов с точностью в один сантиметр.
■ 1
170 180 190 Рост {cm)
Рисунок 18. Пример гистограммы в среде Mathcad
Так, из гистограммы на рис. 18 видно, что больше всего студентов (10 человек) имеют рост 179 см. Но какую-то статистическую закономерность на рис. 18 увидеть довольно сложно. Если бы наша выборка студентов была бы намного больше (студенты всего вуза или всей страны), то следует ожидать, что на гистограмме четко прорисовалась бы некая колоколообразная кривая — кривая нормального распределения: очень мало людей имеют очень низкий или очень высокий рост, а большое количество людей имеют некий средний (среднестатистический) рост.
Из рис. 18 видно, что наша гистограмма немного сдвинута относительно нормального распределения. Но вспомним еще раз второе название нашей статьи: «Есть ложь, наглая ложь и... статистика». В этом изречении скрыт и такой нюанс. Люди при опросах часто вольно или невольно слегка искажают информацию о себе. Можно предположить, что некоторые студенты, называя вслух свой рост, прибавляли к нему 1, 2 или даже 3 см. Скажем несколько мягче. Низкие студенты округляли свой рост до сантиметра в большую сторону (было, например, 156.3 см — стало 157 см), а высокие проводили округление более правильно (было, например, 179.3 см — стало 179 см). В таких опросах (при переписи населения, например) часто лукавят женщины: убавляют себе возраст и вес, приписывают несуществующего мужа и т. д.
В интернете можно найти не очень качественную, но очень интересную фотографию «живой» гистограммы роста людей — рис. 17.
Рисунок 19. Живая гистограмма роста человека
Группа игроков выстроилась на футбольном поле по росту, но не в одну шеренгу, а в колонну с разным числом людей в каждом ряду, который соответствовал определенному росту: 4 фута и 9 дюймов (144.78 см — этот самый левый ряд на рис. 17 пуст), пять футов (152.4 см) и так до 6 футов и 5 дюймов (195.58 см). В таком «опросе» схитрить при ответе о собственном росте будет довольно сложно: завышение этого параметра на 2-3 см будет сразу заметно в соответствующем ряду, если только не встать на цыпочки. В диаграмме, показанной на рис. 17, нет пу-
стых рядов, которые мы наблюдали на рис. 16, потому, что шаги на этих гистограммах разные — один сантиметр и один дюйм (2.52 см).
Если в США и Англии рост человека измеряют в фунтах и дюймах11, то в старой России для этого брали аршины и вершки12, но аршины (а конкретно два аршина) по умолчанию не указывали. Дело в том, что рост взрослого человека очень редко был меньше двух аршин (примерно 142 см) или больше сажени (трех аршин — примерно 213 см). Читаем (и запомним — см. ниже) у Тургенева в повести «Муму»: «Из числа всей ее челяди самым замечательным лицом был дворник Герасим, мужчина двенадцати вершков роста, сложенный богатырем и глухонемой от рожденья». Несложно подсчитать «богатырский» рост Герасима — примерно 196 см. На рис. 20 построена гистограмма роста «авторских» студентов (см. рис. 5) с шагом не в один сантиметр или дюйм, а в один «старый добрый» русский вершок.
Рисунок 20. «Вершковая» гистограмма роста студентов
На «вершковой» гистограмме (рис. 20) четко просматривается некий «статистический колокол» — кривая нормального распределения. Студенты по росту разбиваются на 5 условных групп: низкого роста (4 человека с шестивершковым ростом), скорее низкого, чем среднего роста (9 человек — 7 вершков), среднего роста
11 В одном из романов Агаты Кристи можно прочесть, что рост Эркуля Пуаро равнялся пяти фунтам и пяти дюймах. Такой же рост имел предок Форсайтов, описанный в знаменитой «Саге» Джона Голсу-орси. Знание старинных и иностранных единиц измерения помогает почувствовать «дух» классической литературы.
12Размер сажени был переменной величиной, но к середине XIX века русскую сажень приравняли к семи английским футам (Семь футов под килем!). В сажени 3 аршина, а в аршине 16 вершков (недесятичная система исчислении). «Прелести» такой системы можно сполна почувствовать, читая книгу Сергея Александровича Рачинского «1001 задача для умственного счета». В печатном виде она издавалась в XIX веке. В наше время она возродилась в виде сайта www.1001task.ru с удобными средствами поиска и комментирования задач.
(23 человека — 8 вершков), скорее высокого, чем среднего роста (9 человек — 9 вершков) и высокого роста (2 человека, 10 вершков). Один студент выпал из нашей «статистики» и оказался «герасимовского» роста в 12 вершков: исключение только подтверждает наше статистическое правило. Сравнивая рис. 18 и 20, можно сказать словами персонажа одного мультфильма: «А в вершках-то наша гистограмма красивее и намного «статистичнее»!». Кстати, здесь мы невольно затронули важнейший для обработки данных вопрос об установлении оптимального шага при работе с выборками. Ведь можно было бы группировать данные каким-то иным способом. Например, по признаку попадания значения роста на интервалы длиной в 1 метр, начиная с нуля (тогда гистограмма имела бы явно ложный вид прямоугольника).
На связь роста и веса человека влияет много факторов: на то она и корреляция, а не функция. Один из важных факторов — это возраст человека. Немногим людям удается сохранить свой вес таким, каким он был в молодости13. На вес человека также большое влияние оказывает его генетика (наследственность) и, конечно, образ жизни. У человека, как и у автомобиля (см. ниже), возраст определить несложно. А вот «пробег» человека, конкретно говоря, число шагов, которые он сделал за свою жизнь, тоже хорошо было бы измерять для более точного определения взаимозависимости веса и роста. Пробег замеряет спидометр автомобиля. Некоторые люди на прогулках или пробежках используют шагомер. Если читатель перейдет по ссылке, отмеченный в списке литературы [6], то он может узнать, как в среде Mathcad проводился более сложный регрессионный анализ — как определялась формула, по которой можно оценить стоимость подержанного автомобиля в зависимости от его возраста и пробега (см. таблицу).
Один из авторов статьи в этой связи вспоминает одну история, свидетелем которой он был в Западной Германии в 80-х гг. прошлого века. Один немецкий профессор (автор был на стажировке в одном немецком университете) решил быстро продать свой автомобиль за 2000 марок и выставил его по этой низкой цене на местном автомобильном «блошином» рынке. К нему довольно долго подходили потенциальные покупатели, осматривали машину, задавали вопросы, но все эти контакты оказывались безрезультатными. Но тут к нему подходит один довольно бойкий человек, покупает, как говорится, не глядя машину, меняет на ней... ценник и тут же на глазах изумленного профессора перепродает ее за 4000 марок... Этот перекупщик умел сходу правильно оценить стоимость подержанного автомобиля и понимал, что заниженная цена также плоха для торга, как и завышенная. Проеци-
13 «Когда хорошее настроение, что-то удалось, нет посторонних и помещение подходящее, хочется раскинуть руки, загудеть и побежать по кругу, изображая самолет на вираже... Давно я не совершал такой «полет», а тут «совершил». Жена посмотрела и сказала, что раньше, в молодости я был легкомоторным самолетом, а теперь стал. широкофюзеляжным.» (www.anekdot.ru).
руя этот тезис на наших студентов, можно сказать, что слишком низкий рост также плох для жизни, как и слишком высокий...
Таблица. Данные о стоимости подержанного автомобиля определенной марки
Возраст Пробег Цена Возраст Пробег Цена Возраст Пробег Цена
(лет) (миль) (1Ю$) (лет) (миль) (1Ю$) (лет) (миль) <ир$)
Аде МНеАде Рпсе Аде МНеАде Рпсе Аде МНеАде Рпсе
11,5 88000 1195 13,5 120000 495 7,5 41000 3495
13,5 103000 750 11,5 101000 895 13,5 124000 850
10,5 82000 1295 3,5 39000 4995 10,5 77000 1595
10,5 65000 1495 10,5 78000 1295 6,5 56000 3495
12,5 97000 800 6,5 52000 2695 12,5 83000 895
10,5 70000 1495 9,5 84000 1995 9,5 67000 2495
8,5 51000 2295 4,5 39000 3995 4,5 38000 3990
10,5 80000 1495 4,5 46000 3675 6,5 43000 3400
9,5 79000 1995 12,5 92000 795 13,5 92000 795
6,5 57000 2695 11,5 108000 975 11,5 78000 1295
На рис. 21-23 можно видеть графики, построенные по точкам из таблицы.
Из рис. 21 видно, что корреляция пробега автомобиля от его возраста (некий эксплуатационный клин) — это область, ограниченная двумя условными лучами: минимальная и максимальная средние скорости автомобилей. В нашей выборке они такие: 1.004 км/ч (7.5 лет и 41 000 миль пробега) и 2.046 км/ч (3.5 года и 39 000 миль пробега). Отсюда вывод — если вам предлагают купить подержанный автомобиль, параметры которого не попадают в наш «клин» (1-2 км/ч), то это значит, что машину либо слишком интенсивно эксплуатировали (скорость больше 2 км/ч), либо она непонятно почему простаивала (скорость меньше 1 км/час), либо у нее. подкручен спидометр. Возраст машины скрыть намного трудней.
Цена автомобиля от возраста зависит линейно (а + Ь х, см. рис. 2 и 22). Никакой другой более сложной зависимости из точек, разбросанных на графике рис. 22, выудить невозможно. В нашем случае после каждого года пробега с цены машины нужно будет скидывать в среднем по 393 доллара.
210 180 150 120 90 60 30
5 / 8
о о о ¡Б 39000 гт _ 2 о 3.5 уг / о. / / о о ° о
6 крЬ / / о / 0 о § о 0 ^
/ / О о ^^ о ^^
8 о/-
^^ \ 41000 тт х —1.004 к-рЬ 7.5 уг
1 Возраст автомобиля (лет)
10
15
Рисунок 21. Корреляция пробега автомобиля и его возраста
Рисунок 22. Корреляция цены автомобиля и его возраста
Зависимость цены автомобиля от пробега (рис. 23) более сложная — она меняется экспоненциально: (а + вь х). При такой зависимости (Ь < 0) цена машины никогда не упадет до нуля, чего не скажешь о зависимости цены от возраста (рис. 22), когда в районе 14-15 лет автомобиль становиться «бесценен» в двух смыслах — у него нулевая цена (наш смысл) либо автомобиль становится... «бесценным» музейным экспонатом.
Статистическая обработка данных по цене подержанного автомобиля позволила создать Mathcad-сайт для оценки данной корреляции (рис. 24).
6000
50 100 150 200
Рисунок 23. Корреляция цены автомобиля и его пробега
По графику на рис. 24 можно не только получить оценочную стоимость старого автомобиля, но и отсечь заведомо сомнительные варианты, зафиксированные на рис. 21. Mathcad-документ, по которому велись расчеты (Car_Price.mcd — Mathcad 2001 Premium), можно «скачать» по адресу: ftp://twt.mpei.ac.ru/ochkov/Auto.
Один из авторов статьи два раза приобретал автомобиль по схеме treid-in, когда в стоимость новой машины входит стоимость старой, сдаваемой на комиссию. Так вот, оба раза цена старого автомобиля, предлагаемая салоном, отличалась от той, какая была рассчитана с помощью сайта, показанного на рис. 24, не более, чем на 5%.
А можно ли оценить в рублях не автомобиль, а человека по его возрасту и пробегу, если под пробегом иметь в виду его деловые качества, которые приобретаются в результате «пробега»: учебы, работы по специальности, повышения квалификации, самообразования, решения разных жизненных и производственных (офисных) коллизий.?
Кому-то этот вопрос покажется несколько кощунственным — как можно человека измерять рублями, долларами, фунтами? Но тут можно упомянуть футболистов, многие из которых бегают по полю с «ценниками», где указаны суммы, какие нужно уплатить при переходе футболиста в другой клуб. У специалистов высшей квалификации цена не прописана так явно, но она негласно или гласно существует и фигурирует в переговорах по переходу, например, менеджера из одной фирмы в другую. Если же под стоимостью человека подразумевать его капитал, то тут можно наблюдать разные кривые в координатах «стоимость — возраст» (см. рис. 22).
Одни люди начинают с нуля, а потом создают себе капитал. Другие же растрачивают средства, полученные по наследству, «скатываясь по наклону». Оптимальным же или самым «красивым» считается вариант, когда человек получает какой-
то стартовый капитал и приумножает его, оставляя затем наследникам какую-то разумную сумму и потратив львиную долю на благотворительность, основав, например, новый университет или поддержав свою Alma Mater, свой родной город.
Рисунок 24. Сайт по оценке стоимости подержанного автомобиля
Но если говорить о цене человека, имея в виду не его капитал, а его ценность как специалиста, то тут тоже можно видеть некую горбатую кривую зависимости этого основного параметра человека от его возраста. Пик (максимум) на этой кривой может приходиться на разный возраст у разных специалистов-профессионалов. У упоминавшихся футболистов это 25-30 лет, у руководителя компании 40-60 и т. д.
Но вернемся к гистограммам, показанным на рис. 18-20. Давайте проведем глобальный (мысленный) статистический эксперимент и у всех взрослых людей на планете измерим какой-либо параметр: вес, рост (см. выше), ум и т. д. — словом все то, что можно измерить числом или оценить лингвистическими критериями (гений, талант, очень умный, просто умный... совсем, пардон, дурак, если говорить об умственных способностях человека).
Полученные точки превратим в кривые, где по оси X отложим параметр человека, а по оси У — процент людей с данным параметром. При этом статобработку проведем отдельно для мужчин и женщин. Что мы получим?
Кривые случая 1 (рис. 25) получаются для тех параметров человека, значение которых у мужчин больше, чем у женщин (рост, вес, сила мышц и т. д.). Это связа-
но с эволюцией — если мужская особь крупнее женской, то новые поколения укрупняются. У пауков, например, самцы намного меньше самок — пауки сейчас, слава богу, более мелкие, чем в доисторические времена.
Рост ч&повека
Рисунок 25. Гистограмма роста мужчин и женщин
А так (случай 2, рис. 26) могут выглядеть кривые для параметров, значения которых у человека за последние несколько тысяч лет не менялись — ум, например. Многие вполне обоснованно полагают, что современный человек, если убрать налет образованности и культуры, ненамного умней древнего грека или египтянина. Среди мужчин гениев (нобелевских и прочих лауреатов, великих изобретателей, знаменитых писателей и художников — см. правый край графиков на рис. 26) больше лишь потому, что и. дураков среди мужчин достаточно (левый край графика). Но средняя женщина умнее среднего мужчины — центр «женской» кривой приподнят за счет меньшего разброса по краям: площади фигур под графиками одинаковы — никого не обижая, будем считать, что Господь Бог или Природа (кто как для себя считает) одинаково наделили умом обе половины человечества. С этим можно, конечно, поспорить, но. еще раз взглянем на второе название статьи.
Ум человека
Рисунок 26. Гистограмма интеллектуальных качеств мужчин и женщин
Выводы. Пакет Mathcad (как, впрочем, и другие пакеты) имеет богатый набор встроенных средств для статистической и прочей обработки массивов данных. Их освоение поможет, в частности, отличить ложь (неверные или несколько искаженные исходные данные) от. статистики.
Литература
[1] Богомолова Е. П., Максимова О. В. Проблемы оценивания результатов ЕГЭ по математике // Alma Mater Вестник высшей школы. 2014. № 9. С. 56-60.
[2] Богомолова Е. П., Максимова О. В. Влияние компьютерной поддержки математики на успеваемость студентов технических вузов // Открытое образование. 2014. № 6. С. 6571.
[3] Кузнецов А. А., Бешенков С. А., Захарова Т. Б., Коротенков Ю. Г., Матвеева Н. В. Об информатике, ее подходах и предмете (философия информатики) // Вестник Тамбовского университета. Серия: Естественные и технические науки. 2005. Т. 10. № 3. С. 236-249.
[4] Кузнецов А. А., Левченко И. В., Заславская О. Ю., Гриншкун В. В. Григорьев С. Г. Содержание обучения информатике в основной школе: на пути к фундаментализации // Вестник Московского городского педагогического университета. Серия: Информатика и информатизация образования. 2010. № 20. С. 6-18.
[5] Коротенков Ю. Г. Обучение информатике и математике в рамках единой образовательной области «Математика и информатика» // Информатика и образование. 2014. № 8(257). С. 36-39.
[6] Очков В. Ф. Цена подержанного автомобиля или Путь от корреляции к регрессии в среде Mathcad // КомпьютерПресс. 2001. № 9. (http://twt.mpei.ac.ru/ochkov/car/car.htm)
[7] Теплотехнические этюды с Excel, Mathcad и Интернет / под общ. ред. В. Ф. Очкова. — СПб. : Издательство БХВ-Петербург, 2014. (http://twt.mpei.ac.ru/ochkov/TTMI)
Авторы.
Очков Валерий Федорович — доктор технических наук, профессор, профессор какфедры тепловых электрических станций Национального исследовательского университета «МЭИ»
Богомолова Елена Петровна — кандидат физико-математических наук, доцент, доцент кафедры Высшей математики Национального исследовательского университета «МЭИ»
Interpolation, extrapolation, fitting or «Lies, damned lies, and statistics»
V. F. Ochkov, E. P. Bogomolova
National Research University Moscow "Power Engineering Institute" 14, Krasnokazarmennay street, Moscow, 111250
e-mail: [email protected]
Abstract. The article describes the basic mathematical tools Mathcad, providing statistical processing data sets and the way from correlation to the regression. The proper use of the program Mathcad for correct solving of problems with different types of data dependencies is discussed.
Keywords: Mathcad, data sets, interpolation, approximation, extrapolation, regression, functional dependence, statistical dependence.
Reference
[1] Bogomolova E. P., Maksimova O. V. (2014) Problemy ocenivanija rezul'tatov EGJe po ma-
tema-tike. Alma Mater Vestnik vysshej shkoly, 9, 56-60. (In Rus)
[2] Bogomolova E. P., Maksimova O. V. (2014) Vlijanie komp'juternoj podderzhki matematiki na uspevaemost' studentov tehnicheskih vuzov. Otkrytoe obrazovanie, 6, 65-71. (In Rus)
[3] Kuznecov A. A. et al. (2005) Ob informatike, ee podhodah i predmete (filosofija informatiki).
Vestnik Tambovskogo universiteta. Serija: Estestvennye i tehnicheskie nauki, 10(3), 236-249. (In Rus)
[4] Kuznecov A. A., et al. (2010) Soderzhanie obuchenija informatike v osnovnoj shkole: na puti k fundamentalizacii. Vestnik Moskovskogo gorodskogo pedagogicheskogo universiteta. Serija: Informatika i informatizacija obrazovanija, 20, 6-18. (In Rus)
[5] Korotenkov Ju. G. (2014) Obuchenie informatike i matematike v ramkah edinoj obrazovatel'-noj oblasti «Matematika i informatika». Informatika i obrazovanie, 257, 36-39. (In Rus)
[6] Ochkov V. F. (2001) Cena poderzhannogo avtomobilja ili Put' ot korreljacii k regressii v sre-de Mathcad. KompjuterPress, 9. (http://twt.mpei.ac.ru/ochkov/car/car.htm) (In Rus)
[7] Ochkov V. F. (ed.) (2014) Teplotehnicheskie jetjudy s Excel, Mathcad i Internet. BHV-Peterburg. (http://twt.mpei.ac.ru/ochkov/TTMI) (In Rus)
