寻找序列中的零岛

Question

想象一下你有一个很长的序列。查找序列全为零的间隔（或更准确地说，序列下降到接近零值 abs(X)）的最有效方法是什么：

为简单起见，我们假设以下序列：

sig = [1 1 0 0 0 0 1 1 1 1 1 0 1 0 0 0 1 1 1 1 1 1 1 1 0 0 1 1 1 0];

我试图获取以下信息：

startIndex   EndIndex    Duration
3            6           4
12           12          1
14           16          3
25           26          2
30           30          1

然后使用此信息，我们找到持续时间 >= 某个指定值（例如 3）的间隔，并返回值的索引在所有这些间隔的总和中：

indices = [3 4 5 6 14 15 16];

最后一部分与之前的问题相关：

MATLAB：矢量化数组创建 来自开始/结束索引列表

这是我到目前为止所拥有的：

sig = [1 1 0 0 0 0 1 1 1 1 1 0 1 0 0 0 1 1 1 1 1 1 1 1 0 0 1 1 1 0];
len = length(sig);
thresh = 3;

%# align the signal with itself successively shifted by one
%# v will thus contain 1 in the starting locations of the zero interval
v = true(1,len-thresh+1);
for i=1:thresh
    v = v & ( sig(i:len-thresh+i) == 0 );
end

%# extend the 1's till the end of the intervals
for i=1:thresh-1
    v(find(v)+1) = true;
end

%# get the final indices
v = find(v);

我正在寻求矢量化/优化代码，但我对其他解决方案持开放态度。 我必须强调，空间和时间效率非常重要，因为我正在处理大量的长生物信号。

Answer 1

以下是我以向量化方式解决您的问题所采取的步骤，从给定的向量 sig 开始：

• 首先，对向量进行阈值处理以获得向量 tsig零和一（信号绝对值足够接近零的零，其他地方的零）：

  tsig = (abs(sig) >= eps); %# 使用 eps 作为阈值
  

• 接下来，使用函数 差异 和查找：

  dsig = diff([1 tsig 1]);
  startIndex = find(dsig < 0);
  endIndex = find(dsig > 0)-1;
  持续时间=结束索引-开始索引+1；
  

• 然后，查找持续时间大于或等于某个值（例如 3，在您的示例中）的零字符串：

  stringIndex = (持续时间 >= 3);
  开始索引 = 开始索引(字符串索引);
  endIndex = endIndex(字符串索引);
  

• 最后，使用 我对链接问题的回答生成的方法您的最终索引集：

  索引 = Zeros(1,max(endIndex)+1);
  索引（起始索引）= 1；
  索引(endIndex+1) = 索引(endIndex+1)-1;
  索引 = find(cumsum(索引));
  

Answer 2

您可以通过查找长度为 thresh 的零字符串来解决此问题作为字符串搜索任务（STRFIND 函数非常快）

startIndex = strfind(sig, zeros(1,thresh));

请注意，较长的子字符串将在多个位置进行标记，但一旦我们添加，最终将被连接从 startIndex 开始到 start+thresh-1 结束的间隔的中间位置。

indices = unique( bsxfun(@plus, startIndex', 0:thresh-1) )';

请注意，您始终可以将最后一步与 @gnovice 的 CUMSUM/FIND 解决方案进行交换，来自 链接问题。

Answer 3

function indice=sigvec(sig,thresh)
    %extend sig head and tail to avoid 0 head and 0 tail

    exsig=[1,sig,1];
    %convolution sig with extend sig
    cvexsig=conv(exsig,ones(1,thresh));
    tempsig=double(cvexsig==0);

    indice=find(conv(tempsig,ones(1,thresh)))-thresh;

Answer 4

genovice 的上述答案可以修改为查找向量中非零元素的索引，如下所示：

    tsig = (abs(sig) >= eps);
    dsig = diff([0 tsig 0]);
    startIndex = find(dsig > 0);
    endIndex = find(dsig < 0)-1;
    duration = endIndex-startIndex+1;

Answer 5

我认为最 MATLAB/“矢量化”的方法是使用 [-1 1] 等滤波器计算信号的卷积。您应该查看函数 conv 的文档。然后在 conv 的输出上使用 find 来获取相关索引。

Answer 6

正如 gnovice 所示，我们将进行阈值测试，使“接近零”真正为零：

logcl = abs(sig(:)) >= zero_tolerance;

然后找到累积和不增加的区域：

cs = cumsum(logcl);
islands = cs(1+thresh:end) == cs(1:end-thresh);

记住 gnovice 填写索引范围的好方法

v = Zeros(1,max(endInd)+1); %# 一个零数组
v(起始位置) = 1; %# 将 1 放在间隔的开头
v(endInd+1) = v(endInd+1)-1; %# 在间隔结束后添加 -1 一个索引
索引 = find(cumsum(v)); %# 执行累积和并找到非零条目

我们注意到，我们的 islands 向量已经在 startInd 位置中包含了一个，并且出于我们的目的，endInd 始终是 thresh稍后（较长的运行有岛屿中的运行）

endcap = zeros(thresh,1);
indices = find(cumsum([islands ; endcap] - [endcap ; islands]))

测试

sig = [1 1 0 0 0 0 1 1 1 1 1 0 1 0 0 0 1 1 1 1 1 1 1 1 0 0 1 1 1 0];
logcl = abs(sig(:)) >= .1;
cs = cumsum(logcl);
islands = cs(1+thresh:end) == cs(1:end-thresh);
endcap = zeros(thresh,1);
indices = find(cumsum([islands ; endcap] - [endcap ; islands]))

索引 =

     2
     3
     4
     5
    13
    14
    15