背景

前段时间在知识星球上立了一个Flag，至少写10篇关于 Python，Matlab 和 C# 对比的总结。

这是第 4 篇，对比 Matlab 与 Numpy 中经常用到的各种矩阵，比如零矩阵、单位矩阵、全1矩阵、对角矩阵、均匀分布的随机矩阵、标准正态分布的随机矩阵等。

虽然 Numpy 定义了 matrix 类型，使用该 matrix 类型创建的是矩阵对象。但是由于 NumPy 中同时存在 ndarray 和 matrix 对象，因此用户很容易将两者弄混。这有违 Python 的“显式优于隐式”的原则，因此官方并不推荐在程序中使用 matrix。在这里，我们仍然用 ndarray 来介绍。

1. 零矩阵

零矩阵为矩阵内所有元素为零的矩阵。在具体任务中通常先用零矩阵来确定大小，随后为矩阵元素赋值。

在 Matlab和Python 中都用zeros()函数来实现。

在 Matlab 中，zeros(m,n)表示零矩阵，其中$m$为行数，$n$为列数。若为方阵，则以zeros(n)表示。

在 Numpy 中，numpy.zeros(shape, dtype=None, order='C')函数返回给定形状和类型的零数组。

【例1】求3行4列的零矩阵。

Matlab：

>> A = zeros(3,4);
>> disp(A)0     0     0     00     0     0     00     0     0     0

Python：

import numpy as npA = np.zeros([3, 4])
print(A)
# [[0. 0. 0. 0.]
#  [0. 0. 0. 0.]
#  [0. 0. 0. 0.]]

【例2】求4阶的零矩阵。

Matlab：

>> A = zeros(4);
>> disp(A)0     0     0     00     0     0     00     0     0     00     0     0     0

Python：

import numpy as npA = np.zeros([4, 4])
print(A)
# [[0. 0. 0. 0.]
#  [0. 0. 0. 0.]
#  [0. 0. 0. 0.]
#  [0. 0. 0. 0.]]A = np.zeros(4)
print(A)
# [0. 0. 0. 0.]

np.zeros(4)表示生成 Numpy 中的包含4个零元素的一维数组。

注意：在线性代数里面讲的维数和数组的维数不同，如线代中提到的 n 维行向量在 Numpy 中是一维数组，而线性代数中的 n 维列向量在 Numpy 中是一个shape为(n, 1)的二维数组。

2. 全1矩阵

矩阵中所有元素为1，称为全1矩阵。

在 Matlab和Numpy 中都用ones()函数来表示。

在 Matlab 中，ones(m,n)表示全1矩阵，其中$m$为行数，$n$为列数。若为方阵，则以ones(n)表示。

在 Numpy 中，numpy.ones(shape, dtype=None, order='C')函数返回给定形状和类型的全1数组。

【例1】求4阶全1矩阵。

Matlab：

>> A=ones(4);
>> disp(A)1     1     1     11     1     1     11     1     1     11     1     1     1

Python：

import numpy as npA = np.ones([4, 4])
print(A)
# [[1. 1. 1. 1.]
#  [1. 1. 1. 1.]
#  [1. 1. 1. 1.]
#  [1. 1. 1. 1.]]A = np.ones(4)
print(A)
# [1. 1. 1. 1.]

【例2】求3行4列的全1矩阵。

Matlab：

>> A=ones(3,4);
>> disp(A)1     1     1     11     1     1     11     1     1     1

Python：

import numpy as npA = np.ones([3, 4])
print(A)
# [[1. 1. 1. 1.]
#  [1. 1. 1. 1.]
#  [1. 1. 1. 1.]]

3. 单位矩阵

单位矩阵为主对角线元素全为1，而其它元素全为零的矩阵，数学上以$I$或$E$表示之。

在 Matlab和Numpy 中都用eye()函数来表示，而 Numpy 中还提供了identity()函数来专门生成方阵。

在 Matlab 中，eye(m,n)表示单位矩阵，其中$m$为行数，$n$为列数，若为方阵，则以eye(n)表示。

在 Numpy 中，numpy.eye(N, M=None, k=0, dtype=float, order='C')函数返回一个对角线上为1，其它地方为零的单位数组。numpy.identity(n, dtype=None)函数则返回一个方的单位数组。

【例1】求4行4列的单位矩阵。

Matlab：

>> I=eye(4);
>> disp(I)1     0     0     00     1     0     00     0     1     00     0     0     1

Python：

import numpy as npI = np.eye(4)
print(I)
# [[1. 0. 0. 0.]
#  [0. 1. 0. 0.]
#  [0. 0. 1. 0.]
#  [0. 0. 0. 1.]]A = np.identity(4)
print(A)
# [[1. 0. 0. 0.]
#  [0. 1. 0. 0.]
#  [0. 0. 1. 0.]
#  [0. 0. 0. 1.]]

【例2】求3行4列的单位矩阵。

Matlab：

>> I=eye(3,4);
>> disp(I)1     0     0     00     1     0     00     0     1     0

Python：

import numpy as npI = np.eye(3, 4)
print(I)
# [[1. 0. 0. 0.]
#  [0. 1. 0. 0.]
#  [0. 0. 1. 0.]]

4. 对角矩阵

在 Matlab和Numpy 中都用diag()函数来表示，并且含义相同。

若书写格式为X=diag(v,k)。将向量$v$写入矩阵$X$的主对角线上，而矩阵$X$的其它元素为零。$k$表示上移或下移行数，正表示上移，负表示下移，$k=0$则恰好在主对角线上，当$k=0$时可以默认不写。

另一种书写格式为v=diag(X,k)。它是从矩阵$X$中提取对角线元素到向量$v$上。$k=0$或默认，则提取主对角线元素，否则提取上移$k$行或下移$k$行的对角线元素，正号表示上移，负号表示下移。

【例1】已知行向量v=(1 2 3 4)，将$v$向量元素写入矩阵的主对角线，求对角矩阵，以及上移一行的对角矩阵和下移一行的对角矩阵。

Matlab：

>> v=[1,2,3,4];
>> lambda = diag(v);
>> disp(lambda)1     0     0     00     2     0     00     0     3     00     0     0     4>> lambda = diag(v,1);
>> disp(lambda)0     1     0     0     00     0     2     0     00     0     0     3     00     0     0     0     40     0     0     0     0>> lambda = diag(v,-1);
>> disp(lambda)0     0     0     0     01     0     0     0     00     2     0     0     00     0     3     0     00     0     0     4     0

Python：

import numpy as npv = [1, 2, 3, 4]
lamb = np.diag(v)
# print(lamb)
# [[1 0 0 0]
#  [0 2 0 0]
#  [0 0 3 0]
#  [0 0 0 4]]lamb = np.diag(v, k=1)
print(lamb)
# [[0 1 0 0 0]
#  [0 0 2 0 0]
#  [0 0 0 3 0]
#  [0 0 0 0 4]lamb = np.diag(v, k=-1)
print(lamb)
# [[0 0 0 0 0]
#  [1 0 0 0 0]
#  [0 2 0 0 0]
#  [0 0 3 0 0]
#  [0 0 0 4 0]]

【例2】
已知u=(2 3 5 8 11)，创建Vandermonde矩阵$X$，提取$X$的主对角线赋予向量$v$，主对角线上移一行赋予向量$v1$，下移一行赋予向量$v2$。

>> u=[2,3,5,8,11];
>> X=vander(u);
>> disp(X)16           8           4           2           181          27           9           3           1625         125          25           5           14096         512          64           8           114641        1331         121          11           1>> v=diag(X);
>> disp(v)16272581>> v1=diag(X,1)';
>> disp(v1)8     9     5     1>> v2=diag(X,-1)';
>> disp(v2)81   125    64    11

Python：

import numpy as npdef vander(u):n = len(u)a = np.ones([n, n], dtype=np.int)for i in range(n - 2, -1, -1):a[:, i] = u ** (n - i - 1)return au = np.array([2, 3, 5, 8, 11])
X = vander(u)
print(X)
# [[   16     8     4     2     1]
#  [   81    27     9     3     1]
#  [  625   125    25     5     1]
#  [ 4096   512    64     8     1]
#  [14641  1331   121    11     1]]v = np.diag(X)
print(v)
# [16 27 25  8  1]v1 = np.diag(X, k=1)
print(v1)
# [8 9 5 1]v2 = np.diag(X, k=-1)
print(v2)
# [ 81 125  64  11]

5. 均匀分布的随机矩阵

在 Matlab和Numpy 中都可以生成均匀分布的随机矩阵。

在 Matlab 中，有rand()，rands()，randi()三种方式：

• rand(n)得到一个$n\times n$的矩阵。
• rands(n)得到一个$n\times 1$的列向量。
• rand(m,n)得到一个由均匀分布在(0,1)区间内元素组成的$m\times n$的矩阵。
• rands(m,n)得到一个由均匀分布在(-1,1)区间内元素组成的$m\times n$的矩阵。
• randi(iMax，n)是一个由均匀分布在[1，iMax]区间内整型元素组成的$n\times n$的矩阵。
• randi(iMax，m，n)是一个由均匀分布在[1，iMax]区间内整型元素组成的$m\times n$的矩阵。
• randi([iMin,iMax],m,n)是一个由均匀分布在[iMin，iMax]区间内整型元素组成的$m\times n$矩阵。

在 Numpy 中，有uniform()，rand()，randint()三种方式：

• numpy.random.uniform(low=0.0, high=1.0, size=None) 从均匀分布中抽取样本，区间为[low, high) 。
• numpy.random.rand() 作为uniform(d0, d1, ..., dn)的特列，可以得到[0,1)之间的均匀分布的随机数。
• numpy.random.randint(low, high=None, size=None, dtype='l') 作为uniform的另一特例，可以得到[low,high)之间均匀分布的随机整数，如果high=None，那么区间为[0,low)。

【例1】求两个4阶在(0,1)范围内的均匀分布的随机矩阵，两个4阶在[1,100]范围内的整型均匀分布的随机数，两个4×1的在(-1,1)范围内的列向量，分别赋值给A1，A2。

Matlab：

>> A1=rand(4);
>> A2=rand(4);
>> disp(A1)0.2232    0.6585    0.3061    0.94720.8876    0.8102    0.2829    0.67460.7311    0.2015    0.4056    0.70160.3146    0.1130    0.8241    0.6380>> disp(A2)0.2420    0.6504    0.2952    0.63280.9741    0.0814    0.3742    0.21570.8539    0.8477    0.7121    0.41150.7950    0.1764    0.6574    0.4196>> A1 = randi(100,4);
>> A2 = randi(100,4);
>> disp(A1)47    56    32    6165    25    94    9383     5    19    9225    76    55     4>> disp(A2)69    86    47    414    84    98     351    80    59     930    81    73    59>> A1=rands(4);
>> A2=rands(4);
>> disp(A1)0.27960.40790.5119-0.9015>> disp(A2)0.7674-0.40590.3743-0.3374

Python：

import numpy as npA1 = np.random.rand(4, 4)
A2 = np.random.rand(4, 4)
print(A1)
# [[0.26928205 0.36481389 0.66676956 0.86820657]
#  [0.91229641 0.04239373 0.88579269 0.41331319]
#  [0.45072909 0.16651979 0.09096757 0.63638447]
#  [0.96256553 0.12203043 0.88213325 0.70119275]]print(A2)
# [[0.99001364 0.59498628 0.18862191 0.24256341]
#  [0.1211043  0.32129381 0.9460552  0.05122717]
#  [0.63987355 0.61879931 0.31145657 0.32688994]
#  [0.73492717 0.43796817 0.9866393  0.93774641]]A1 = np.random.randint(1, 100, [4, 4])
A2 = np.random.randint(1, 100, [4, 4])
print(A1)
# [[ 9 44 72 10]
#  [33 89 55 47]
#  [77 59 84 83]
#  [43 23 86 59]]print(A2)
# [[ 5 89 33 11]
#  [21 95 73 92]
#  [ 2 55 30 25]
#  [54 40 65 83]]A1 = np.random.uniform(-1, 1, 4)
A2 = np.random.uniform(-1, 1, 4)
print(A1.reshape(4, 1))
# [[ 0.93055573]
# #  [-0.86806273]
# #  [-0.8755103 ]
# #  [ 0.54197419]]print(A2.reshape(4, 1))
# [[-0.32703134]
#  [ 0.80505337]
#  [-0.83637509]
#  [ 0.69088863]]

【例2】求三个$3\times 4$阶均匀分布的随机矩阵，范围在(0,1)，(-1,1)和[10,100]，分别赋值给A1，A2，A3。

Matlab：

>> A1=rand(3,4);
>> disp(A1)0.7850    0.5495    0.0590    0.17400.0589    0.5529    0.1722    0.84260.3510    0.1447    0.4489    0.4642>> A2=rands(3,4);
>> disp(A2)-0.3603   -0.5592    0.2269    0.0017-0.4025   -0.7993   -0.8550    0.62550.2765    0.6801    0.4448   -0.2550>> A3 = randi([10,100],3,4);
>> disp(A3)59    66    63    8131    55    50    6315    88    98    34

Python：

import numpy as npA1 = np.random.rand(3, 4)
print(A1)
# [[0.39787012 0.13709285 0.43893989 0.96172473]
#  [0.43530643 0.05443109 0.36316048 0.94816304]
#  [0.55265213 0.60320904 0.24952952 0.83926361]]A2 = np.random.uniform(-1, 1, [3, 4])
print(A2)
# [[ 0.32297696 -0.78124371 -0.39123924  0.24991417]
#  [ 0.85884802 -0.23765347  0.56611597  0.90097203]
#  [-0.72628246 -0.10791575 -0.11244534  0.40308801]]A3 = np.random.randint(10, 100, [3, 4])
print(A3)
# [[29 53 65 89]
#  [10 77 44 33]
#  [28 85 35 42]]

【例3】使用种子使得每次随机数相同。

Matlab：

>> s=rng(20200821); % 保存随机数seed
>> disp(s)Type: 'twister'Seed: 20200821State: [625×1 uint32]>> rng(s) % 调用随机数seed
>> A1=rand(4);
>> disp(A1)0.0172    0.6571    0.8932    0.31100.6420    0.8373    0.8906    0.38150.0978    0.8899    0.6041    0.88750.4353    0.7270    0.5803    0.3294>> rng(s)
>> A2=rand(4);
>> disp(A2)0.0172    0.6571    0.8932    0.31100.6420    0.8373    0.8906    0.38150.0978    0.8899    0.6041    0.88750.4353    0.7270    0.5803    0.3294

Python：

import numpy as npnp.random.seed(20200821)
A1 = np.random.rand(4, 4)
print(A1)
# [[0.32452405 0.49496274 0.88643592 0.75054701]
#  [0.88940824 0.4339118  0.45194462 0.80973677]
#  [0.39462383 0.72057367 0.21335562 0.87752275]
#  [0.42963472 0.92526431 0.22639798 0.83541325]]np.random.seed(20200821)
A2 = np.random.rand(4, 4)
print(A2)
# [[0.32452405 0.49496274 0.88643592 0.75054701]
#  [0.88940824 0.4339118  0.45194462 0.80973677]
#  [0.39462383 0.72057367 0.21335562 0.87752275]
#  [0.42963472 0.92526431 0.22639798 0.83541325]]

6. 正态分布的随机矩阵

在 Matlab和Numpy 中都用randn()函数来表示抽取标准正态分布的样本。

注：既然是随机数，无论产生多少个，都只能算是样本，所谓期望，是对样本数量足够多的情况而言的。一般而言，样本数量越大，其均值越接近0，但并不是说就一定为0。

【例1】列出4阶正态分布的随机矩阵。

Matlab：

randn(m,n)产生一个期望为0，标准差为1且服从正态分布的随机矩阵，其中$m$为行数，$n$为列数，若为方阵，则以randn(n)表示。

>> A = randn(4);
>> disp(A)0.5611   -0.6045   -0.0025   -1.1567-1.1985   -0.1034    0.1916    0.6170-1.0359    1.1587   -0.4981    0.0227-0.5967   -1.8604    0.9566    0.8366>> A=randn(100000,1);
>> histogram(A)

Python：

• numpy.random.randn(d0, d1, ..., dn)从均值为0，标准差为1的标准正态分布中抽取样本。

import numpy as np
import matplotlib.pyplot as pltA = np.random.randn(4, 4)
print(A)
# [[-1.6652191   2.57166227 -2.12784124  1.52365121]
#  [-1.03609184 -0.32386994  0.45424631  0.58238759]
#  [-0.02458478 -1.69242129  0.01512539  0.8165994 ]
#  [-0.39558196 -1.62699892 -0.5198507   2.29518816]]A = np.random.randn(100000, 1)
plt.hist(A, bins=20)
plt.show()

总结

以上总结的不一定全，但先有个框架等后面在实践的过程中慢慢补充。今天就到这里吧。See You！

后台回复「搜搜搜」，随机获取电子资源！
欢迎关注，请扫描二维码：