Чичкарёв Е.А. Компьютерная математика с Maxima

Подождите немного. Документ загружается.

8.3. Пример построения статистической модели 211

-10 -5 0 5 10 15 20 25

Рис. 8.2. Траектория броуновского движения модельной частицы

8.2.3 Броуновское движение

Наличие генератора случайных чисел дает возможность моделировать движение броуновской

частицы.

Эйнштейн первый рассчитал параметры броуновского движения, показав, что нерегулярное пе-

ремещение частиц, взвешенных в жидкости, вызвано случайными ударами соседних молекул, совер-

шающих тепловое движение. В соответствии с теорией Смолуховского-Эйнштейна, среднее значение

квадрата смещения броуновской частицы (s

) за время t прямо пропорционально температуре Т и

обратно пропорционально вязкости жидкости h, размеру частицы r и постоянной Авогадро N

2RT t

6πhrN

где R — газовая постоянная. Броуновские частицы имеют размер порядка 0,1–1 мкм, т.е. от одной

тысячной до одной десятитысячной доли миллиметра.

Построим несколько упрощённую модель броуновского движения, предполагая, что смещение

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

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

необходимые функции (использован генератор random_normal).

(%i1) load("distrib");

(%o1) /usr/share/maxima/5.17.1/share/contrib/distrib/distrib.mac

(%i2) x:0$ y:0$ xy:[[0,0]]$ m:0$ s:1$

Nmax:500$ for i:1 thru Nmax do (x:x+random_normal(m,s),

y:y+random_normal(m,s), xy:append(xy,[[x,y]]));

(%o8) done

(%i9) plot2d([discrete,xy]);

Результат построения графика приведен на рис. 8.2.

8.3 Пример построения статистической модели

Рассмотрим построение задачи с практическим содержанием.

В таблице 8.1 п риведены данные (взяты из статьи В.Ф. Очкова - http://twt.mpei.ac.ru/ochkov/)

о зависимости цены подержанного автомобиля от его пробега и "возраста"(времени использования).

В статье-первоисточнике задача исследования этой зависимости решалась средствами MathCad.

Рассмотрим её решение средствами Maxima.

212 Глава 8. Решение физических и математических задач с Maxima

Возраст

(лет)

Пробег

(миль)

Цена ($) Возраст

(лет)

Пробег

(миль)

Цена ($)

11.5 88000 1195 13.5 103000 750

10.5 82000 1295 10.5 65000 1495

12.5 97000 800 10.5 70000 1495

8.5 51000 2295 10.5 80000 1495

9.5 79000 1995 6.5 57000 2695

13.5 120000 495 11.5 101000 895

3.5 39000 4995 10.5 78000 1295

6.5 52000 2695 9.5 84000 1995

4.5 39000 3995 4.5 46000 3675

12.5 92000 795 11.5 108000 975

7.5 41000 3495 13.5 124000 850

10.5 77000 1595 6.5 56000 3495

12.5 83000 895 9.5 67000 2495

4.5 38000 3990 6.5 43000 3400

13.5 92000 795 11.5 78000 1295

Таблица 8.1. Стоимость подержанного автомобиля в зависимости от его возраста и пробега

В дальнейшем предполагается, что исходные данные для решения подготовлены в виде файла

cars.txt. Для считывания используем пакет numericalio. В памяти данные представляются матрицей,

а для построения отдельных графиков - списками (переменные age, mile,price - см. ниже).

(%i1) load("draw");

(%o1) /usr/share/maxima/5.13.0/share/draw/draw.lisp

(%i2) load("numericalio");

(%o2) /usr/share/maxima/5.13.0/share/contrib/numericalio/numericalio.lisp

(%i3) data:read_matrix("cars1.txt")$

(%i4) age:makelist(data[k,1], k, 1, 30)$

(%i5) mile:makelist(data[k,2], k, 1, 30)$

(%i6) price:makelist(data[k,3], k, 1, 30)$

Простейшую линейную регрессию можно построить, используя функцию simple_linear_regression

(пакет stats). Построим зависимость цены автомобиля от его стоимости и пробега:

(%i21) xy:makelist([age[k],price[k]], k, 1, 30)$(%i22) simple_linear_regression(xy);

8.3. Пример построения статистической модели 213

1000

2000

3000

4000

5000

0 2 4 6 8 10 12 14

Price

Age

Table

price=f(age)

Рис. 8.3. Зависимость цены подержанного автомобиля от его возраста

(%o22)







SIMP LELINEARREGRESSION

model = 5757.594446543255 − 392.7181715149224 x

correlation = −0.96881779424672

p_value = 0.0







Фактический вывод результатов несколько сокращён.

Построим аналогичную зависимость цены автомобиля от пробега, но не в линейной, а в экспо-

ненциальной форме:

(%i26) xy:makelist([mile[k],log(price[k])], k, 1, 30)$(%i27) simple_linear_regression(xy);

(%o27)







SIMP LELINEARREGRESSION

model = 9.174960600286802 − 2.3747715120748162 10

−5

correlation = −0.93011255642444

p_value = 5.9285909514983359 10

−14







Полученные зависимости представлены в виде выражений Maxima:

(%i28) fun1:5757.6-392.7*x$(%i29) exp(9.175);

(%o29) 9652.768071616591

(%i30) fun2:9653*exp(-2.375*10^(-5)*x);

(%o30) 9653 e

−2.3750000000000001 10

−5

Проиллюстрируем полученные результаты графически (рис. 8.3 и рис. 8.4 ):

(%i34) draw2d(terminal=eps,key="Table",xlabel="Age",ylabel="Price",point_size = 3,point_type=3,points(age,price),key="price=f(age)",explicit(f

(%i41) draw2d(terminal=eps,key="Table",xlabel="Mile",ylabel="Price",point_size = 3,point_type=3,points(mile,price),key="price=f(mile)",explicit(

214 Глава 8. Решение физических и математических задач с Maxima

1000

2000

3000

4000

5000

6000

7000

8000

9000

0 20000 40000 60000 80000 100000 120000

Price

Mile

Table

price=f(mile)

Рис. 8.4. Зависимость цены подержанного автомобиля от его пробега

20000

40000

60000

80000

100000

120000

140000

1000

2000

3000

4000

5000

6000

7000

Price

Table

price=f(age,mile)

Age

Mile

Price

Рис. 8.5. Иллюстрация зависимости отклика (цены подержанного автомобиля) от двух независимых

факторов (возраста и пробега автомобиля)

Для построения модели в виде зависимости цены автомобиля от пробега и возраста одновре-

менно целесообразно использовать более сложную функцию lsquares_estimates (пакет lsquares).

Искомая модель была представлена уравнением:

P rice = a + b ∗Age + c ∗ Mile + d ∗M ile

)

Необходимые команды Maxima:

(%i5) lsquares_estimates (data, [x,y,z], z=a+b*x+c*y+d*y^2, [a,b,c,d]);

(%o5)

[[a =

36712000090549571

5117101479342

, b = −

80056614985946

284283415519

, c = −

194393701258481

3411400986228000

, d =

2937180994967

10234202958684000000

]]

(%i6) float(%);

(%o6)

[[a = 7174.37405506961, b = −281.6084604857839, c = −0.056983539033745, d = 2.8699655525931528 10

−7

]]

Следует отметить, что сильно нелинейные задачи решаются при помощи lsquares_estimate мед-

ленно, поэтому результаты построения модели сильно зависят от обоснованности постановки задачи

оценивания. ллюстрация Графическая иллюстрация представлена на рис.

Приложения

Таблица 1: Сокращённый список основных функций

Maxima

Функция или

переменная

Краткое описание

’, , %

простейшие команды , см. стр. 17

addcol Функция добавляет столбец к матрице , см. стр.

30-35

addrow Функция добавляет строку к матрице , см. стр.

30-35

algsys Функция решает полиномиальные системы урав-

нений. Допускаются системы из одного уравне-

ния с одной неизвестной. Кроме того, допускают-

ся недоопределенные системы , см. стр. 63-67

allroots Функция, которая находит и печатает все (в том

числе и комплексн ые) корни полиномиального

уравнения с действительными либо комплексны-

ми коэффициентами , см. стр. 63-67

antidiﬀ Функция выполняет интегрирование выражений с

произвольными функциями, перед ее первым вы-

зовом следует загрузить пакет "antid см. стр. 98-

102

append Функция позволяет склеивать два списка , см.

стр. 22

arrayinfo Функция печатает информацию о массиве — его

вид, число индексов, размер , см. стр. 26

arrays Переменная содержит список имен массивов пер-

вого и второго видов, определенных на данный

момент , см. стр. 29

array Функция определяет массив с данным именем,

определенным количеством индексов и заданным

размером , см. стр. 26

assume Функция вводит информацию о переменной в базу

данных , см. стр. 36-40

216 Глава 8. Приложения

atom Функция возвращает "true если аргумент не име-

ет структуры, т.е. составных частей (например,

число или переменная не имеют структуры) , см.

стр. ??

atvalue Функция позволяет задать значение функции и ее

производных при некоторых значениях аргумен-

тов , см. стр. 126

at Функция вычисляет значение выражения в задан-

ной точке с учетом свойства , см. стр. 126

augmented_lagrangian_methodФункция осуществляет минимизацию ФНП с

ограничениями, см. стр. 161

batch Функция запускает файл с программой. Операто-

ры выполняются один за другим либо до конца

файла, либо до синтаксической ошибки, либо до

некорректной операции, см. стр. 156

bс2 Функция позволяет учесть краевые условия в ре-

шениях дифференциальных уравнений второго

порядка , см. стр. 116

cabs Функция возвращает модуль комплексного выра-

жения, см. стр. 56-58

carg Функция возвращает фазу комплексного выраже-

ния, см. стр. 56-58

cfdisrep Функция преобразует список (как правило ре-

зультат выполнения функции cf) в собственно

цепную дробь , см. стр. 114

cf Функция Создает цепную дробь, аппроксимирую-

щую данное выражение. Выражение должно со-

стоять из целых чисел, квадратных корней целых

чисел и знаков арифметических операций. Воз-

вращаемый результат - список , см. стр. 114

cfdirep Функция преобразует список в собственно цепную

дробь, см. стр. 114

changevar реализует замену переменных в интеграле , см.

стр. 98-102

charpoly Функция является до некоторой степени избыточ-

ной — она вычисляет характеристический поли-

ном матрицы (корни этого полинома — собствен-

ные значения матрицы), см. стр. 58- 63

closeﬁle Функция прекращает вывод в файл , см. стр. 154

COl Функция выделяет заданный столбец матрицы ,

см. стр. 30-35

combine Функция объединяет слагаемые с идентичным

знаменателем , см. стр. 36-40

8.3. Пример построения статистической модели 217

compile Функция сначала транслирует функцию Maxima

на язык LI SP, а затем компилирует эту функцию

LISP’a до двоичных кодов и загружает их в па-

мять , см. стр. 152

conjugate Функция для вычисления комплексно-

сопряжённых выражений , см. стр. 56-58

cons Функция позволяет добавлять элемент в начало

списка , см. стр. 23

contrib_ode Функция решает дифференциальные уравнения

(больше возможностей, чем у ode2), см. стр. 128-

131

copylist Функция создаёт копию списка , см. стр. 21

create_list Функция создаёт список , см. стр. 22

copymatrix Функция Создаёт копию матрицы , см. стр. 30-35

cspline Функция строит сплайн-интерполяцию , см. стр.

159

deﬁne Функция позволяет преобразовать выражение в

функцию , см. стр. 46

demoivre Функция заменяет все экспоненты с мнимы-

ми показателями на соответствующие три-

гонометрические функции , см. стр. 56 - 58

denom Функция выделяет знаменатель , см. стр. 40

depends Функция позволяет декларировать, что перемен-

ная зависит от одной или нескольких других пе-

ременных , см. стр. 94

desolve Функция решает дифференциальные уравнения и

системы дифференциальных уравнений методом

преобразования Лапласа, см. стр. 127

determinant Функция вычисляет детерминант матрицы, см.

стр. 58- 63

diﬀ Функция выполняет дифференцирование , см.

стр. 79

display2d Переменная включает или выключает "дву-

мерное"рисование дробей, степеней, и т.п. Из-

начально установлено значение "true см. стр. 153

display Функция печатает значения своих аргументов

вместе с их именем, каждое в отде льной строке,

см. стр. 153

disp Функция печатает значения своих аргументов,

причем каждое значение печатается в отдельной

строке? см. стр. 153

218 Глава 8. Приложения

divide Функция позволяет вычислить частное и остаток

от деления одного многочлена на другой, см. стр.

36 - 40

draw2d строит двумерные графики, см. стр. 170-172

draw3d строит трёхмерные графики , см. стр. 170-172

echelon Функция преобразует матрицу к верхней тре-

угольной, см. стр. 58- 63

eigenvalues Функция аналитически вычисляет собственные

значения матрицы, см. стр. 58- 63

eigenvectors Функция аналитически вычисляет собственные

значения и собственные вектора матрицы, если

это возможно , см. стр. 58- 63

eliminate Функция и сключает из системы уравнений ука-

занные переменные. Оставшиеся уравнения при-

водятся к виду с нулевой правой частью, которая

опускается , см. стр. 67

endcons Функция позволяет добавлять элемент в конец

списка, см. стр. 23

ev Функция является основной функцией, обрабаты-

вающей выражения , см. стр. см. стр. 36-40

expand Функция раскрывает скобки, см. стр. 36-40

exponentialize Функция приводит комплексное выражение к экс-

поненнциальной форме , см. стр. 56-58

express Функция преобразует дифференциальные опера-

торы в выражения, см. стр. 95

factor Функция представляет в виде произведения неко-

торых сомножителей заданное выражение, см.

стр. 36 - 40

factorsum Функция факторизует отдельные слагаемые в вы-

ражении , м. стр. 36 - 40

ﬁllarray Функция позволяет заполнять одноиндексные

массивы третьего вида из списка , см. стр. 29

ﬁnd_root Функция находит корень уравнения на заданном

интервале методом деления отрезка пополам , см.

стр. 156

ﬁrst Функция выделяет первый элемент списка , см.

стр. 23

ﬂoat Функция конвертирует любые числа в выражени-

ях в числа машинной точности , см. стр. 18

fourier Функция позволяет вычислить коэффициенты ря-

да Фурье см. стр. ??

foursimp Функция позволяет упростить коэффициенты ря-

да Фурье см. стр. ??

8.3. Пример построения статистической модели 219

fullratsimp Функция вызывает функцию "ratsimp"до тех пор,

пока выражение не перестанет меняться , 41-43

genmatrix Функция возвращает матрицу заданной размер-

ности, составленную из элементов индексного

массива , см. стр. 30-35

gfactorsum Функция представляет в виде сомножителей сла-

гаемые выражения с комплексными числами, см.

стр. 36-40

gfactor Функция предс тавляет в виде сомножителей вы-

ражение с комплексными числами, см. стр. 36-40

gradef Функция определяет результат дифференцирова-

ния функции по своим аргументам , см. стр. 94

gramschmidt Функция вычисляет ортонормированную систему

векторов, см. стр. 58-63

ic1 Функция п озволяет учесть начальное условие в

решениях дифференциальных уравнений первого

порядка , см. стр. 116

ic2 Функция позволяет учесть начальные условия в

решениях дифференциальных уравнений второго

порядка , см. стр. 116

ident Функция возвращает единичную матрицу задан-

ной размерности , см. стр. 30-35

ilt Функция реализует обратное преобразование Ла-

пласа , см. стр. 102-103

imagpart Функция возвращает действительную часть вы-

ражения , см. стр. 56-58

integrate Функция выполняет интегрирование заданного

выражения по указанной переменной (неопреде-

ленная константа не добавляется). Можно также

указать пределы интегрирования — в этом случае

вычисляется определенный интеграл , см. стр. 98-

102

invert функция выполняет обращение матрицы , см. стр.

30 - 35

join функция выполняет компоновку списков , см. стр.

kill Функция уничтожает всю информацию (как свой-

ства, так и присвоенное значение) об объекте или

нескольких объектах , см. стр. 19

lagrange Функция строит интерполяцию полиномом Ла-

гранжа, см. стр. 158

220 Глава 8. Приложения

lambda создает лямбда-выражение (безымянную функ -

цию). Лямбда-выражение может использоваться

в некоторых случаях как обычная функция, см.

стр. ??

laplace Функция реализует прямое преобразование Ла-

пласа , см. стр. 102-103

last Функция выделяет последний элемент списка , см.

стр. 23

lbfgs Функция осуществляет минимизацию ФНП, см.

стр. 160

ldisplay Функция печатает значения своих аргументов

вместе с их именем и метками "%t см. стр. 153

ldisp Функция печатает значения своих аргументов

вместе с метками "%t см. стр. 153

length Функция возвращает длину списка , см. стр. 21

lhs Функция выделяет левую часть уравнения, см.

стр. 63-67

limit функция осуществляет вычисление пределов , см.

стр. 70

linearinterpol Функция строит линейную интерполяцию , см.

стр. 158

linsolve Функция решает системы линейных и полиноми-

альных уравнений. Допускаются недоопределен-

ные системы, см. стр. 63-67

listarray Функция печатает содержимое массивов первого

и второго видов , см. стр. 26

load Функция загружает тот или иной файл :

load(someﬁle); Тип загрузки зависит от типа фай-

ла (макрос Maxima, программа на Lisp, бинарный

файл) , см. стр. 155

logcontract Функция компактифицирует логарифмы в дан-

ном выражении , см. стр. 45-45

make_array Функция оздает массивы третьего вида, содержи-

мое которых печатается автоматически , см. стр.

makelist Функция позволяет создавать списки , см. стр. 21

map Функция применяет заданную функцию к каждо-

му элементу списка , см. стр. 25

matrix Функция возвращает матрицу, заданную поэле-

ментно , см. стр. 30-35

matrixmap Функция для заполнения матрицы значениями

некоторой функции , см. стр. 30-35