我有两个向量
data vector: A = [1 2 2 1 2 6; 2 3 2 3 3 5]
label vector: B = [1 2 1 2 3 NaN]
我想拿以具有相同标签的所有列的平均值,并将其输出为矩阵,该矩阵由标签编号排序,忽略NAN。因此,在此示例中,我想要:
labelmean(A,B) = [1.5 1.5 2; 2 3 3]
这可以用这样的前面进行。
function out = labelmean(data,label)
out=[];
for i=unique(label)
if isnan(i); continue; end
out = [out, mean(data(:,label==i),2)];
end
但是,我正在处理包含许多数据点和标签的大型阵列。此外,此代码段将经常执行。我想知道是否有更有效的方法可以在不循环浏览每个单独标签的情况下进行此操作。
这是一种方法:
- 获取不包含
NaN
s的标签的索引。 - 创建一个零零矩阵和乘以
A
的零矩阵将提供所需的行总和。 - 将矩阵除以每列的总和,以使总和变为平均值。
- 应用矩阵乘法以获取结果,然后转换为完整的矩阵。
代码:
I = find(~isnan(B)); % step 1
t = sparse(I, B(I), 1, size(A,2), max(B(I))); % step 2
t = bsxfun(@rdivide, t, sum(t,1)); % step 3
result = full(A*t); % step 4
这将是使用accumarray
的好情况。将accumarray
视为微型MapReduce范式。有键和值,因此accumarray
的作业是将共享相同键的所有值分组在一起,然后对这些值进行一些操作。就您而言,键将是B
中的元素,但是对于B
中相应值所需的值是行的值。基本上,对于B
中的每个值,B
中的位置告诉您需要在A
中访问哪一行。因此,我们只需要抓住所有映射到同一ID的行位置,访问A
的行,然后在所有行上找到平均值即可。我们需要小心,因为我们忽略了NaN
的值。我们可以在调用accumarray
之前过滤掉。传统上您在accumarray
中所做的"某件事"应该输出一个数字,但实际上我们正在为每个标签输出一个列向量。因此,诀窍是将输出包裹到单元格数组中,然后使用与逗号分隔列表结合的cat
将输出转换为矩阵。
这样,类似的事情应该起作用:
% Sample data
A = [1 2 2 1 2 6; 2 3 2 3 3 5];
B = [1 2 1 2 3 NaN];
% Find non-NaN locations
mask = ~isnan(B);
% Generate row locations that are not NaN as well as the labels
ind = 1 : numel(B);
Bf = B(mask).';
ind = ind(mask).';
% Find label-wise means
C = accumarray(Bf, ind, [], @(x) {mean(A(:,x), 2)});
% Convert to numeric matrix
out = cat(2, C{:});
如果您不喜欢使用临时变量来查找这些非NaN
值,我们可以以更少的代码行进行此操作,但是您仍然需要行索引向量来确定我们需要在哪里进行采样来自:
% Sample data
A = [1 2 2 1 2 6; 2 3 2 3 3 5];
B = [1 2 1 2 3 NaN];
% Solution
ind = 1 : numel(B);
C = accumarray(B(~isnan(B)).', ind(~isnan(B)).', [], @(x) {mean(A(:,x), 2)});
out = cat(2, C{:});
使用您的数据,我们得到:
>> out
out =
1.5000 1.5000 2.0000
2.0000 3.0000 3.0000
这个答案不是一种新方法,而是给定答案的基准,因为如果您谈论性能,您始终必须基准测试。
clear all;
% I tried to make a real-life dataset (the original author may provide a
% better one)
A = [1:3e4; 1:10:3e5; 1:100:3e6]; % large dataset
B = repmat(1:1e3, 1, 3e1); % large number of labels
labelmean(A,B);
labelmeanLuisMendoA(A,B);
labelmeanLuisMendoB(A,B);
labelmeanRayryeng(A,B);
function out = labelmean(data,label)
tic
out=[];
for i=unique(label)
if isnan(i); continue; end
out = [out, mean(data(:,label==i),2)];
end
toc
end
function out = labelmeanLuisMendoA(A,B)
tic
B2 = B(~isnan(B)); % remove NaN's
t = full(sparse(1:numel(B2),B2,1,size(A,2),max(B2))); % template matrix
out = A*t; % sum of columns that share a label
out = bsxfun(@rdivide, out, sum(t,1)); % convert sum into mean
toc
end
function out = labelmeanLuisMendoB(A,B)
tic
B2 = B(~isnan(B)); % step 1
t = sparse(1:numel(B2), B2, 1, size(A,2), max(B2)); % step 2
t = bsxfun(@rdivide, t, sum(t,1)); % step 3
out = full(A*t); % step 4
toc
end
function out = labelmeanRayryeng(A,B)
tic
ind = 1 : numel(B);
C = accumarray(B(~isnan(B)).', ind(~isnan(B)).', [], @(x) {mean(A(:,x), 2)});
out = cat(2, C{:});
toc
end
输出为:
Elapsed time is 0.080415 seconds. % original
Elapsed time is 0.088427 seconds. % LuisMendo original answer
Elapsed time is 0.004223 seconds. % LuisMendo optimised version
Elapsed time is 0.037347 seconds. % rayryeng answer
对于此数据集Luismendo优化版本是明显的赢家,而他的第一个版本比原始版本慢。
=>不要忘记基于您的性能!
编辑:测试平台规格
- MATLAB R2016B
- Ubuntu 64位
- 15.6 Gib Ram
- Intel®Core™i7-5600U CPU @ 2.60GHz×4