Python求矩阵的内积、外积、克罗内克直积、Khatri-Rao积
文章目录
- 矩阵乘法
- 内积和外积
- 直积
- Khatri-Rao积
矩阵乘法
线性代数研究的核心对象是矩阵，所谓矩阵就是由 m m m行 n n n列的数组成的一个举行的数阵，从编程的角度理解，就是二维数组。
作为一种特殊的代数结构，需要为矩阵引入特殊的计算方法，其中矩阵加法、减法相对来说比较直观，而乘法比较复杂，矩阵乘法的定义为
( A B ) i j = ∑ a i k b k j (AB)_{ij}=\sum a_{ik}b_{kj} (AB)ij=∑aikbkj
在python中可直接通过
A
@
B
A@B
A@B来实现，而常用的乘号*表示元素之间一对一的乘法。
此外，numpy中提供了np.matmul进行矩阵乘法。
内积和外积
对于向量 a , b a,b a,b，其内积表示为
a ⋅ b = ∑ i a i b i a\cdot b=\sum_i a_ib_i a⋅b=i∑aibi
外积则得到一个矩阵
( A B ) i j = a i b j (AB)_{ij}=a_ib_j (AB)ij=aibj
而矩阵乘法则可表示为向量外积之和。
在numpy中，提供了内积inner和外积outer，令
a
=
[
1
2
3
4
]
,
b
=
[
0
3
2
1
]
a=\begin{bmatrix}1&2\\3&4\end{bmatrix}, b=\begin{bmatrix}0&3\\2&1\end{bmatrix}
a=[1324],b=[0231]，则计算结果如下
import numpy as np
a = np.array([1,2,3,1]).reshape(2,2)
b = np.array([0,3,2,1]).reshape(2,2)
print(np.inner(a,b))
'''
[[6 4]
[3 7]]
'''
print(np.outer(a,b))
'''
[[0 3 2 1]
[0 6 4 2]
[0 9 6 3]
[0 3 2 1]]
'''
直积
numpy和scipy.linalg中都实现了矩阵直积，又称克罗内克积
A ⊗ B = [ a 11 B a 12 B ⋯ a 1 n B a 21 B a 22 B ⋯ a 2 n B ⋯ ⋯ ⋯ ⋯ a m 1 B a m 2 B ⋯ a m n B ] A\otimes B=\begin{bmatrix} a_{11}B&a_{12}B&\cdots&a_{1n}B\\ a_{21}B&a_{22}B&\cdots&a_{2n}B\\ \cdots&\cdots&\cdots&\cdots\\ a_{m1}B&a_{m2}B&\cdots&a_{mn}B\\ \end{bmatrix} A⊗B= a11Ba21B⋯am1Ba12Ba22B⋯am2B⋯⋯⋯⋯a1nBa2nB⋯amnB
import scipy.linalg as sl
c = sl.kron(a,b)
print(c)
'''
[[0 3 0 6]
[2 1 4 2]
[0 9 0 3]
[6 3 2 1]]
'''
print(np.kron(a,b))
'''
[[0 3 0 6]
[2 1 4 2]
[0 9 0 3]
[6 3 2 1]]
'''
下面来看一下二者的速度
>>> import timeit
>>> y = np.random.rand(30,30)
>>> x = np.random.rand(30,30)
>>> timeit.timeit(lambda :np.kron(x,y), number=10)
0.023692000002483837
>>> timeit.timeit(lambda :sl.kron(x,y), number=10)
0.05632819999300409
看来这次是scipy输了。
Khatri-Rao积
Khatri-Rao积是两个具有相同列数的矩阵 A , B A,B A,B对应列向量的克罗内克积排列而成的，可表示为
A ⊙ B = [ a 1 ⊗ b 1 a 2 ⊗ b 2 ⋯ a n ⊗ b n ] A\odot B = \begin{bmatrix} a_1\otimes b_1&a_2\otimes b_2&\cdots&a_n\otimes b_n \end{bmatrix} A⊙B=[a1⊗b1a2⊗b2⋯an⊗bn]
Khatri-Rao积满足结合律，在scipy.linalg中，用khatri_rao可求矩阵的Khatri-Rao积。
d = sl.khatri_rao(a,b)
print(d)
'''
[[0 6]
[2 2]
[0 3]
[6 1]]
'''