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

Линейная алгебра

Раздел математики, связанный с математическими структурами, закрытыми операциями сложения и скалярного умножения и включающий теорию систем линейных уравнений, матриц, определителей, векторных пространств. , и линейные преобразования - Источник: Merriam Webster

Если вы еще не прошли - Знакомство с Numpy, обязательно пройдите его! Нам понадобится Numpy для реализации большинства математических концепций!

СОДЕРЖАНИЕ

  • Умножение матриц
  • Полномочия матрицы
  • Транспонировать
  • Обратный
  • Трассировка
  • Определитель
  • Характеристические многочлены и теорема Кэли-Гамильтона
  • Прогнозы
  • Линейные системы
  • Исключение Гаусса
  • Элементарные операции со строками
  • Обратный
  • Решение системы

Начиная

Библиотеки, которые вам нужно будет импортировать для этого руководства - Numpy, Matplotlib и Scipy

import numpy as np
import matplotlib.pyplot as plt
import scipy.linalg as la
%matplotlib inline

Вычислим 2I + 3A − AB (I - единичная матрица размера 2) для -

A = np.array([[1,3],[-1,7]])
print(A)
Out:
[[ 1  3]
 [-1  7]]
B = np.array([[5,2],[1,2]])
print(B)
Out:
[[5 2]
 [1 2]]
I = np.eye(2)
print(I)
Out:
[[1. 0.]
 [0. 1.]]
2*I + 3*A - A@B
Out:
array([[-3.,  1.],
       [-5., 11.]])

Матричные силы

Для степеней матрицы нет символа, поэтому мы должны импортировать функцию matrix_power из подпакета numpy.linalg.

from numpy.linalg import matrix_power as mpow
M = np.array([[3,4],[-1,5]])
print(M)
Out:
[[ 3  4]
 [-1  5]]
mpow(M,2)
Out:
array([[ 5, 32],
       [-8, 21]])
mpow(M,5)
Out:
array([[-1525,  3236],
       [ -809,    93]])

Сравните с оператором матричного умножения:

M @ M @ M @ M @ M
Out:
array([[-1525,  3236],
       [ -809,    93]])

Транспонировать

Мы можем выполнить транспонирование с помощью .T attribute:

print(M)
Out:
[[ 3  4]
 [-1  5]]
print(M.T)
Out:
[[ 3 -1]
 [ 4  5]]

Обратите внимание, что M.MT - это симметричная матрица:

M @ M.T
Out:
array([[25, 17],
       [17, 26]])

Обратный

Обратное можно найти с помощью функции scipy.linalg.inv:

A = np.array([[1,2],[3,4]])
print(A)
Out:
[[1 2]
 [3 4]]
la.inv(A)
Out:
array([[-2. ,  1. ],
       [ 1.5, -0.5]])

След

Мы можем найти след матрицы с помощью функции numpy.trace:

np.trace(A)
Out:
5

Детерминант

Находим определитель с помощью функции scipy.linalg.det:

A = np.array([[1,2],[3,4]])
print(A)
Out:
[[1 2]
 [3 4]]
la.det(A)
Out:
-2.0

Характеристические многочлены и теорема Кэли-Гамильтона

Характеристический многочлен квадратной матрицы A 2 на 2 равен

Теорема Кэли-Гамильтона утверждает, что любая квадратная матрица удовлетворяет своему характеристическому многочлену.

Для матрицы A размера 2 это означает, что

Давайте проверим теорему Кэли-Гамильтона для нескольких разных матриц -

print(A)
Out:
[[1 2]
 [3 4]]
trace_A = np.trace(A)
det_A = la.det(A)
I = np.eye(2)
A @ A - trace_A * A + det_A * I
Out:
array([[0., 0.],
       [0., 0.]])

Прогнозы

Формула для проецирования вектора v на вектор w:

Давайте представим функцию proj , которая вычисляет проекцию v на w:

def proj(v,w):
    '''Project vector v onto w.'''
    v = np.array(v)
    w = np.array(w)
    return np.sum(v * w)/np.sum(w * w) * w  # or (v @ w)/(w @ w) * w
proj([1,2,3],[1,1,1])
Out:
array([2., 2., 2.])

Решение линейных систем

Линейные системы

линейная система уравнений - это набор линейных уравнений

В матричной записи линейной системой является Ax = b, где

Гауссово исключение

Общая процедура решения линейной системы уравнений называется исключением Гаусса. Идея состоит в том, чтобы выполнить элементарные операции со строками, чтобы привести систему к ее эшелонированной форме, а затем решить.

Элементарные операции со строками

Элементарные операции со строками включают:

  1. Добавить k раз строку j в строку i
  2. Умножить строку i на скаляр k
  3. Поменять местами строки i и j

Каждая из элементарных операций со строкой является результатом умножения матрицы на элементарную матрицу. Чтобы добавить k раз строку i к строке j в матрице A, мы умножаем A на матрицу E, где E равно единичной матрице , за исключением i, j запись - Ei, j = k.

Например, если A равно 3 на 3, и мы хотим добавить 3 раза строку 2 к строке 0 (с использованием 0 индексации), тогда -

Проверим расчет:

A = np.array([[1,1,2],[-1,3,1],[0,5,2]])
print(A)
Out:
[[ 1  1  2]
 [-1  3  1]
 [ 0  5  2]]
E1 = np.array([[1,0,3],[0,1,0],[0,0,1]])
print(E1)
Out:
[[1 0 3]
 [0 1 0]
 [0 0 1]]
E1 @ A
Out:
array([[ 1, 16,  8],
       [-1,  3,  1],
       [ 0,  5,  2]])

Чтобы умножить k раз строку i в матрице A, мы умножаем A на матрицу E, где E равно единичной матрице , за исключением i, j : Ei, i = k.

Например, если A равно 3 на 3, и мы хотим умножить строку 1 на -2, тогда -

Проверим расчет:

E2 = np.array([[1,0,0],[0,-2,0],[0,0,1]])
print(E2)
Out:
[[ 1  0  0]
 [ 0 -2  0]
 [ 0  0  1]]
E2 @ A
Out:
array([[ 1,  1,  2],
       [ 2, -6, -2],
       [ 0,  5,  2]])

Наконец, чтобы поменять местами строку i и строку j в матрице A, мы умножаем A на матрицу E, где E равно единичной матрице , за исключением Ei, i = 0, Ej, j = 0, Ei, j = 1 и Ej, i = 1.

Например, если A равно 3 на 3, и мы хотим поменять местами строку 1 и строку 2, тогда -

Проверим расчет:

E3 = np.array([[1,0,0],[0,0,1],[0,1,0]])
print(E3)
Out:
[[1 0 0]
 [0 0 1]
 [0 1 0]]
E3 @ A
Out:
array([[ 1,  1,  2],
       [ 0,  5,  2],
       [-1,  3,  1]])

Реализация

Теперь напишем функции для реализации операций с элементарной строкой. Прежде всего, давайте напишем функцию с именем add_rows, которая принимает входные параметры A, k, i и j и возвращает массив NumPy, полученный в результате добавления k раз строки j к строке i в матрице A.

Если i = j, то предположим, что функция масштабирует строку i на k + 1, поскольку это будет результатом k умножить строку, которую я добавил к строке i.

def add_row(A,k,i,j):
    "Add k times row j to row i in matrix A."
    n = A.shape[0]
    E = np.eye(n)
    if i == j:
        E[i,i] = k + 1
    else:
        E[i,j] = k
    return E @ A

Давайте проверим нашу функцию:

M = np.array([[1,1],[3,2]])
print(M)
Out:
[[1 1]
 [3 2]]
add_row(M,2,0,1)
Out:
array([[7., 5.],
       [3., 2.]])

Давайте напишем функцию с именем scale_row, которая принимает 3 входных параметра A, k и i и возвращает матрицу, которая получается из k раз строки i в матрице A.

def scale_row(A,k,i):
    "Multiply row i by k in matrix A."
    n = A.shape[0]
    E = np.eye(n)
    E[i,i] = k
    return E @ A

Давайте проверим нашу функцию:

M = np.array([[3,1],[-2,7]])
print(M)
Out:
[[ 3  1]
 [-2  7]]
scale_row(M,3,1)
Out:
array([[ 3.,  1.],
       [-6., 21.]])

Давайте напишем функцию с именем switch_rows, которая принимает 3 входных параметра A, i и j и возвращает матрицу, полученную в результате переключения строк i и j в матрице A.

def switch_rows(A,i,j):
    "Switch rows i and j in matrix A."
    n = A.shape[0]
    E = np.eye(n)
    E[i,i] = 0
    E[j,j] = 0
    E[i,j] = 1
    E[j,i] = 1
    return E @ A

Давайте проверим нашу функцию:

A = np.array([[1,1,1],[1,-1,0]])
print(A)
Out:
[[ 1  1  1]
 [ 1 -1  0]]
switch_rows(A,0,1)
Out:
array([[ 1., -1.,  0.],
       [ 1.,  1.,  1.]])

Примеры

Найдите обратное

Применим наши функции к расширенной матрице [M | I], чтобы найти обратную матрицу M:

M = np.array([[5,4,2],[-1,2,1],[1,1,1]])
print(M)
Out:
[[ 5  4  2]
 [-1  2  1]
 [ 1  1  1]]
A = np.hstack([M,np.eye(3)])
print(A)
Out:
[[ 5.  4.  2.  1.  0.  0.]
 [-1.  2.  1.  0.  1.  0.]
 [ 1.  1.  1.  0.  0.  1.]]
A1 = switch_rows(A,0,2)
print(A1)
Out:
[[ 1.  1.  1.  0.  0.  1.]
 [-1.  2.  1.  0.  1.  0.]
 [ 5.  4.  2.  1.  0.  0.]]
A2 = add_row(A1,1,1,0)
print(A2)
Out:
[[1. 1. 1. 0. 0. 1.]
 [0. 3. 2. 0. 1. 1.]
 [5. 4. 2. 1. 0. 0.]]
A3 = add_row(A2,-5,2,0)
print(A3)
Out:
[[ 1.  1.  1.  0.  0.  1.]
 [ 0.  3.  2.  0.  1.  1.]
 [ 0. -1. -3.  1.  0. -5.]]
A4 = switch_rows(A3,1,2)
print(A4)
Out:
[[ 1.  1.  1.  0.  0.  1.]
 [ 0. -1. -3.  1.  0. -5.]
 [ 0.  3.  2.  0.  1.  1.]]
A5 = scale_row(A4,-1,1)
print(A5)
Out:
[[ 1.  1.  1.  0.  0.  1.]
 [ 0.  1.  3. -1.  0.  5.]
 [ 0.  3.  2.  0.  1.  1.]]
A6 = add_row(A5,-3,2,1)
print(A6)
Out:
[[  1.   1.   1.   0.   0.   1.]
 [  0.   1.   3.  -1.   0.   5.]
 [  0.   0.  -7.   3.   1. -14.]]
A7 = scale_row(A6,-1/7,2)
print(A7)
Out:
[[ 1.          1.          1.          0.          0.          1.        ]
 [ 0.          1.          3.         -1.          0.          5.        ]
 [ 0.          0.          1.         -0.42857143 -0.14285714  2.        ]]
A8 = add_row(A7,-3,1,2)
print(A8)
Out:
[[ 1.          1.          1.          0.          0.          1.        ]
 [ 0.          1.          0.          0.28571429  0.42857143 -1.        ]
 [ 0.          0.          1.         -0.42857143 -0.14285714  2.        ]]
A9 = add_row(A8,-1,0,2)
print(A9)
Out:
[[ 1.          1.          0.          0.42857143  0.14285714 -1.        ]
 [ 0.          1.          0.          0.28571429  0.42857143 -1.        ]
 [ 0.          0.          1.         -0.42857143 -0.14285714  2.        ]]
A10 = add_row(A9,-1,0,1)
print(A10)
Out:
[[ 1.          0.          0.          0.14285714 -0.28571429  0.        ]
 [ 0.          1.          0.          0.28571429  0.42857143 -1.        ]
 [ 0.          0.          1.         -0.42857143 -0.14285714  2.        ]]

Давайте проверим, правильно ли мы нашли обратное M:

Minv = A10[:,3:]
print(Minv)
Out:
[[ 0.14285714 -0.28571429  0.        ]
 [ 0.28571429  0.42857143 -1.        ]
 [-0.42857143 -0.14285714  2.        ]]
result = Minv @ M
print(result)
Out:
[[ 1.00000000e+00  4.44089210e-16  2.22044605e-16]
 [-6.66133815e-16  1.00000000e+00 -2.22044605e-16]
 [ 0.00000000e+00  0.00000000e+00  1.00000000e+00]]

Успех! Мы сможем увидеть результат более четко, если округлим до 15 десятичных знаков:

np.round(result,15)
Out:
array([[ 1.e+00,  0.e+00,  0.e+00],
       [-1.e-15,  1.e+00, -0.e+00],
       [ 0.e+00,  0.e+00,  1.e+00]])

Решить систему

Давайте воспользуемся нашими функциями для выполнения исключения Гаусса и решения простой системы уравнений Ax = b.

A = np.array([[6,15,1],[8,7,12],[2,7,8]])
print(A)
Out:
[[ 6 15  1]
 [ 8  7 12]
 [ 2  7  8]]
b = np.array([[2],[14],[10]])
print(b)
Out:
[[ 2]
 [14]
 [10]]

Сформируйте расширенную матрицу M:

M = np.hstack([A,b])
print(M)
[[ 6 15  1  2]
 [ 8  7 12 14]
 [ 2  7  8 10]]

Выполнение операций со строками:

M1 = scale_row(M,1/6,0)
print(M1)
Out:
[[ 1.          2.5         0.16666667  0.33333333]
 [ 8.          7.         12.         14.        ]
 [ 2.          7.          8.         10.        ]]
M2 = add_row(M1,-8,1,0)
print(M2)
Out:
[[  1.           2.5          0.16666667   0.33333333]
 [  0.         -13.          10.66666667  11.33333333]
 [  2.           7.           8.          10.        ]]
M3 = add_row(M2,-2,2,0)
print(M3)
Out:
[[  1.           2.5          0.16666667   0.33333333]
 [  0.         -13.          10.66666667  11.33333333]
 [  0.           2.           7.66666667   9.33333333]]
M4 = scale_row(M3,-1/13,1)
print(M4)
Out:
[[ 1.          2.5         0.16666667  0.33333333]
 [ 0.          1.         -0.82051282 -0.87179487]
 [ 0.          2.          7.66666667  9.33333333]]
M5 = add_row(M4,-2,2,1)
print(M5)
Out:
[[ 1.          2.5         0.16666667  0.33333333]
 [ 0.          1.         -0.82051282 -0.87179487]
 [ 0.          0.          9.30769231 11.07692308]]
M6 = scale_row(M5,1/M5[2,2],2)
print(M6)
Out:
[[ 1.          2.5         0.16666667  0.33333333]
 [ 0.          1.         -0.82051282 -0.87179487]
 [ 0.          0.          1.          1.19008264]]
M7 = add_row(M6,-M6[1,2],1,2)
print(M7)
Out:
[[1.         2.5        0.16666667 0.33333333]
 [0.         1.         0.         0.1046832 ]
 [0.         0.         1.         1.19008264]]
M8 = add_row(M7,-M7[0,2],0,2)
print(M8)
Out:
[[1.         2.5        0.         0.13498623]
 [0.         1.         0.         0.1046832 ]
 [0.         0.         1.         1.19008264]]
M9 = add_row(M8,-M8[0,1],0,1)
print(M9)
Out:
[[ 1.          0.          0.         -0.12672176]
 [ 0.          1.          0.          0.1046832 ]
 [ 0.          0.          1.          1.19008264]]

Успех! Решение Ax = b:

x = M9[:,3].reshape(3,1)
print(x)
Out:
[[-0.12672176]
 [ 0.1046832 ]
 [ 1.19008264]]

Или мы можем сделать это простым способом…

x = la.solve(A,b)
print(x)
Out:
[[-0.12672176]
 [ 0.1046832 ]
 [ 1.19008264]]

scipy.linalg.solve

Нас больше всего интересуют линейные системы Ax = b, в которых существует единственное решение x. Это тот случай, когда A - квадратная матрица (m = n) и det (A) ≠ 0. Для решения такой системы мы можем использовать функцию scipy.linalg.solve.

Функция возвращает решение системы уравнений Ax = b. Например:

A = np.array([[1,1],[1,-1]])
print(A)
Out:
[[ 1  1]
 [ 1 -1]]
b1 = np.array([2,0])
print(b1)
Out:
[2 0]

И решаем:

x1 = la.solve(A,b1)
print(x1)
Out:
[1. 1.]

Обратите внимание, что выходные данные x возвращаются как 1D Numpy-массив, когда вектор b (правая часть) вводится как 1D Numpy-массив. Если мы вводим b как 2D-массив Numpy, тогда на выходе будет 2D-массив Numpy.

Например:

A = np.array([[1,1],[1,-1]])
b2 = np.array([2,0]).reshape(2,1)
x2 = la.solve(A,b2)
print(x2)
Out:
[[1.]
 [1.]]

Наконец, если правая часть b представляет собой матрицу, то на выходе будет матрица того же размера. Это решение Ax = b, когда b - матрица.

Например:

A = np.array([[1,1],[1,-1]])
b3 = np.array([[2,2],[0,1]])
x3 = la.solve(A,b3)
print(x3)
Out:
[[1.  1.5]
 [1.  0.5]]

Простой пример

Давайте вычислим решение системы уравнений

Создайте матрицу коэффициентов:

A = np.array([[2,1],[1,1]])
print(A)
Out:
[[2 1]
 [1 1]]

И вектор b:

b = np.array([1,-1]).reshape(2,1)
print(b)
Out:
[[ 1]
 [-1]]

И решаем:

x = la.solve(A,b)
print(x)
Out:
[[ 2.]
 [-3.]]

Мы можем проверить решение, вычислив обратное к A:

Ainv = la.inv(A)
print(Ainv)
Out:
[[ 1. -1.]
 [-1.  2.]]

И умножьте обратное значение A и b, чтобы найти x:

x = Ainv @ b
print(x)
Out:
[[ 2.]
 [-3.]]

Получаем тот же результат. Успех!

Обратить или решить

Плохая идея использовать обратное значение A для решения Ax = b, если A большое. Это слишком затратно с точки зрения вычислений.

Давайте создадим большую случайную матрицу A и вектор b и вычислим решение x двумя способами:

N = 1000
A = np.random.rand(N,N)
b = np.random.rand(N,1)

Проверьте первые записи A:

A[:3,:3]
Out:
array([[0.35754719, 0.63135432, 0.6572258 ],
       [0.18450506, 0.14639832, 0.23528745],
       [0.27576474, 0.46264005, 0.26589724]])

А для b:

b[:4,:]
Out:
array([[0.82726751],
       [0.96946096],
       [0.31351176],
       [0.63757837]])

Теперь сравним скорость scipy.linalg.solve с scipy.linalg.inv:

%%timeit
x = la.solve(A,b)
Out:
2.77 s ± 509 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)
%%timeit
x = la.inv(A) @ b
Out:
4.46 s ± 2.04 s per loop (mean ± std. dev. of 7 runs, 1 loop each)

Решение с scipy.linalg.solve примерно в два раза быстрее!

Это все, что касается Части 1! Я знаю, что сразу нужно принять во внимание многое! Но ты дожил до конца! Престижность на этом! Не забудьте ознакомиться с другими частями этой статьи - Часть 2, Часть 3, Часть 4 и Часть 5 !

Получите доступ к экспертному обзору - Подпишитесь на DDI Intel

Дополнительные ресурсы и ссылки

Есть много других хороших ресурсов, если вы все еще заинтересованы в получении максимальной отдачи от этой темы -





Для полной реализации ознакомьтесь с моим репозиторием GitHub -