首先推荐使用matlab 2006a版本,该版本优点很多(不过有一个小bug,就是通过GUI自动生成的m文件居然一大堆warning,希望在已经发布了的2006b版本中有改善),其中对于编程人员来说比较突出的一个就是编辑窗口的自动语法检查功能。这可以在一定程度上避免使用没有被定义或赋值的变量,另外,也可以帮助你优化代码,【例1】的【方案3】就是因为我看到matlab编辑窗口的warning而得到的启发。顺便提一下,虽然matlab不像其他语言那样,对变量采用“先定义,后使用”的规则,但是,我的经验是,在使用一个变量之前,最好先对它进行“定义”,这里的“定义”是指为它分配空间,这样不但可以提高运行的速度(这在matlab的帮助中也提到,详见Preallocating Arrays一节),而且还可以减少出错的几率,特别是在循环赋值、且变量大小不固定的时候(对此可参阅这个帖子:http://vib.hit.edu.cn/forum/thread-23732-1-9.html)。
下面说说如何对matlab提速的问题,我会使用两个例子来说明。
【例1】任务描述:根据A的取值使用imshow函数显示矩阵B
A = randn(100, 100);
B = zeros(size(A));
【方案1】
[X,Y] = find(A > 0.6);
For i = 1:length(X)
B(X(i),Y(i)) = 1;
End
【方案2】
B = zeros(size(A));
X = find(A > 0.6);
B(X) = 1;
【方案3】:
B = zeros(size(A));
X = logical(A > 0.6);
B(X) = 1;
事实上,【方案1】到【方案2】的改进在“再谈Matlab的多维数组问题”一文中已经提及过,但是没想到自己在矩阵输入的时候注意到了,但是在矩阵输出的时候就忘记了,看来程序是需要不断修改、优化的,技巧也是需要不断巩固的。至于【方案3】,是得益于matlab的warning提示。
然而,这并不是表示所有类似的地方都可以用logical代替find,当遇到循环次数与X有关时,用find会更有效,对此可以参考【例2】的【方案2】,这是用logical实现不到的。
【例2】任务描述:有一个四维矩阵A(大小为61*73*61*210),其中前三维表示一个包含大脑结构的立方体,最后一维表示大脑中每个点对应的一个长度为210的时间序列。另外有一个三维矩阵Mask(大小是61*73*61),B是二值的,其中1表示该点是前景点(大脑),0表示该点是背景点。任务是对A中属于前景点的时间序列进行EMD处理,从而判断该前景点是否属于激活区。
显然,这个问题要利用循环来完成。在本人写的“再谈Matlab的多维数组问题”一文中已经提到可以把多维数组转化为一维数组来处理,在这里,也要利用这个思想。而且,从下文可以看到,正是因为 “多维变一维”的出现,才令程序得到更进一步的提速。
首先利用reshape函数把四维矩阵A变成二维矩阵B,把三维矩阵Mask变成一维矩阵C:
B = reshape(A, 61*73*61, 210);
C = Mask(:);
t = 1:210;
【方案1】
iTotalVoxel = 61*73*61;
for k = 1:iTotalVoxel
if C(k) == 1
temp = B(k,: );
imf = emd(t, temp);
…
end
end
【方案2】
D = find(C);
iTotalBrainVoxel = length(D);
for k = 1: iTotalBrainVoxel
temp = B(D(k),: );
imf = emd(t, temp);
…
end
【方案1】明显是基于C语言的套路,而【方案2】则充分避免了matlab的弱点――循环,经过改进以后(“多维转一维”为此提供了保证,多维的话,恐怕要使用形如A(X(k),Y(k),Z(k),:)的形式了),由于循环次数的降低(大约降低为原来的1/3),故运行时间大致上减少了一半。
由此可见,在matlab中,想加快运行速度,不但要减少循环的层数,而且,还要减少循环的次数。
以下是对【例2】(实现大脑激活区检测)的不同实现的结果比较:(核心工作是对73000个前景点相应的时间序列进行EMD处理,生成10个左右的imf,其中抽取每个imf时平均迭代大概10次左右)
版本1:完全用matlab编写——运行时间大约是200分钟,也就是说,要在matlab做73000*100次循环,最头痛的是迭代时的需要产生样条包络,默认的spline函数耗时相当多;
版本2:用matcom转换成cpp文件后在Borland C++ Builder中运行,完全脱离matlab——由于spline函数耗时太多,因此转换前改用了折线包络,而非样条包络,运行时间为15分钟左右,不过这个结果是对仿真数据,非实际数据而言的,因为我仍未解决matrix.h和matlib.h不能共存的问题,故无法对实际数据进行测试,不过一般来说实际数据比仿真数据运行速度更慢;
版本3:把核心代码(基于样条包络的EMD算法)做成dll文件后在matlab中调用——整个程序,从数据输入到数据输出,只在一个地方使用了一层for循环(就是【例2】中【方案2】的循环),且结合上述优化方案,对于实际数据大概5分钟就能得出结果。
小结:要学好matlab,有效地使用matlab,一定要摆脱C++的思想。能不用循环的地方,尽量不要使用(例如求极值点等这些算法基本上可以不使用循环便可实现),逐渐抛弃C++的逐点运算的思想,多从矩阵的整体(或分块)上考虑。虽然matlab 2006a版对循环已经改进不少,但是,循环的确是造成程序运行速度降低的主要原因。matlab提供的远远不止在调用函数上的方便(例如在C++中编写fft、dwt等函数可能需要几十甚至几百行,而在matlab中只需要一两个语句),运行速度慢或许是没有使用好它,让它发挥出所长所致的。想matlab更高效地为你服务,那就需要不断修改、优化你的代码吧(我的程序编写大概用了一个星期,而修改、优化的时间就用了两个多星期,呵呵)。最后,套用某人的一句话来作结:与其抱怨matlab运行速度慢,不如先改进一下你的算法吧。