Error and anomaly detection for intra-participant time-series data

Figures & data

Figure 1. Ankle plantar-dorsi flexion angle for one participant on a standard treadmill running at 3.35 m/s: (a) before, (b), after stage 1 (b; α1 = 0.0001) and (c) after stage 2 (α2 = 0.01; b = 1) of a two-stage outlier removal method. The data were the right ankle from heel strike (0%) to heel strike (100%). The stance phase is approximately 0–40%. A potential spatial (*) and spatial–temporal (**) outliers are indicated, which are removed, and the number of cycles at each stage are indicated (k). The cycles that were deleted are listed here for stage 1 (3, 4, 5, 10, 26, 27, 41, 42) and stage 2 (11), which correspond to the column numbers in the supplementary material used to create this figure, to facilitate comparison with other outlier detection methods.

Table 1. Reduction in number of cycles, variability and normal distribution of lower limb angle data before (Raw) and after a two-stage time-series outlier removal method (stage 1 and stage 2).

Stage	Settings			Cycles (count)		SD (°)		Normal (count)
	α1	α2	b	Mean		Mean		Mean
Raw				38.2		11.9		22.6
				Reduction	Sig	Reduction	Sig	Reduction	Sig
Stage 1	0.01			12.3	Y	10.3	N	16.0	Y
	0.001			6.0	Y	10.0	N	13.0	Y
	0.0001			3.5	Y	9.7	N	10.6	Y
	3 sets				3		0		3
Stage 1 to 2	0.01	0.01	0	0.5	Y	0.2	Y	0.4	N
	0.001	0.001	1	<0.1	N	<0.1	Y	0.1	N
	0.0001	0.01	1	2.6	Y	0.4	Y	1.3	Y
	0.0001	0.0001	3	0	N	<0.1	N	0	N
	36 sets				17		22		5
Raw to 2	36 sets				36		0		36

Notes: The data are for six participants × 3 conditions × 2 sides × 2 angles; stage 1 has three sets of statistical significance levels (α1 = 0.01, 0.001 and 0.0001), all presented; stage 2 has 36 sets of α1, α2 (0.01, 0.001 and 0.0001) and moving-window size (b = 0, 1, 2, and 3) but only four sets covering the least to most stringent combinations are presented; statistical significance (Sig) indicated as no (N; p > 0.05) or yes (Y; p < 0.05), and the numbers indicate the total of the three sets (stage 1) or 36 combinations of sets (stage 1 to stage 2, and raw to stage 2) that were statistically significant.

Table 2. Reduction in number of cycles, variability and normal distribution of knee vertical displacement data before (Raw) and after a two-stage time-series outlier removal method (stage 1 and stage 2).

Stage	Settings			Cycles (count)		SD (°)		Normal (count)
	α1	α2	b	Mean		Mean		Mean
Raw				38.4		0.01		19.1
				Reduction	Sig	Reduction	Sig	Reduction	Sig
Stage 1	0.01			11.3	Y	0.003	Y	13.9	Y
	0.001			5.3	Y	0.002	Y	10.4	Y
	0.0001			3.0	Y	0.002	Y	8.1	Y
	3 sets				3		3		3
Stage 1 to 2	0.01	0.01	0	0.3	Y	<0.001	Y	0.1	N
	0.001	0.001	1	<0.1	N	<0.001	N	0	N
	0.0001	0.01	1	2.2	Y	<0.001	Y	1.3	Y
	0.0001	0.0001	3	0	N	0	N	0	N
	36 sets				12		12		6
Raw to 2	36 sets				36		36		36

Notes: These data are for six participants × 3 conditions × 2 sides; stage 1 has three sets of statistical significance levels (α1  =  0.01, 0.001 and 0.0001), all presented; stage 2 has 36 sets of α1, α2 (0.01, 0.001 and 0.0001) and moving-window size (b  =  0, 1, 2 and 3) but only four sets covering the least to most stringent combinations are presented; statistical significance (Sig) indicated as no (N; p > 0.05) or yes (Y; p  <  0.05), and the numbers indicate the total of the three sets (stage 1) or 36 combinations of sets (stage 1 to stage 2, and raw to stage 2) that were statistically significant.

Table

function [x2D,trialskept1,trialskept2]=ts2Doutlier(xraw,b,alpha1,alpha2)

% input: xraw data over time in rows (all equal length); k trials in columns

% input: b is moving-window size, e.g. 0 or 1 or 2 or 3

% input: alpha1,alpha2 p values for 1D,2D stages (e.g. 0.0001,0.01)

% output: x2D is x with 1D and 2D outliers removed

% output: trialskept1 & 2 is list of trial numbers kept after stages 1 & both 1and 2

% toolboxes required: image (padarray; std2); statistics (mad; tinv)

% Mullineaux DR, Irwin G. 2017. Error and anomaly detection for intra-participant time-series data. % Int Biomech. 4:28–35. doi:https://doi.org/10.1080/23335432.2017.1348913

[n,k]=size(xraw); %n number time points; k number of trials

trialskept2 = 1:k; %numbered list of k trials, to keep record of trials kept

x1D = xraw; %create x after 1D outliers removed (temporary as copy of xraw)

%==========================================================================

% 1D outlier: identify and remove columns/trials with any 1D outliers

t1 = abs(tinv(alpha1/2,k-1)); %2-tailed t-statistic for alpha1 and k-1 df

for p = 1:n % loop through each point p

xloop = x1D(p,:); % temporary x for the loop

medianpj = nanmedian(xloop); % median at point p for all j = 1:k trials

MADCI = t1*1.4826*mad(xloop,1); %MAD CI (1.4826 makes MAD~ = 1SD)

lower = medianpj-MADCI; % lower and upper limits

upper = medianpj+MADCI; % below find outliers outside limits …

x1D(p,xloop<=lower | xloop>=upper)=NaN; % … and replace with NaN

end

nanCt = sum(isnan(x1D),1); % sum NaNs in each column

trialskept2(:,nanCt>0)=[]; %reduce list of trials to trial numbers kept

trialskept1 = trialskept2; %reduce list of trials to trial numbers kept

x1D(:,nanCt>0)=[]; % removes columns with any NaNs (outliers) to leave x1D

%==========================================================================

% 2D outlier detection; start with padding data by b samples at each end

xpad=padarray(x1D,b,’symmetric’); %temporary x padded using reflection

[nb,k2]=size(xpad); %rows is pb = n + 2*b; columns is k2<=k (trials for x1D)

% ———- % detrend data (by removing mean cycle)

xpadmean=mean(xpad,2); % create mean cycle

xpadmean=repmat(xpadmean,1,k2); % make matrix same size as xpad

xpaddetrended=xpad-xpadmean; % remove mean cycle from all cycles

% ———- % create moving SD for all (i.e. 1 column)

xpadmwSD=nan(nb,1); % create correct size matrix with NaNs

for m=b+1:n+b % loop to incorporate moving-window index m (j index in std2)

startm=m-b; % start row of index m

endm=m+b; % end row of index m

xpadmwSD(m,1)=std2(xpaddetrended(startm:endm,:)); % calculate moving SD

end %m index for rows in loop; j index for columns accounted by std2

% ———- % remove padding (before removing outliers in next section)

startr=1+b; % start row of original data (then add n-1 to get end)

xmean=xpadmean(startr:startr+n-1,:); % mean of x1D

mwSD=xpadmwSD(startr:startr+n-1,:); % movingSD for x1D

mwSD=repmat(mwSD,1,k2); % make mwSD matrix same size as x1D

x2D=x1D; %create x after 2D outliers removed (temporary as copy of x1D)

% ———- % remove 2D outliers outside lower/upper limits

t2=abs(tinv(alpha2/2,k2-1)); %2-tailed t-statistic for alpha2 and k2-1 df

mwCI=t2*mwSD; %scale SD to give moving-window CI

lower=xmean-mwCI; upper=xmean+mwCI;% lower and upper limits

x2D(x2D<=lower | x2D>=upper)=NaN; % find outliers and replace with NaN

nanCt = sum(isnan(x2D),1); % sum NaNs in each column

trialskept2(:,nanCt > 0)=[]; %reduce list of trials to trial numbers kept

x2D(:,nanCt>0)=[]; % removes columns with any NaNs (outliers) to leave x2D

%==========================================================================