Построение эмпирических формул методом наименьших квадратов

КУРСОВАЯ РАБОТА

Построение эмпирических формул методом наименьших квадратов

Задание

Требуется исследовать эмпирическую зависимость содержания ионов Cl ^- (с, мг-экв/л) от плотности воды (с, кг/мі) поступающей в скважину вместе с нефтью. Используя МНК, построить и исследовать «специализированное» уравнение. Используя МНК, построить и исследовать линейное уравнение общего вида Вычислить по этим уравнениям прогнозные значения величины содержания ионов Cl^-, если плотность пластовой воды равна с_{прогнозн}. В качестве с_{прогнозн} взять величину, равную с_max-0,1 (с_max — с_min), где с_max и с_min — максимальное и минимальное значения плотности с в таблице исходных данных.

Целью курсовой работы является углубление знаний по информатике, развитие и закрепление навыков работы с табличным процессором МS Ехсе1, и применение их для решения задач из предметной области, связанной с исследованиями. В задании формулируются условия задачи, исходные данные, форма выдачи результатов, указываются основные математические зависимости для решения задачи. В соответствии с методом решения задачи разрабатывается алгоритм решения. Поставленные задачи анализа данных можно решить с помощью инструментария MS Excel. Контрольный расчёт в среде МаthCAD позволяет убедиться в правильности решения задач.

Основными этапами курсовой работы являются:

1) формализация поставленной задачи;

2) выбор, обоснование и изложение метода решения задачи;

3) решение задачи с помощью инструментария;

4) оформление пояснительной записки и защита отчёта.

Оформление пояснительной записки должно начинаться с титульного листа, содержать теоретические сведения и описание средств решения задач, постановку задачи.

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

Построение эмпирических формул методом наименьших квадратов

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

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

Таблица 1

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

Если между величинами х и у существует некоторая функциональная зависимость, но её аналитический вид неизвестен, то возникает практическая задача — найти эмпирическую формулу

y^T =F (x, a₁, a₂,…, a_m), (2)

где а₁, а₂,…, а_т — коэффициенты. Вид функции и значения коэффициентов а₁, а₂,…, а_т подбираются таким образом, чтобы значения = F(x_i, a_l, a₂,…, a_m), вычисленные по эмпирической формуле при различных значениях x_i, как можно меньше отличались бы от опытных значений у_i.

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

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

Для аппроксимации вначале определяют класс функций, из которых выбирается аппроксимирующая функция F (x, a_l, a₁,…, a_m), и далее отыскивают наилучшие значения коэффициентов.

Чаще всего для аппроксимации используют метод наименьших квадратов (МНК). Поясним геометрический смысл этого метода. Каждая пара чисел (x_i, y_i) из исходной таблицы определяет точку М_i на плоскости XOY. Используя формулу (2) с различными значениями коэффициентов а₁, а₂,…, а_т, можно построить множество кривых, которые будут являться графиками теоретических функций F (x, a_l, a₂,…, a_m). Величина = F(x_i, a₁, а₂,…, а_т) называется теоретическим значением функции в точке х_i. Разность (- у_i) называется отклонением или остатком и представляет собой расстояние по вертикали от точки М, до графика эмпирической функции.

Рис. 1. Геометрический смысл метода наименьших квадратов Согласно методу наименьших квадратов, наилучшими коэффициентами a1, a2,…, a_m считаются те, для которых сумма квадратов отклонений найденных теоретических значений функции от заданных эмпирических значений будет минимальной. Следовательно, задача состоит в определении коэффициентов а₁, а₂,…, а_т таким образом, чтобы сумма квадратов отклонений была наименьшей.

S(a_l, a₂,…, a_m) = (3)

Построение эмпирических формул состоит из двух этапов: выяснение общего вида этой формулы и определение её наилучших параметров.

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

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

Для того, чтобы найти набор коэффициентов а₁, а₂,…, а_т, при которых достигается минимум функции S(a₁, a₂,…, a_m), определяемой формулой (3), используем необходимое условие экстремума функции нескольких переменных — равенство нулю частных производных. В результате получим нормальную систему дли определения коэффициентов a₁, а₂,…, а_т:

(4)

Таким образом, нахождение коэффициентов a₁, а₂,…, а_т сводится к решению системы (4). Эта система упрощается, если эмпирическая формула (2) линейна относительно параметров a₁, а₂,…, а_т, тогда система (4) будет линейной.

Конкретный вид системы (4) зависит от того, из какого класса эмпирических формул мы ищем зависимость (2). В случае линейной зависимости у = а₁ +а₂х система (4) примет вид:

(5)

где a₁ и а₂ — неизвестные, а суммы (); () и т. д. дают конкретные значения коэффициентов и свободных членов в системе линейных уравнений (5). Эта линейная система может быть решена любым известным методом (с помощью обратной матрицы, методом Гаусса, простых итераций, по формулам Крамера, и т. д.). В случае квадратичной зависимости у = а₁+ а₂х + а₃х² система (4) примет вид

(6)

Линеаризация нелинейных зависимостей

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

у = а₁х^а2. (7)

экспоненциальная зависимость:

(8)

и показательная зависимость:

(9)

В приведенных выше зависимостях а₁ и а₂ являются коэффициентами, которые необходимо определить численно.

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

В случае степенной зависимости линеаризацию выполним путём логарифмирования уравнения (7). В результате чего получим соотношение:

(10)

Обозначим ln y, ln х и ln a₁, соответственно через z, t и b, тогда зависимость (9) может быть записана в виде z = b + a₂t, что позволяет применить формулы (5) с заменой a₁, на b и пересчетом исходных данных z_i = ln y_i, а t_i = ln x_i. После вычисления b определяем значение коэффициента a₁ исходной зависимости по формуле а₁= е^b.

Линеаризацию экспоненциальной зависимости выполняем путём логарифмирования равенства (8), после чего получаем соотношение

(11)

Обозначим ln y и ln а₁ соответственно через Z и с, тогда зависимость (6.4) может быть записана в виде z=c+а₂· х, что позволяет применить формулы для вычисления коэффициентов линейной зависимости (с заменой а на с и у_i на z_i).

Линеаризующие преобразования для различных видов функций приведены в таблице 2.

Таблица 2

Специальный вид линейной зависимости

В традиции некоторых разделов науки зависимость содержания различных ионов (минерализации) от плотности пластовой воды принято представлять в виде

Y = с· (Х-1000), (12)

где X — плотность пластовой воды; Y — степень минерализации; с — некоторый параметр, зависящий от типа растворённых ионов (минералов).

Такой вид зависимости имеет ряд преимуществ по сравнению с традиционным. Во-первых, он обеспечивает нулевое значение степени минерализации для дистиллированной воды, которая имеет плотность, равную 1000, и при этом значении X, и, очевидно, Y равны нулю. Во-вторых, параметр с имеет простой смысл — он показывает, на сколько повышается степень минерализации для пластовой воды, если её плотность возрастает на единицу по сравнению с плотностью дистиллированной воды.

В качестве примера приведём эмпирическую зависимость содержания ионов от плотности воды, поступающей в скважину вместе с нефтью. Для этого был проанализирован состав вод по скважинам горизонта Д1 Ромашкинского месторождения, эксплуатирующимся в течение длительного времени. Для основных ионов пластовой воды: Na+, К+, Сa₂+, Mg₂+, Cl- — указанная зависимость в пределах изменения плотности пластовых вод 1030−1185 кг/мі (типичная зависимость) приведена на рис. 2. Очевидно, что она носит линейный характер и хорошо аппроксимируется уравнением

с (с) = с · (с -1000), (13)

где с — плотность воды, поступающей в скважину; с — постоянная величина для данного вида иона, характеризующая концентрацию.

Рис. 2. Зависимости содержания ионов от плотности воды, поступающей в скважину вместе с нефтью

Полученная закономерность подтверждена результатами обработки данных об изменении состава вод при заводнении продуктивных пластов девонского горизонта и верхнего карбона на 121-м месторождении Урало-Поволжья и Западной Сибири.

Чтобы определить коэффициент с в зависимости (12), достаточно сделать замену

Z = (X — 1000) (14)

Тогда зависимость (12) примет вид

Y = c· Z (15)

Чтобы найти коэффициент с с помощью МНК, воспользуемся следующими соотношениями. Цель МНК — найти такое значение с, чтобы сумма квадратов отклонений S(c) была минимальной. Сумма квадратов отклонений в данном случае равна Необходимым условием экстремума является соотношение

(16)

Соотношение (16) может быть преобразовано к виду:

Отсюда можно определить с:

(17)

Элементы теории корреляции

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

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

Коэффициент корреляции вычисляется по формуле:

(18)

где — среднеарифметические значения по x и y соответственно.

Коэффициент корреляции по абсолютной величине не превосходит 1. Чем ближе |r| к 1, тем теснее линейная связь между х и у, и тем более целесообразна аппроксимация таблично заданной функции линейной зависимостью.

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

Чтобы определить, насколько хорошо построенная зависимость отображает эмпирические данные, водится ещё одна характеристика — коэффициент детерминированности R².

Пусть S_ocm — сумма квадратов отклонений теоретических значений функции от эмпирических данных:

(19)

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

Введём понятие регрессионной суммы квадратов:

(20)

Эта величина характеризует разброс теоретических данных относительно среднего значения.

Для линейной зависимости справедливо следующее соотношение:

(21)

Обозначим:, тогда для линейной зависимости справедливо следующее равенство:

(22)

Коэффициент детерминированности R² определяют по формуле:

(23)

Поскольку S_noлн ? 0, S_ocm ?0, и из формулы (23) следует, S_noлн ? S_ocm, то из формулы (15) следует

0?R²?1. (24)

Чем меньше остаточная сумма квадратов S_ocm по сравнению с общей суммой квадратов S_nom, тем больше значение коэффициента детерминированности R². Коэффициент детерминированности R² показывает, насколько хорошо полученная теоретическая функция описывает взаимосвязь между эмпирическими данными. Если этот коэффициент равен 1, то имеет место полное совпадение выбранной теоретической модели с фактическими данными. В противоположном случае, если коэффициент детерминированности близок к нулю, то выбранная эмпирическая формула неудачна, и она не может использоваться для вычисления значений функции.

Коэффициент детерминированности служит показателем тесноты связи между фактором х и откликом у, описываемой данным уравнением.

Иногда показателям тесноты связи можно дать качественную оценку (шкала Чеддока) (табл. 3).

Таблица 3

Решение задачи в MS Excel

Расчеты для определения коэффициентов зависимостей (1) и (13) с использованием MS Excel приведены на рис. 3−15.

На рис. 3−4 в ячейке С35 содержится вычисленный коэффициент корреляции r, он примерно равен 0,993. Это значение близко к единице, что позволяет сделать вывод о том, что содержание ионов Cl — (с, мг-экв/л) и плотность воды (с, кг/мі), поступающей в скважину вместе с нефтью, связаны линейной зависимостью.

На рис. 3−4 также приведены расчеты по определению коэффициента с для «специализированного» уравнения (12).

Коэффициент с для иона Cl-, характеризующий концентрацию, находится в ячейке G34 и равен 25,76. Итак, искомое «специализированное» уравнение имеет вид:

с© = 25,77 · (с — 1000) (25)

Для определения качества уравнения вычислим коэффициент детерминированности R². (рис. 5−6). Значение, равное 0,9864, приведено в ячейке J35. Такая величина (близкая к единице), позволяет сделать вывод, что «специализированное» уравнение хорошо описывает эмпирические данные.

Рис. 3. Решение задачи. Расчеты для «специализированного» уравнения на листе MS Excel в режиме отображения данных (начало) Рис. 4. Решение задачи. Расчеты для «специализированного» уравнения на листе MS Excel в режиме отображения формул (начало) Рис. 5. Решение задачи. Расчеты для «специализированного» уравнения на листе MS Excel в режиме отображения данных (окончание) Рис. 6. Решение задачи. Расчеты для «специализированного» уравнения на листе MS Excel в режиме отображения формул (окончание) Для определения коэффициентов зависимости (1), т. е. «общего» уравнения с использованием MS Excel приведены на рис. 7−10. При этом коэффициенты системы (5) определены в интервале ячеек B66: D66 (рис. 7−8). Система нормальных уравнений примет вид:

Коэффициенты матрицы системы содержатся в интервале ячеек B70: C71; вектор правых частей — в интервале ячеек B74: C75. Вектор решения содержится в интервале ячеек E74: E75. Таким образом, а₁ и а₂ равны -25 828 и 25,822 соответственно, и уравнение общего вида может быть записано как с© = 25,822· с — 25 828 (26)

Рис. 7. Решение задачи. Расчеты для «общего» уравнения на листе MS Excel в режиме отображения данных (начало) Рис. 8. Решение задачи. Расчеты для «общего» уравнения на листе MS Excel в режиме отображения формул (начало) Для определения качества «общего» уравнения вычисляем коэффициент детерминированности R² (рис. 9−10). Значение, приведенное в ячейке I68, равно 1. Такая величина (близкая к единице) позволяет сделать вывод, что полученное «общее» уравнение хорошо описывает эмпирические данные и делает это лучше, чем «специализированное» уравнение, у которого = 0,99. Это улучшение невелико, поскольку разница в коэффициентах мала и составляет 0,01.

Рис. 9. Решение задачи. Расчеты для «общего» уравнения на листе MS Excel в режиме отображения данных (окончание) Рис. 10. Решение задачи. Расчеты для «общего» уравнения на листе MS Excel в режиме отображения формул (окончание) Рис. 11. Решение задачи. Расчеты для определения прогнозного значения на листе MS Excel в режиме отображения данных Рис. 12. Решение задачи. Расчеты для определения прогнозного значения на листе MS Excel в режиме отображения данных Расчёты по этим уравнениям прогнозных значений величины содержания ионов Сl-, если плотность пластовой воды равна заданному значению с_{прогнозн}, приведены на рис. 11−12. Для определения величины с_{прогнозн} в ячейку F79 вводим формулу, соответствующую выражению с_тах — 0.1· (с_max — с_min), где с_тax и с_min — максимальное и минимальное значения плотности с в таблице исходных данных (интервал ячеек Е7: ЕЗ1). В результате получим значение с_{прогноз}, равное 1157,8. Для вычисления прогнозных значений величины содержания ионов Сl - подставим это значение в найденные уравнения. Полученные прогнозные значения величин содержания ионов для «специализированного» и «общего» уравнений находятся в ячейках G83 и G84, соответственно. Абсолютная и относительная разности прогнозных значений приведены в ячейках H83 и I83.

Коэффициент с зависимости (13) также может быть найден с помощью MS Excel, если воспользоваться специальной возможностью средства «Диаграмма-Тренд». Для этого на вкладке «Параметры» нужно установить отметку (v) в элементе управления с надписью «пересечение кривой с осью Y в точке 0» (рис. 13).

Рис. 13. Вкладка «Параметры» средства «Диаграмма-Тренд»

Результат применения этого средства приведен на рис. 14. Заметим что по оси абцисс откладываются не истинные значения плотности, а их отклонения от плотности дистиллированной воды, равной 1000. Коэффициент с и коэффициент детерминированности получились равными 25,822 и 0,9864, что полностью совпало со значением, вычисленным по расчетным формулам.

Рис. 14. Решение задачи. Построение средством «Диаграмма-Тренд» зависимости для «специализированного» уравнения на листе MS Excel

Коэффициенты а₁ и а₂ уравнения «общего» вида (1) также определены с помощью MS Excel «Диаграмма-Тренд» на диаграмме рис. 15. Все полученные величины совпали со значениями, полученными по расчётным формулам.

Рис. 15. Решение задачи. Построение средством «Диаграмма-Тренд» зависимости для «специализированного» уравнения на листе MS Excel.

Решение задачи в MathCad

Расчёты для определения коэффициентов зависимости (1) и (13) с использованием Mathcad приведены на рис. 16−20.

Исходные данные запишем в текстовый файл. Для этого можно скопировать исходные данные из интервала ячеек В7: С31 (рис. 4) в буфер обмена и перенести в текстовый файл, созданный с помощью «Блокнота». Далее (при необходимости) установить символ «точка» в качестве разделителя между целой и дробной частью чисел, после чего сохранить полученный файл с именем mnk_oil.txt. Вид текстового файла приведен на рис. 16. Этот текстовый файл будет использован для ввода данных как при расчетах в Mathcad, так и в MATLAB.

Все вычисления в Mathcad достаточно прозрачны. Для ввода данных используется команда > Вставка > Компонента из главного меню Mathcad.

Будет запущен Мастер (специальная программа), который обеспечит выполнение всей процедуры ввода, результат ввода будет находиться в двумерном массиве W. Далее эти данные (по столбцам) присваиваем одномерным массивам x и y, с которыми выполняем все дальнейшие операции.

Рис. 16. Решение задачи. Подготовка данных и запуск Мастера Компонентов для ввода данных в системе MathCad

Рис. 17. Решение задачи 1 в системе Mathcad (начало)

Рис. 18. Решение задачи 1 в системе Mathcad (продолжение)

Рис. 19. Решение задачи 1 в системе Mathcad (продолжение)

Рис. 20. Решение задачи 1 в системе Mathcad (продолжение)

Рис. 21. Решение задачи 1 в системе Mathcad (продолжение)

Рис. 22. Решение задачи 1 в системе Mathcad (окончание)

Основные выводы по решению задачи

Все вычисленные величины и коэффициенты зависимостей в системе Mathcad полностью совпали с соответствующими значениями, вычисленными в MS Excel.

Эмпирические зависимости содержания ионов Cl ^- (с, мгэкв/л) от плотности воды с (кг/мі) могут быть записаны в виде:

— «специальное» уравнение;

— уравнение «общего» вида.

Сопоставим величины коэффициентов детерминации и. Можно заметить, что, а значит уравнение «общего» вида формально хуже описывает экспериментальные данные, чем «специализированное» уравнение. Стоит заметить, что отличие незначительно, поскольку объяснённая доля дисперсии величины «с» «общим» уравнением всего на 0,01% меньше объяснённой доли дисперсии величины «с» «специализированным» уравнением. Таким образом, принимая во внимание, что весьма близко к единице и это уравнение умеет простой и ясный физический смысл, приходим к выводу, что для описания зависимости содержания ионов Cl ^- от плотности воды, поступающей в скважину вместе с нефтью, можно использовать «специализированное» уравнение. Другим выводом в пользу использования «специализированного» уравнения является тот факт, что графики этих зависимостей, приведённые на рис. 21, отличаются весьма незначительно.

Прогнозные значения содержания ионов Cl^-, вычисленные по этим уравнениям при равны 4066,5 и 4068,1, соответственно. Отклонение этих величин незначительно и равно 0,04%. Полученные прогнозные значения величины для «специализированного» и «общего» уравнений имеют расхождение в единичных знаках. Это объясняется различной настройкой точности вычислений Worksheet Options в MS Excel и в MathCAD.

1. Мановян А. К. Технология первичной переработки нефти и природного газа. М.: Химия. 2001, 568 с.

2. Очков В. Ф. MathCAD 14 для студентов и инженеров: русская версия. СПб: BHV — Петербург. 2009, 512 с.

3. Ломтадзе В. Д. Словарь по инженерной геологии. СПб.: СПГГИ (ТУ), 1999, 424 с.