在MATLAB中索引向量效率低下吗?

背景

我的问题是由简单的观察所激发的,这些观察有点破坏了经验丰富的MATLAB用户常常持有/做出的信仰/假设:

  • MATLAB在内置函数和基本语言function(如索引向量和matrix)方面进行了很好的优化。
  • MATLAB中的循环很慢(尽pipe是JIT),如果algorithm可以用本地“vector化”方式表示,通常应该避免。

底线:核心的MATLABfunction是高效的,试图超越它使用MATLAB代码是困难的,即使不是不可能的。

调查向量索引的性能

下面显示的示例代码与基本相同:我为所有向量条目分配一个标量值。 首先,我分配一个空的向量x

 tic; x = zeros(1e8,1); toc Elapsed time is 0.260525 seconds. 

有了x我想将其所有条目设置为相同的值。 在实践中,你会以不同的方式做,例如, x = value*ones(1e8,1) ,但这里的要点是调查向量索引的性能。 最简单的方法是写:

 tic; x(:) = 1; toc Elapsed time is 0.094316 seconds. 

我们把它称为方法1(从分配给x的值)。 这似乎是非常快(至less比内存分配更快)。 因为我在这里做的唯一的事情就是在内存上运行,所以我可以通过计算获得的有效内存带宽并将其与计算机的硬件内存带宽进行比较来估计此代码的效率:

 eff_bandwidth = numel(x) * 8 bytes per double * 2 / time 

在上面,我乘以2因为除非使用SSEstream,存储器中的设置值要求向量既被读取也被写入存储器。 在上面的例子中:

 eff_bandwidth(1) = 1e8*8*2/0.094316 = 17 Gb/s 

我的电脑的基于STREAM标准的内存带宽约为17.9 Gb / s,确实如此 – 在这种情况下,MATLAB提供了接近峰值的性能! 到现在为止还挺好。

如果要将所有vector元素设置为某个值,则方法1是合适的。 但是,如果要访问元素的每个step条目,则需要用例如1:step:endreplace: 以下是与方法1的直接速度比较:

 tic; x(1:end) = 2; toc Elapsed time is 0.496476 seconds. 

虽然你不会期望它有任何不同,但方法2显然是很大的麻烦:因素5无故减速。 我怀疑在这种情况下,MATLAB显式地分配索引向量( 1:end )。 这是通过使用显式向量大小而不是end确认的:

 tic; x(1:1e8) = 3; toc Elapsed time is 0.482083 seconds. 

方法2和3同样不好。

另一种可能是显式创build索引向量id并使用它来索引x 。 这为您提供了最灵活的索引function。 在我们的情况下:

 tic; id = 1:1e8; % colon(1,1e8); x(id) = 4; toc Elapsed time is 1.208419 seconds. 

现在,这真的是一些东西 – 比方法1减less了12倍! 我知道它应该比方法1更差,因为用于id的额外内存,但是为什么比方法2和3差得多呢?

让我们尝试循环尝试 – 尽可能听起来毫无希望。

 tic; for i=1:numel(x) x(i) = 5; end toc Elapsed time is 0.788944 seconds. 

一个很大的惊喜 – 一个循环击败了vectorized方法4,但仍然比方法1,方法2和方法3慢。事实certificate,在这种情况下,你可以做得更好:

 tic; for i=1:1e8 x(i) = 6; end toc Elapsed time is 0.321246 seconds. 

这也许是这项研究最奇怪的结果 – 一个MATLAB写的循环显着优于本地向量索引 。 那当然不是这样的。 请注意,JIT'ed循环仍然比方法1几乎获得的理论峰值慢3倍。所以还有很大的改进空间。 这是令人惊讶的(一个更强的词会更适合)通常的“向量化”索引( 1:end )更慢。

问题

  • 在MATLAB中简单的索引非常低效(方法2,3和4比方法1慢),还是我错过了什么?
  • 为什么方法4 (太多)比方法2和3
  • 为什么使用1e8而不是numel(x)作为一个循环界限加快了代码的系数2?

编辑阅读Jonas的评论后,这里是另一种使用逻辑索引的方法:

 tic; id = logical(ones(1, 1e8)); x(id) = 7; toc Elapsed time is 0.613363 seconds. 

比方法4好得多。

为了方便:

 function test tic; x = zeros(1,1e8); toc tic; x(:) = 1; toc tic; x(1:end) = 2; toc tic; x(1:1e8) = 3; toc tic; id = 1:1e8; % colon(1,1e8); x(id) = 4; toc tic; for i=1:numel(x) x(i) = 5; end toc tic; for i=1:1e8 x(i) = 6; end toc end 

当然,我只能推测。 但是,当我运行JIT编译器启用与禁用testing时,我得到以下结果:

  % with JIT no JIT 0.1677 0.0011 %# init 0.0974 0.0936 %# #1 I added an assigment before this line to avoid issues with deferring 0.4005 0.4028 %# #2 0.4047 0.4005 %# #3 1.1160 1.1180 %# #4 0.8221 48.3239 %# #5 This is where "don't use loops in Matlab" comes from 0.3232 48.2197 %# #6 0.5464 %# logical indexing 

划分显示我们在哪里有任何速度增加:

 % withoutJit./withJit 0.0067 %# w/o JIT, the memory allocation is deferred 0.9614 %# no JIT 1.0057 %# no JIT 0.9897 %# no JIT 1.0018 %# no JIT 58.7792 %# numel 149.2010 %# no numel 

在初始化时出现明显的加速,因为在JITclosures的情况下,似乎MATLAB会延迟内存分配,直到它被使用,所以x =零(…)并没有做任何事情。 (谢谢@angainor)。

方法1到4似乎没有受益于JIT。 我猜#4可能会很慢,因为subsref额外inputtesting,以确保input是正确的forms。

numel结果可能与编译器难以处理不确定的迭代次数有关,或者由于检查循环的边界是否正常而导致一些开销(认为,没有JITtesting表明只有0.1s那)

令人惊讶的是,在我的机器上的R2012b上,逻辑索引似乎比#4慢。

我认为这再一次certificate了MathWorks在加速代码方面做了很多工作,如果您想要获得最快的执行时间,那么“不使用循环”并不总是最好的在这一刻)。 不过,我发现向量化一般来说是一个很好的方法,因为(a)JIT在更复杂的循环上失败,(b)学习vector化使得你更好地理解Matlab。

结论:如果您想要速度,请使用分析器,如果切换Matlab版本,请重新configuration。


作为参考,我使用了稍微修改过的testing函数

 function tt = speedTest tt = zeros(8,1); tic; x = zeros(1,1e8); tt(1)=toc; x(:) = 2; tic; x(:) = 1; tt(2)=toc; tic; x(1:end) = 2; tt(3)=toc; tic; x(1:1e8) = 3; tt(4)=toc; tic; id = 1:1e8; % colon(1,1e8); x(id) = 4; tt(5)=toc; tic; for i=1:numel(x) x(i) = 5; end tt(6)=toc; tic; for i=1:1e8 x(i) = 6; end tt(7)=toc; %# logical indexing tic; id = true(1e8,1)); x(id)=7; tt(8)=toc; 

我没有对所有问题的答案,但是我对方法2,3和4有一些完善的推测。

关于方法2和3.看起来似乎MATLAB为向量索引分配了内存,并用从11e8值来填充1e8 。 要了解它,让我们看看发生了什么。 默认情况下,MATLAB使用double作为其数据types。 分配索引数组需要与分配x相同的时间

 tic; x = zeros(1e8,1); toc Elapsed time is 0.260525 seconds. 

目前,索引数组只包含零。 如方法1中那样,以最佳方式将值赋值给x向量需要0.094316秒。 现在,索引vector必须从内存中读取,以便它可以用于索引。 这是额外的0.094316/2秒。 回想一下,在x(:)=1向量x必须同时从内存中读取和写入。 所以只能读一半的时间 假设这是在x(1:end)=value ,那么方法2和3的总时间应该是

 t = 0.260525+0.094316+0.094316/2 = 0.402 

这几乎是正确的,但不完全。 我只能推测,但用值填充索引向量可能是作为一个额外的步骤,并需要额外的0.094316秒。 因此, t=0.4963 ,或多或less符合方法2和3的时间。

这些只是猜测,但他们似乎确认,MATLAB 明确创build索引向量时做本地向量索引。 就个人而言,我认为这是一个性能错误。 MATLAB的JIT编译器应该足够聪明,以了解这个微不足道的构造,并将其转换为调用正确的内部函数。 就目前而言,在今天的内存带宽上,有限的体系结构索引在理论峰值的20%处执行。

所以,如果你关心性能,你将不得不将x(1:step:end)作为MEX函数来实现

 set_value(x, 1, step, 1e8, value); 

现在这在MATLAB中显然是非法的,因为你不能在现场修改MEX文件中的数组。

编辑关于方法4,可以尝试分析各个步骤的性能如下:

 tic; id = 1:1e8; % colon(1,1e8); toc tic x(id) = 4; toc Elapsed time is 0.475243 seconds. Elapsed time is 0.763450 seconds. 

第一步,分配和填充索引向量的值与单独的方法2和3相同。 看起来太多了 – 分配内存和设置值( 0.260525s+0.094316s = 0.3548s )最多需要的时间,所以在某个地方有额外的开销0.12秒,其中I不能理解。 第二部分( x(id) = 4 )看起来效率也很低:需要花费时间设置x的值,并读取id向量( 0.094316s+0.094316/2s = 0.1415s )加上一些错误检查在id值上。 用C编程,这两个步骤是:

 create id 0.214259 x(id) = 4 0.219768 

使用的代码检查double索引实际上表示一个整数,并且它符合x的大小:

 tic(); id = malloc(sizeof(double)*n); for(i=0; i<n; i++) id[i] = i; toc("create id"); tic(); for(i=0; i<n; i++) { long iid = (long)id[i]; if(iid>=0 && iid<n && (double)iid==id[i]){ x[iid] = 4; } else break; } toc("x(id) = 4"); 

第二步需要超过预期的0.1415s – 这是由于对id值进行错误检查的必要性。 开销似乎对我来说太大了 – 也许它可以写得更好。 但是,所需的时间是0.4340s ,而不是1.208419s 。 MATLAB做了什么 – 我不知道。 也许有必要这样做,我只是没有看到它。

当然,使用doubles作为索引引入了两个额外的开销水平:

  • 尺寸是uint32尺寸的两倍。 回想一下,内存带宽是这里的限制因素。
  • 双打需要转换成整数索引

方法4可以使用整数索引在MATLAB中编写:

 tic; id = uint32(1):1e8; toc tic x(id) = 8; toc Elapsed time is 0.327704 seconds. Elapsed time is 0.561121 seconds. 

这明显提高了30%的性能,并certificate应该用整数作为向量指标。 但是,开销仍然存在。

正如我现在所看到的,我们无法做任何事情来改善在MATLAB框架内工作的情况,我们必须等到Mathworks解决这些问题。