在两个条件下移动 window - Matlab

Question

我有一个时间序列数据，值的幅度范围从 1 到 100。我需要使用移动 window（持续 10 秒）和 flag/highlight 数据段在 91 到 97 的阈值范围内下降至少 10 秒或更长时间。可用数据如下link

https://www.dropbox.com/s/fx7z9qzg8gxb4x3/data.mat?dl=0

非常感谢这方面的任何帮助。

我使用了以下代码：

`x= load('data.mat');
time = 1:length(x);
window = 1000; % 10 second window - data was sampled at 100Hz.
[TF,L,U,C] = isoutlier(x); % I tried to find outlier
figure;
plot(time,x,time(TF),x(TF),'*')
legend('Data','Outlier')
figure;
plot(time,x,time(TF),x(TF),'x',time,93*ones(1,length(x)),time,
97*ones(1,length(x)))`

得到下图[https://i.stack.imgur.com/fZa9m.jpg] 但不确定。如何将其用作阈值 windows.

提前致谢

Answer 1

可能有一种方法可以通过实际的滑动来做到这一点 window，但是由于我不得不以不同的方式对您执行类似的操作，所以我将重用几年前使用的函数这非常适合用另一种方式解决这个问题。

基本上流程是：

找到所有验证条件 (x>91 & x<97) 的索引。
找到在上述步骤中找到的每个连续 group/interval 的开始和停止索引。
计算每个间隔的长度（上面给出的开始和停止很容易）
丢弃间隔太短

现在你有一个间隔列表（开始和停止索引）验证你的所有条件（信号水平和持续时间）。有了它，您就可以：

重建有效索引的单个向量
创建一个仅包含标记数据的辅助向量（x 的数据验证了所有条件）。
显示 :-)

在代码中它看起来像这样：

%% Your inputs
load('data.mat');
time = 1:length(x);

%% Setup
winSize = 1000 ;    % Fs * Duration = 100 * 10 = 1000 points
lvlmin = 91 ;       % Minimum level to flag
lvlmax = 97 ;       % Maximum level to flag

%% Find which interval to flag
% find all the indices where the condition is true
idx = ( x>lvlmin ) & ( x<lvlmax ) ;
% get the start and stop index of each group of consecutive indices
itBounds = get_interval_boundaries( idx ) ;
% get the length of each interval/group
itLenght = diff(itBounds,1,2)+1 ; 
% only consider intervals >= winSize
it2flag = itLenght >= winSize ; 
nint = sum(it2flag) ;       % found 241 valid intervals out of 596.

%% Clear [idx] of the short intervals
itbad = itBounds( ~it2flag , : ) ; % keep only the intervals to discard
for k=1:size(itbad,1)
    idx(itbad(k,1):itbad(k,2)) = false ;
end

%% Display
flaggedTime = time(idx) ;
flaggedData = x(idx) ;

figure
plot(time,x)
hold on
plot(flaggedTime,flaggedData,'.')

lx = [time(1) time(end)] ;
plot( lx , [lvlmin lvlmin], '-.k')
plot( lx , [lvlmax lvlmax], '-.k')

%% OR, alternatively, keep vectors the same lenght by adding NaNs
flaggedData = x ;
flaggedData(~idx) = NaN ;

figure
plot(time,x)
hold on
plot(time,flaggedData)

以及如何标记数据的预览：

您需要 get_interval_boundaries.m 的代码。我本可以只用更少的代码编写所需的功能，但由于这是可用的并且可以完美运行，因此无需重新发明轮子：

function itbound = get_interval_boundaries(vec)
% function itbound = get_interval_boundaries(vec)
%
% This function takes a vector of index as input (or a logical index array)
% It returns a nx2 table containing on each line the first and last index
% of consecutive intervals of indexes.
% ex:
% A = [1 2 3 5 7 8 9] ;
% [itbound] = get_interval_boundaries(A)
% itbound =
%      1     3
%      5     5
%      7     9
% 
% 09-Oct-2011 - Hoki:  creation
% 15-Sep-2012 - Hoki: Corrected last index special case
%                       (return last idx instead of 0)
% 01-Sep-2014 - Hoki: Corrected first index special case
%                       (return [1 1] in case vec is a scalar)

%% Check vec type (logical array or direct indexing)

% Return empty vector if input is empty
itbound = [] ;
if isempty(vec)
    return
end

% Check the type of input vector
if islogical(vec)
    idxDirect = find(vec) ;

elseif isnumeric(vec)
    idxDirect = vec ;
else
    errordlg('bad type for ''vec''. Variable should be numeric or logical',mfilename,'modal')
    return
end

%% Detect intervals
Npts = length(idxDirect) ;

% return [] in case vec is all [0]
if Npts == 0 ; return ; end

itbound(1,1) = idxDirect(1) ;
% return [x x] in case vec is a scalar value [x]
if Npts == 1
    itbound(1,2) = idxDirect(1) ;
    return
end

j=1 ;
for k = 2:Npts
    if idxDirect(k)==idxDirect(k-1)+1
        if k~=Npts
            % Cycle the loop
            continue
        else
            % Last point: Assign closing boundary of last interval
            itbound(j,2) = idxDirect(k) ;
        end
    else
        % Assign closing boundary of current interval
        itbound(j,2) = idxDirect(k-1) ;
        % Assign opening boundary of next interval
        j = j + 1 ;
        itbound(j,1) = idxDirect(k) ;
        % If we're on the very last index, close the interval.
        if k==Npts
            itbound(j,2) = idxDirect(k) ;
        end
    end
end

在两个条件下移动 window - Matlab

Moving window with two conditions - Matlab

matlab

window

highlight

threshold