面向MATLAB用户的NumPy#

简介#

MATLAB®和NumPy有很多共同之处,但NumPy是为与Python一起使用而创建的,而不是作为MATLAB的克隆。本指南将帮助MATLAB用户开始使用NumPy。

一些关键区别#

在MATLAB中,基本类型,即使对于标量,也是多维数组。MATLAB中的数组赋值存储为双精度浮点数的二维数组,除非您指定维度和类型。这些数组的二维实例上的操作基于线性代数中的矩阵运算。

在NumPy中,基本类型是多维array。NumPy中的数组赋值通常存储为n维数组,其类型足以保存序列中的对象,除非您指定维度和类型。NumPy逐元素执行运算,因此使用*对二维数组进行乘法不是矩阵乘法——而是逐元素乘法。(自Python 3.5起可用的@运算符可用于传统的矩阵乘法。)

MATLAB从1开始对索引编号;a(1)是第一个元素。参见注释 INDEXING

NumPy与Python一样,从0开始对索引编号;a[0]是第一个元素。

MATLAB的脚本语言是为线性代数而创建的,因此某些数组操作的语法比NumPy的更简洁。另一方面,添加GUI和创建完整应用程序的API或多或少是事后才考虑的。

NumPy基于Python,这是一种通用语言。NumPy的优势是可以访问Python库,包括:SciPyMatplotlibPandasOpenCV等等。此外,Python通常嵌入到其他软件中作为脚本语言,从而也可以在那里使用NumPy。

MATLAB数组切片使用按值传递语义,并采用延迟复制写入方案,以防止在需要之前创建副本。切片操作会复制数组的部分内容。

NumPy数组切片使用按引用传递,不复制参数。切片操作是数组的视图。

粗略等价物#

下表给出了某些常用MATLAB表达式的粗略等价物。这些是相似的表达式,而不是等价物。有关详细信息,请参阅文档

在下表中,假定您已在Python中执行了以下命令

import numpy as np
from scipy import io, integrate, linalg, signal
from scipy.sparse.linalg import cg, eigs

另请在下文中假设,如果注释谈论“矩阵”,则参数是二维实体。

通用等价物#

MATLAB

NumPy

注释

help func

info(func)help(func)func?(在IPython中)

获取函数 *func* 的帮助

which func

参见注释 HELP

找出 *func* 的定义位置

type func

np.source(func)func??(在IPython中)

打印 *func* 的源代码(如果不是原生函数)

% comment

# comment

使用文本 comment 注释一行代码

for i=1:3
    fprintf('%i\n',i)
end
for i in range(1, 4):
   print(i)

使用range编写for循环打印数字1、2和3

a && b

a and b

短路逻辑与运算符(Python原生运算符);仅限标量参数

a || b

a or b

短路逻辑或运算符(Python原生运算符);仅限标量参数

>> 4 == 4
ans = 1
>> 4 == 5
ans = 0
>>> 4 == 4
True
>>> 4 == 5
False

Python中的布尔对象TrueFalse,而不是MATLAB的逻辑类型10

a=4
if a==4
    fprintf('a = 4\n')
elseif a==5
    fprintf('a = 5\n')
end
a = 4
if a == 4:
    print('a = 4')
elif a == 5:
    print('a = 5')

创建一个if-else语句来检查a是否为4或5并打印结果

1*i1*j1i1j

1j

复数

eps

np.finfo(float).epsnp.spacing(1)

双精度浮点数中,从1到下一个可表示的实数的距离

load data.mat

io.loadmat('data.mat')

加载保存到文件data.mat中的MATLAB变量。(注意:在MATLAB/Octave中将数组保存到data.mat时,请使用最新的二进制格式。scipy.io.loadmat 将创建一个包含保存的数组和更多信息的字典。)

ode45

integrate.solve_ivp(f)

使用Runge-Kutta 4,5积分ODE

ode15s

integrate.solve_ivp(f, method='BDF')

使用BDF方法积分ODE

线性代数等价物#

MATLAB

NumPy

注释

ndims(a)

np.ndim(a)a.ndim

数组a的维度数

numel(a)

np.size(a)a.size

数组a的元素个数

size(a)

np.shape(a)a.shape

数组a的“大小”

size(a,n)

a.shape[n-1]

获取数组a的第n维的元素个数。(请注意,MATLAB使用基于1的索引,而Python使用基于0的索引,参见注释INDEXING

[ 1 2 3; 4 5 6 ]

np.array([[1., 2., 3.], [4., 5., 6.]])

定义一个2x3的二维数组

[ a b; c d ]

np.block([[a, b], [c, d]])

从块abcd构造矩阵

a(end)

a[-1]

访问MATLAB向量(1xn或nx1)或一维NumPy数组a(长度n)中的最后一个元素

a(2,5)

a[1, 4]

访问二维数组a中第二行第五列的元素

a(2,:)

a[1]a[1, :]

二维数组a的整行第二行

a(1:5,:)

a[0:5]a[:5]a[0:5, :]

二维数组a的前5行

a(end-4:end,:)

a[-5:]

二维数组a的最后5行

a(1:3,5:9)

a[0:3, 4:9]

二维数组 a 的第一到第三行和第五到第九列。

a([2,4,5],[1,3])

a[np.ix_([1, 3, 4], [0, 2])]

第2、4和5行以及第1和3列。这允许修改矩阵,并且不需要规则切片。

a(3:2:21,:)

a[2:21:2,:]

数组 a 中每隔一行,从第三行开始到第二十一行

a(1:2:end,:)

a[::2, :]

从第一行开始,每隔一行取a 的数据

a(end:-1:1,:)flipud(a)

a[::-1,:]

a 的行按逆序排列

a([1:end 1],:)

a[np.r_[:len(a),0]]

a 的末尾添加第一行数据的副本

a.'

a.transpose()a.T

a 的转置

a'

a.conj().transpose()a.conj().T

a 的共轭转置

a * b

a @ b

矩阵乘法

a .* b

a * b

元素级乘法

a./b

a/b

元素级除法

a.^3

a**3

元素级幂运算

(a > 0.5)

(a > 0.5)

一个矩阵,其第 i,j 个元素为 (a_ij > 0.5)。MATLAB 结果是一个由逻辑值 0 和 1 组成的数组。NumPy 结果是一个由布尔值 FalseTrue 组成的数组。

find(a > 0.5)

np.nonzero(a > 0.5)

查找满足(a > 0.5) 条件的索引

a(:,find(v > 0.5))

a[:,np.nonzero(v > 0.5)[0]]

提取向量 v > 0.5 的对应列在a 中的列

a(:,find(v>0.5))

a[:, v.T > 0.5]

提取列向量 v > 0.5 的对应列在a 中的列

a(a<0.5)=0

a[a < 0.5]=0

a 中小于 0.5 的元素置零

a .* (a>0.5)

a * (a > 0.5)

a 中小于 0.5 的元素置零

a(:) = 3

a[:] = 3

将所有元素设置为相同的标量值

y=x

y = x.copy()

NumPy 使用引用赋值

y=x(2,:)

y = x[1, :].copy()

NumPy 切片使用引用赋值

y=x(:)

y = x.flatten()

将数组转换为向量(注意,这会强制复制)。为了获得与 MATLAB 相同的数据顺序,请使用 x.flatten('F')

1:10

np.arange(1., 11.)np.r_[1.:11.]np.r_[1:10:10j]

创建一个递增向量(参见注释 RANGES

0:9

np.arange(10.)np.r_[:10.]np.r_[:9:10j]

创建一个递增向量(参见注释 RANGES

[1:10]'

np.arange(1.,11.)[:, np.newaxis]

创建一个列向量

zeros(3,4)

np.zeros((3, 4))

一个 3x4 的二维数组,所有元素均为 64 位浮点型零

zeros(3,4,5)

np.zeros((3, 4, 5))

一个 3x4x5 的三维数组,所有元素均为 64 位浮点型零

ones(3,4)

np.ones((3, 4))

一个 3x4 的二维数组,所有元素均为 64 位浮点型一

eye(3)

np.eye(3)

一个 3x3 的单位矩阵

diag(a)

np.diag(a)

返回二维数组 a 的对角线元素组成的向量

diag(v,0)

np.diag(v, 0)

返回一个方阵,其对角线元素为向量 v 的元素,其他元素为零

rng(42,'twister')
rand(3,4)
from numpy.random import default_rng
rng = default_rng(42)
rng.random((3, 4))

或旧版本:random.rand((3, 4))

生成一个 3x4 的随机数组,使用默认随机数生成器和种子 = 42

linspace(1,3,4)

np.linspace(1,3,4)

在 1 和 3 之间生成 4 个等间距的样本,包含 1 和 3

[x,y]=meshgrid(0:8,0:5)

np.mgrid[0:9.,0:6.]np.meshgrid(r_[0:9.],r_[0:6.])

两个二维数组:一个包含 x 值,另一个包含 y 值

ogrid[0:9.,0:6.]np.ix_(np.r_[0:9.],np.r_[0:6.]

在网格上计算函数的最佳方法

[x,y]=meshgrid([1,2,4],[2,4,5])

np.meshgrid([1,2,4],[2,4,5])

np.ix_([1,2,4],[2,4,5])

在网格上计算函数的最佳方法

repmat(a, m, n)

np.tile(a, (m, n))

创建 m 行 n 列的 a 的副本

[a b]

np.concatenate((a,b),1)np.hstack((a,b))np.column_stack((a,b))np.c_[a,b]

连接 ab 的列

[a; b]

np.concatenate((a,b))np.vstack((a,b))np.r_[a,b]

连接 ab 的行

max(max(a))

a.max()np.nanmax(a)

a 的最大元素(对于 MATLAB,ndims(a)<=2;如果存在 NaN,nanmax 将忽略这些 NaN 并返回最大值)

max(a)

a.max(0)

数组 a 每列的最大元素

max(a,[],2)

a.max(1)

数组 a 每行的最大元素

max(a,b)

np.maximum(a, b)

逐元素比较 ab,并返回每对元素中的最大值

norm(v)

np.sqrt(v @ v)np.linalg.norm(v)

向量 v 的 L2 范数

a & b

logical_and(a,b)

逐元素 AND 运算符(NumPy ufunc)参见注释 LOGICOPS

a | b

np.logical_or(a,b)

逐元素 OR 运算符(NumPy ufunc)参见注释 LOGICOPS

bitand(a,b)

a & b

按位 AND 运算符(Python 原生和 NumPy ufunc)

bitor(a,b)

a | b

按位 OR 运算符(Python 原生和 NumPy ufunc)

inv(a)

linalg.inv(a)

二维方阵 a 的逆矩阵

pinv(a)

linalg.pinv(a)

二维数组 a 的伪逆矩阵

rank(a)

np.linalg.matrix_rank(a)

二维数组 a 的矩阵秩

a\b

如果 a 为方阵,则使用 linalg.solve(a, b);否则使用 linalg.lstsq(a, b)

求解 a x = b 中的 x

b/a

改为求解 a.T x.T = b.T

求解 x a = b 中的 x

[U,S,V]=svd(a)

U, S, Vh = linalg.svd(a); V = Vh.T

a 的奇异值分解

chol(a)

linalg.cholesky(a)

二维数组的 Cholesky 分解

[V,D]=eig(a)

D,V = linalg.eig(a)

a 的特征值 \(\lambda\) 和特征向量 \(v\),其中 \(\mathbf{a} v = \lambda v\)

[V,D]=eig(a,b)

D,V = linalg.eig(a, b)

a, b 的特征值 \(\lambda\) 和特征向量 \(v\),其中 \(\mathbf{a} v = \lambda \mathbf{b} v\)

[V,D]=eigs(a,3)

D,V = eigs(a, k=3)

查找二维数组 ak=3 个最大特征值和特征向量

[Q,R]=qr(a,0)

Q,R = linalg.qr(a)

QR 分解

[L,U,P]=lu(a),其中 a==P'*L*U

P,L,U = linalg.lu(a),其中a == P@L@U

LU 分解,带部分主元(注意:P(MATLAB) == P(NumPy).transpose())

共轭梯度法

cg

共轭梯度求解器

fft(a)

np.fft.fft(a)

a 的傅里叶变换

ifft(a)

np.fft.ifft(a)

a 的傅里叶逆变换

sort(a)

np.sort(a)a.sort(axis=0)

对二维数组 a 的每一列进行排序

sort(a, 2)

np.sort(a, axis=1)a.sort(axis=1)

对二维数组 a 的每一行进行排序

[b,I]=sortrows(a,1)

I = np.argsort(a[:, 0]); b = a[I,:]

将数组 a 保存为数组 b,按第一列排序

x = Z\y

x = linalg.lstsq(Z, y)

执行形式为 \(\mathbf{Zx}=\mathbf{y}\) 的线性回归

decimate(x, q)

signal.resample(x, np.ceil(len(x)/q))

带低通滤波的降采样

unique(a)

np.unique(a)

数组 a 中唯一值的向量

squeeze(a)

a.squeeze()

去除数组 a 的单元素维度。注意,MATLAB 始终返回二维或更高维度的数组,而 NumPy 则返回 0 维或更高维度的数组。

备注#

子矩阵:可以使用索引列表通过 ix_ 命令对子矩阵进行赋值。例如,对于二维数组 a,可以执行以下操作:ind=[1, 3]; a[np.ix_(ind, ind)] += 100

帮助:没有 MATLAB 的 which 命令的直接等效命令,但是命令 help 通常会列出函数所在的文件名。Python 还具有一个 inspect 模块(执行 import inspect),它提供一个 getfile 函数,该函数通常有效。

索引:MATLAB 使用基于 1 的索引,因此序列的初始元素的索引为 1。Python 使用基于 0 的索引,因此序列的初始元素的索引为 0。混淆和争论出现是因为两者各有优缺点。基于 1 的索引与常见的语言习惯一致,其中序列的“第一个”元素的索引为 1。基于 0 的索引 简化了索引。另见 Edsger W. Dijkstra 教授的一篇文章

范围:在 MATLAB 中,0:5 可用作范围字面量和“切片”索引(在括号内);但是,在 Python 中,像 0:5 这样的构造 *只能* 用作切片索引(在方括号内)。因此,创建了有点古怪的 r_ 对象,以允许 NumPy 具有类似的简短范围构造机制。请注意,r_ 不是像函数或构造函数那样被调用,而是使用方括号进行 *索引*,这允许在参数中使用 Python 的切片语法。

逻辑运算符:NumPy 中的 &| 是按位与/或,而在 MATLAB 中 &| 是逻辑与/或。两者看起来可能相同,但存在重要区别。如果您使用了 MATLAB 的 &| 运算符,则应使用 NumPy 的 ufunc 函数 logical_and/logical_or。MATLAB 和 NumPy 的 &| 运算符之间显著的区别在于

  • 非逻辑 {0,1} 输入:NumPy 的输出是输入的按位与。MATLAB 将任何非零值视为 1 并返回逻辑与。例如,NumPy 中的 (3 & 4)0,而在 MATLAB 中,34 都被视为逻辑真,(3 & 4) 返回 1

  • 优先级:NumPy 的 & 运算符优先级高于逻辑运算符(如 <>);MATLAB 的则相反。

如果您知道您有布尔参数,则可以使用 NumPy 的按位运算符,但要小心使用括号,如下所示:z = (x > 1) & (x < 2)logical_andlogical_or 的 NumPy 运算符形式的缺失是 Python 设计的一个不幸的结果。

重塑和线性索引:MATLAB 始终允许使用标量或线性索引访问多维数组,而 NumPy 则不允许。线性索引在 MATLAB 程序中很常见,例如,矩阵上的 find() 返回它们,而 NumPy 的 find 行为不同。转换 MATLAB 代码时,可能需要先将矩阵重塑为线性序列,执行一些索引操作,然后重新塑形。由于重塑(通常)会产生指向相同存储区的视图,因此应该可以相当有效地做到这一点。请注意,NumPy 中重塑使用的扫描顺序默认为“C”顺序,而 MATLAB 使用 Fortran 顺序。如果您只是将序列转换为线性序列然后再转换回来,这无关紧要。但是,如果您正在转换依赖于扫描顺序的 MATLAB 代码中的重塑,那么此 MATLAB 代码:z = reshape(x,3,4); 在 NumPy 中应该变为 z = x.reshape(3,4,order='F').copy()

“数组”还是“矩阵”?我应该使用哪个?#

历史上,NumPy 提供了一个特殊的矩阵类型 np.matrix,它是 ndarray 的子类,它使二元运算成为线性代数运算。您可能会在一些现有代码中看到它代替 np.array 使用。那么,应该使用哪个呢?

简短答案#

使用数组.

  • 它们支持 MATLAB 中支持的多维数组代数

  • 它们是 NumPy 的标准向量/矩阵/张量类型。许多 NumPy 函数返回数组,而不是矩阵。

  • 元素级运算和线性代数运算之间有明确的区别。

  • 如果您愿意,您可以拥有标准向量或行/列向量。

在 Python 3.5 之前,使用数组类型的唯一缺点是您必须使用 dot 而不是 * 来乘以(约简)两个张量(标量积、矩阵向量乘法等)。从 Python 3.5 开始,您可以使用矩阵乘法 @ 运算符。

鉴于上述情况,我们打算最终弃用 matrix

详细答案#

NumPy 包含 array 类和 matrix 类。array 类旨在成为用于各种数值计算的通用 n 维数组,而 matrix 旨在专门促进线性代数计算。实际上,两者之间只有少数几个关键区别。

  • 运算符 *@,函数 dot()multiply()

    • 对于 array`*` 表示逐元素乘法,而`@` 表示矩阵乘法;它们对应的函数分别是 multiply()dot()。(在 Python 3.5 之前,`@` 运算符不存在,必须使用 dot() 进行矩阵乘法)。

    • 对于 matrix`*` 表示矩阵乘法,而对于逐元素乘法,则必须使用 multiply() 函数。

  • 向量(一维数组)的处理

    • 对于 array向量形状 1xN、Nx1 和 N 都是不同的。诸如 A[:,1] 这样的操作返回形状为 N 的一维数组,而不是形状为 Nx1 的二维数组。对一维 array 进行转置不会有任何作用。

    • 对于 matrix一维数组总是会被向上转换为 1xN 或 Nx1 矩阵(行向量或列向量)。A[:,1] 返回形状为 Nx1 的二维矩阵。

  • 高维数组 (ndim > 2) 的处理

    • array 对象可以具有大于 2 的维度

    • matrix 对象始终只有两个维度

  • 便捷属性

    • array 具有 .T 属性,它返回数据的转置。

    • matrix 还具有 .H、.I 和 .A 属性,它们分别返回矩阵的共轭转置、逆矩阵和 asarray()

  • 便捷构造函数

    • array 构造函数接受(嵌套的)Python 序列作为初始化器。例如,array([[1,2,3],[4,5,6]])

    • matrix 构造函数还可以接受方便的字符串初始化器。例如 matrix("[1 2 3; 4 5 6]")

两者都有优缺点

  • array

    • :) 逐元素乘法很容易:A*B

    • :( 你必须记住矩阵乘法有它自己的运算符 `@`。

    • :) 你可以将一维数组视为行向量列向量A @ vv 视为列向量,而 v @ Av 视为行向量。这可以节省你很多转置操作。

    • :) array 是 NumPy 的“默认”类型,因此它得到了最多的测试,并且是最有可能被使用 NumPy 的第三方代码返回的类型。

    • :) 它非常擅长处理任何数量维度的数。

    • :) 如果您熟悉张量代数,则其语义更接近张量代数。

    • :) 所有运算符(*/+- 等)都是逐元素的。

    • :( 来自 scipy.sparse 的稀疏矩阵与数组的交互性不太好。

  • matrix

    • :\\ 行为更类似于 MATLAB 矩阵。

    • <:( 最大二维。要保存三维数据,你需要 array 或可能是 matrix 的 Python 列表。

    • <:( 最小二维。你不能有向量。它们必须转换为单列或单行矩阵。

    • <:( 由于 array 是 NumPy 中的默认类型,因此某些函数即使你给它们提供 matrix 作为参数,也可能会返回 array。对于 NumPy 函数来说,这种情况不应该发生(如果发生,则这是一个错误),但是基于 NumPy 的第三方代码可能不会像 NumPy 那样遵守类型保持。

    • :) A*B 是矩阵乘法,所以它看起来就像你在线性代数中写的那样(对于 Python >= 3.5,普通数组使用 `@` 运算符具有相同的便利性)。

    • <:( 逐元素乘法需要调用函数 multiply(A,B)

    • <:( 运算符重载的使用有点不合理:* 不执行逐元素运算,但 / 执行。

    • scipy.sparse 的交互更简洁。

因此,强烈建议使用 array。事实上,我们打算最终弃用 matrix

自定义您的环境#

在 MATLAB 中,您可以用来自定义环境的主要工具是使用您最喜欢的函数的位置修改搜索路径。您可以将此类自定义内容放入 MATLAB 启动时将运行的启动脚本中。

NumPy,或者更确切地说 Python,具有类似的功能。

  • 要修改您的 Python 搜索路径以包含您自己的模块的位置,请定义 PYTHONPATH 环境变量。

  • 要在启动交互式 Python 解释器时执行特定的脚本文件,请将 PYTHONSTARTUP 环境变量定义为包含启动脚本名称。

与 MATLAB 不同,在 MATLAB 中,路径上的任何内容都可以立即调用,而在 Python 中,您需要首先执行“import”语句才能访问特定文件中 的函数。

例如,您可以创建一个类似这样的启动脚本(注意:这只是一個示例,不是“最佳实践”的陈述)

# Make all numpy available via shorter 'np' prefix
import numpy as np
#
# Make the SciPy linear algebra functions available as linalg.func()
# e.g. linalg.lu, linalg.eig (for general l*B@u==A@u solution)
from scipy import linalg
#
# Define a Hermitian function
def hermitian(A, **kwargs):
    return np.conj(A,**kwargs).T
# Make a shortcut for hermitian:
#    hermitian(A) --> H(A)
H = hermitian

使用已弃用的 matrix 和其他 matlib 函数

# Make all matlib functions accessible at the top level via M.func()
import numpy.matlib as M
# Make some matlib functions accessible directly at the top level via, e.g. rand(3,3)
from numpy.matlib import matrix,rand,zeros,ones,empty,eye