How to mark step changes of various sizes

This is my code for part of a project im doing:

time_run2 = transpose((1:numel(second_run_data(:, 1))') / 60);

figure(3);

plot(time_run2, second_run_data(:, 1))

whish is used to plot a section of data from a larger set. This data follows a step pattern however each step is brocken up into a bunch of little steps and then a section of flat followed by another step and more flat. some of the steps are slightly curved so i cant simply use the ischange() function to find them becasue it wil find multiple changes in each step. My task is too mark each step with a circle or an x whatever but im only allowed to have 1 mark per step. Is there any way of doing this becasue I have no idea

Answer 1

[William Rose]

@Justin,

You may want to differentiate your signal with diff(), and then use findpeaks(). Differentiation amplifies high frequency noise, so it can produce spurious peaks. Therefore some smoothing may be required as well. You can use Matlab's smooth(), or some other lowpass filter. If the smoothing is by a linear filter (and smooth(), with the default method, is linear), then it does not matter whether you smooth then differentiate, or vice versa. I have used this approach to identify changes in heart rate in noisy data.

You say "some of the steps are curved" which makes me wonder if there is a bit of judgement involved in the definition of a step. It is OK if there is - that can be incorporated in code - but you have to give some examples for guidance. It would help if you would post the data for which you want to identify changes - maybe with a hand drawing of where you think the steps are.

Good luck.

[Justin]

I have attatched the data file if that helps, was hard to describe sorry about that.

[William Rose]

@Justin,

d=importdata('FT02_run_2 (1).csv');
t=d.data(:,1);                  % time [min]
Q=d.data(:,2);                  % flow [mL/min]
fs=(length(t)-1)/(t(end)-t(1)); % sampling rate [1/min]
fprintf('Sampling rate=%.1f samples/min\n',fs);

Sampling rate=60.0 samples/min

% plot the data
figure
subplot(311), plot(t,Q,'-r.')
grid on; ylabel('Flow (mL/min')
subplot(312), plot(t,Q,'-r.'); 
ylim([9,14])  % same plot, expanded y-scale
grid on; xlabel('Time (min)'); ylabel('Flow (mL/min')
subplot(313), plot(t,Q,'-r.'); 
xlim([34,40]); ylim([12,13.5])  % same plot, expanded y-scale
grid on; xlabel('Time (min)'); ylabel('Flow (mL/min')

Do you want to identify the times when there is a change? Let's ignore the first and last few tenths of a minute, when the system starts at, and ends at, zero. There are a variety of ways to attempt it. Let's try some and see what happens.

1. Find "change times", Tc, where mean(Q) for a time Tw before Tc is different from mean(Q) for a time Tw after Tc, by at least some threshold flow difference, Qt. This is just a start.

[Justin]

perfect thank you so much really appreciate it

[William Rose]

@Justin, you are welcome.

d=importdata('FT02_run_2 (1).csv');
t=d.data(:,1);                  % time [min]
q=d.data(:,2);                  % flow [mL/min]
N=length(t);                    % number of samples
fs=(N-1)/(t(end)-t(1)); % sampling rate [1/min]

Find transition times

% Define constants
Tw=0.9;             % segment width [min]
Ts=0.03;            % separation between before and after segments [min]
nw=round(Tw*fs);    % segment width [samples]
ns=round(Ts*fs);    % segment separation [samples]
deltaQmin=0.1;     % minimum flow change [L/min]
sdmax=0.05;         % max st dev of flow within a segment [L/min]

Look for transitions

j=0;  % number of transitions found so far
for i=ns+nw:N-nw
    s1=std(q(i-ns-nw+1:i-ns));    % seg 1 SD
    s2=std(q(i+1:i+nw));          % seg 2 SD
    m1=mean(q(i-ns-nw+1:i-ns));    % seg 1 mean
    m2=mean(q(i+1:i+nw));          % seg 2 mean
    if s1<sdmax && s2<sdmax && abs(m1-m2)>deltaQmin
        j=j+1;
        tc(j)=i/fs;
    end
end
ntrans=length(tc);  % number of transitions found

Display and plot results

fprintf('Found %d transitions.\n',ntrans)

Found 26 transitions.

% plot the data
figure
plot(t,q,'-r.')
grid on; xlabel('Time (min)'); ylabel('Flow (mL/min'); ylim([9,14])
hold on
xline(tc,'b')

It is far from perfect, but it is an example of something you could build on.

Good luck.