默认情况下,所有用于计算相关性或协方差的内置函数都返回一个矩阵。我正试图写一个有效的函数来计算种子区域和其他各种区域之间的相关性,但我不需要其他区域之间的关联。因此,我认为计算全相关矩阵将是低效的。
相反,我可以计算每个区域和种子区域之间的相关矩阵,选择一个非对角点并存储它,但我觉得在这种情况下循环也是低效的。
更具体地说,我的三维空间中的每个点都有一个时间维度。我正试图计算给定点和给定半径内空间中所有点之间的平均相关性。我想重复这个过程数十万次,对于许多不同的半径长度,等等,所以我希望这个过程尽可能有效。
那么,在不计算我会忽略的相关性的情况下,计算单个向量和其他几个向量之间的相关性的最佳方法是什么?
谢谢,Chris
编辑:这是我的代码。。。
function [corrMap] = TIME_meanCorrMap(A,radius)
% Even though the variable is "radius", we work with cubes for simplicity...
% So, the radius is the distance (in voxels) from the center of the cube an edge.
denom = ((radius*2)^3)-1;
dim = size(A);
corrMap = zeros(dim(1:3));
for x = radius+1:dim(1)-radius
rx = [x-radius : x+radius];
for y = radius+1:dim(2)-radius
ry = [y-radius : y+radius];
for z = radius+1:dim(3)-radius
rz = [z-radius : z+radius];
corrCoefs = zeros(1,denom);
seed = A(x,y,z,:);
i=0;
for xx = rx
for yy = ry
for zz = rz
if ~all([x y z] == [xx yy zz])
i = i + 1;
temp = corrcoef(seed,A(xx,yy,zz,:));
corrCoeffs(i) = temp(1,2);
end
end
end
end
corrMap = mean(corrCoeffs);
end
end
end
编辑:以下是一些补充公认答案的更多时间。使用bsxfun()进行归一化,并使用矩阵乘法计算相关性:
tic; for i=1:10000
x=rand(100);
xz = bsxfun(@rdivide,bsxfun(@minus,x,mean(x)),std(x));
cc = xz(:,2:end)' * xz(:,1) ./ 99;
end; toc
Elapsed time is 6.928251 seconds.
使用zscore()进行归一化,矩阵乘法计算相关性:
tic; for i=1:10000
x=rand(100);
xz = zscore(x);
cc = xz(:,2:end)' * xz(:,1) ./ 99;
end; toc
Elapsed time is 7.040677 seconds.
使用bsxfun()进行归一化,使用corr()计算相关性。
tic; for i=1:10000
x=rand(100);
xz = bsxfun(@rdivide,bsxfun(@minus,x,mean(x)),std(x));
cc = corr(x(:,1),x(:,2:end));
end; toc
Elapsed time is 11.385707 seconds.
当然可以改进您当前使用的for
循环。如果您有足够的RAM,可以使用矩阵乘法对相关计算进行并行化。但是,这将要求您将4维数据矩阵A展开为不同的形状。最有可能的是,你正在处理三维体素fMRI数据,在这种情况下,你必须将[x y z time]矩阵重塑为[index time]矩阵。我认为你可以处理这种重塑。一旦您准备好了seed
时间进程[Time by 1]和target
时间进程[Timeby NumTargets],您就可以执行一些更高效的计算。
有效计算所需相关性的一种快速方法是在MATLAB中使用corr
函数。该函数将接受2个矩阵自变量,并且它将非常有效地计算自变量1的列和自变量2的列之间的所有成对相关性,例如
T = 200; %time samples
N = 20; %number of other voxels
seed = randn(T,1); %data from seed voxel
targets = randn(T,N); %data from target voxels
%here is the for loop method
tic
for n = 1:N
tmp = corrcoef(seed, targets(:,n));
tmpcc = tmp(1,2);
end
looptime = toc;
%here is the parallel method
tic
cc = corr(seed, targets);
matrixtime = toc;
在我的机器上,corr
中的并行运算比循环法快一个与T*N成比例的因子。
如果你愿意自己完成底层的矩阵运算,那么可以比corr
函数快一点,而且在任何情况下都值得知道它们是什么。两个向量之间的相关性基本上是一个标准化的点积,因此使用上面的约定,您可以用以下方式计算相关性
zseed = zscore(seed); %normalize the seed timecourse by z-scoring
ztargets= zscore(targets); %normalize the target timecourses by z-scoring
ztargets = ztargets'; %flip columns and rows for convenience
cc2 = ztargets*zseed./(T-1); %compute many dot products with one matrix multiplication
上面的代码基本上就是corr
函数要做的,这就是为什么它比循环快得多。请注意,大多数操作时间都在zscore
操作中,如果使用bsxfun
命令有效地计算zscore
,则可以提高corr
函数的性能。目前,我希望这能为您提供一些方向,让您了解如何计算种子时间进程和许多目标时间进程之间的相关性,而不必循环并分别计算每一个时间进程。