Linear segmentation of noisy data

Hi,

I would like to divide my noisy data into a set of linear segments.

As seen in the above figure. I can measure the red noisy signal. The signal in its nature should be composed of linear segments.

Is there a robust way to perform this segmentation? I.e. convert the red noisy data into the black lined segments.

Thanks!

Regards,

Omar

Answer 1

[Star Strider]

This would be easier with your data. Lacking them, I created my own.

Try this:

x = 0:32;
y = [-0.5*(0:11)+10+randi([-1 1], 1, 12), 4.5*(12:14)-50+randi([-1 1], 1, 3), 0*(15:28)+13+randi([-1 1], 1, 14), -4*(29:32)+129+randi([-1 1], 1, 4)];

Lv = ischange(y, 'linear', 'Threshold',5);
Idx = [1 find(Lv), numel(y)];

for k = 1:numel(Idx)-1
    P(:,k) = [x(Idx(k)) 1; x(Idx(k+1)) 1] \ y(Idx(k:k+1)).';
end
Slopes = P(1,:)
Intercepts = P(2,:)

figure
plot(x, y, '-o')
hold on
plot(x(Idx), y(Idx), 'p')
hold off
grid
axis([0  35    0  15])

In this run, that produced this plot:

and these data:

Slopes =
   -0.6667    3.6667         0   -3.0000

Intercepts =
   11.0000  -41.0000   14.0000   98.0000

You will likely need to adjust the 'Threshold' value (here 5) to work with your data. It should produce the appropriate slopes and intercepts for the regression lines. The ischange function was introduced in R2017b. A similar function, findchangepts in the Signal Processing Toolbox (with significantly different arguments and outputs) was introduced in R2016a.

[Omar Alahmad]

Hi Star Strider,

Thanks for going through my question.

What I noticed is that the solution you offered works around half of the time (depending on how the signal looks like).

I ran your code several times. Sometimes I get the correct segmentation, and at other times I do not (see below images).

Here, only the solution in the middle image works well for me.

Is it possible to make the solution more robust, such that computing the correct segments has a much higher probability?

Many thanks in advance!

Best regards,

Omar

[Star Strider]

I have no idea what your signals are or what you are doing with them.

There are several options, one being to use some sort of filtering, such as the Savitsky-Golay filter or a frequency-selective filter, to smoothe the irregularities and enhance the inflection points first. Another is to define an independent variable with greater resolution, then interpolate. That should enhance the inflection points enough to provide reliable detection of the inflection points. There is also the option of exploring the other name-value pair options in ischange, or experimenting with findcchangepts. It all depends on what you want to do with the data.

You need to experiment to see what works best.

[Omar Alahmad]

My signals come in the form as shown in the first figure I posted.

The idea is to detect the change points and construct linear segments in real-time as data is continously being updated.

So I might make it work for one set of data, but perhaps not all.

In any case, I will have to play around with the function parameters and, as you suggested, look at applying some filters first to the data.

Thanks for your help.

Best regards,

Omar

[Star Strider]

As always, my pleasure!

[Image Analyst]

Omar, you don't seem to be getting the hint, so let me make it very clear: attach any data that you are having trouble with in a .mat file with the paper clip icon.

Answer 2

[Bruno Luong]

Similar topic discussed here

Using my solution of BSFK

% Generate random data 
N = 5; % number of linear segments + 1
breaks = cumsum([0, 1+rand(1,N)]);
yb = rand(size(breaks));
coefs = zeros(N,2);
for k=1:N
    coefs(k,:)= polyfit(breaks([k,k+1])-breaks(k),yb([k,k+1]),1);
end
pp = struct('form', 'pp',...
    'breaks', breaks, ...
    'pieces', N, ...
    'coefs', coefs, ...
     'order', 2,...
       'dim', 1);
n = 1000; % number of data point 
sigma = 0.01; % Gaussian noise std
x = linspace(min(breaks),max(breaks),n);
y = ppval(pp,x) + 0.05*randn(size(x));

%%
close all
BSFK(x,y,2,[],[],struct('annimation',1)); % FEX