Any suggestion for the method I should use for calculating the sum of the event which has zero values and its duration with matlab?

Hi! I need a suggestion about which method I should use in doing this task. Before, I use Excel VBA to do this and it worked well. Then, I have to change my method to Matlab and I got confused as I am really new in this programming software. I have a data in the form of table, consist of datetime and flow (numeric). The flow is dynamic (like the plot in peak analysis), as it has both zero and non zero values in it. I need to calculate the total events of the condition with zero values (dry period) and the duration of each events. Do you have any suggestion for me to do it in matlab? Here is what I do in Excel VBA (picture attached)

1. I add one column contain sequence number from the first row until the last row beside the datetime data and flow data.
2. I filter the data based on the flow. If it is not zero, I delete the flow and other data in its row.
3. Now, the sequence number may not be in order form because some have been deleted. From this sequence number, I then make a column more beside sequence number column for unique number. The first row would be'1', and if it is still in sequenced form, it remains 1, but if it jumps to three number above or more, it should be '2' and this pattern continue until the last row.
4. After I got all the unique number which represents the event I need, I use pivot table to find maximum and minimum datetime of each event, then I substract it. That way I can find the duration of each event.

Do you have any suggestion for me to do it in matlab (with matlab code?) or should I also do the same as in ExcelVBA? In matlab I also have the table form consist of the data. I really appreciate for your help. Thank you so much in advance.

Answer 1

[Guillaume]

Looks like you use excel macros, not VBA. In any case, no you should not do the same method at all in matlab. Matlab is a lot more powerful but requires a different way of thinking.

Giving us a file we can play with would have been more useful than screenshots. Anyway, first thing you need to do is import your data in matlab. For that use readtable with appropriate options (possibly none). Something like:

flowtable = readtable('c:\somewhere\somefile.xlsx');

If necessary rename the columns of the table to something useful. In the following I assume the flow column is called Flow. You can then detect dry periods, their beginnings and ends:

dryperiod = flowtable.Flow == 0;
transitions = diff([false, dryperiod, false]); %is 1 at the start of a dry period, -1 after the end
drystarts = find(transitions == 1);
dryends = find(transitions ==  -1) - 1;

You can then build your result table. I'm assuming that your flowtable has a single column named DateTime of type datetime

drytable = table(flowtable.DateTime(drystarts), ...
                 flowtable.DateTime(dryends), ...
                 flowtable.DateTime(dryends) - flowtable.DateTime(drystarts), ...
                 'VariableNames', {'Start', 'End', 'Duration'})

[Kasih Ditaningtyas Sari Pratiwi]

Hi, Guillaume. Thank you for your reply. I try to use your code above :

%%Dry period analysis
flowtable = finalCSVnew(:,[1,7:8]);

dryperiod = flowtable{:,2} == 0;
transitions = diff([false, dryperiod, false]);                              %is 1 at the start of a dry period, -1 after the end
drystarts = find(transitions == 1);
dryends = find(transitions ==  -1) - 1;

drytable = table(flowtable.DateAndTime(drystarts), flowtable.DateAndTime(dryends), ...
                 flowtable.DateAndTime(dryends) - flowtable.DateAndTime(drystarts), ...
                 'VariableNames', {'Start', 'End', 'Duration'});

However I got an error : "Error using horzcat Dimensions of matrices being concatenated are not consistent". the source of the table for flowtable is also in the table form. Here we do not use matrices, or do you think it is because of my data source which is in the form of table? Or Does it has something to do with my datetime format? I use "dd.mm.uuuu HH:mm:ss" here

[Guillaume]

Sorry my mistake, it should have been

transitions = diff([false; dryperiod; false]);

i.e. a vertical concatenation

[Kasih Ditaningtyas Sari Pratiwi]

Thank You Guillame. This code works perfectly. I have to study more because at some part, I still did not understand. For example, I want to change the threshold, for dry period it should be < 2. What should I change in this code?

dryperiod = flowtable{:,2} == 0;
transitions = diff([false; dryperiod; false]);                              %is 1 at the start of a dry period, -1 after the end
drystarts = find(transitions == 1);
dryends = find(transitions ==  -1) - 1;

Should it be like this?

dryperiod = flowtable{:,2} < 2;
transitions = diff([false; dryperiod; false]);                              %is 1 at the start of a dry period, -1 after the end
drystarts = find(transitions == 1);
dryends = find(transitions ==  -1) - 1;

And I still did not get why dryends should be "find(transitions == -1) - 1;" and "transitions=diff([false; dryperiod; false]);", the comma is for horizontal matrices and the semicolon is for vertical matrices? could you please explain what is the function of those two command?

Thank you so much Guillame for your help. I was thinking about doing the same as in my excel, but this is shorter and really fast. Thank you again :)

[Guillaume]

You're basically creating a logical vector of true (1) and false (0) values that tells you when the dry periods are. You can create this any way you want. If indeed it's when the flow rate is below 2, then flowtable{:,2} < 2 will do that.

The beginning and end of the dry period is when that logical vectors transitions from false to true and true to false respectively. So you need to detect these transitions. The only times such a transition wouldn't happen is when you start or end with a dry period since you start or end directly with a true. To ensure we still get a transition in that case, we simply surround the logical vector by false.

There are several ways you can detect the transitions:

• If you take the diff as I have done, transitions from false (0) to true (1) will be 1 and transitions from true to false will be -1. You can search either with find. Because of the added false at the beginning, find returns the index of second element in the sequence [false; true} (for starts) and [true; false] (for ends). The actual end is the first element hence the subtraction in that case
• Because the vector is surrounded by false, you're guaranteed that the non-zeros values in the diff will always be in pairs ordered as 1, -1. So, you if you get the locations of the non-zero values, odd indices of this vector are the starts, end even indices are the ends:

transitions = find(diff([false; dryperiod; false]));
drystarts = transition(1:2:end); 
dryends = transitions(2:2:end);

• You could also look directly for these [false; true] and [true; false] transitions with strfind (which despite its name also find patterns in numeric vectors):

drystarts = strfind([false; dryperiod]);
dryends = strfind([dryperiod; false]);

[Kasih Ditaningtyas Sari Pratiwi]

Thank you so much for the clear explanation, Guillaume. I now understand about the code that you wrote.

I am sorry I have one more question, just to check whether I got it right or not. If for example I try to do vice versa. This time, I want to find the value above zero (I want to eliminate zero), the duration for each event as in the dry period, and the max flow for each event. Is this right? :

    flowtable = finalCSVnew(:,[1,7:8]);                                         

    peakflowperiod = flowtable{:,2} > 0 ;                                           
    transitions = diff([false; peakflowperiod; false]);                              

      pfstarts = find(transitions == 1);
      pfends = find(transitions == -1) -1;
      location= flowtable.Flow(pfstarts):flowtable.Flow(pfends);
      maxflow = max(location);

      peakflowtable = table(flowtable.DateAndTime(pfstarts), flowtable.DateAndTime(pfends), ...
                       flowtable.DateAndTime(pfends) - flowtable.DateAndTime(pfstarts), ...
                       flowtable.Volume(pfends)-flowtable.Volume(pfstarts), ...
                       maxflow, ...
                     'VariableNames', {'Start', 'End', 'Duration','ChangeInVolume','MaxFlow'});

    numMaxflow=length(peakflowtable.Duration);                                     % calculate the number of max flow

Is it correct for finding the max value in that range? I actually also use findpeaks but the drawbacks of using findpeaks is that they will have many max flow in one event and thus the duration for each event is not also correct. So I consider using function max. Thank you again for your help Guillaume. I am sorry I am asking too much >.<

[Guillaume]

The bit where you find the max flow is not going to work. You can't pass vectors to the colon operator (well, you can but it's only going to use the first element). One way to do what you want:

maxflow = arrayfun(@(s, e) max(flowtable.Flow(s:e)), pfstarts, pfends);

I'd change the last line to:

numMaxflow = height(peakflowtable);

I never use length. I prefer numel for vectors or size for matrices. But for tables you have the height and width methods.

[Kasih Ditaningtyas Sari Pratiwi]

Awesome. It works perfectly. I also check the result with my excel and it is also correct, but with shorter code and fast. Many many thanks Guillaume. I just check arrayfun function in the mathwork description, here we have @(s:e), what does it mean? Because in the example, it only has one data. Thank you again and I am sorry for asking many things >.<

[Guillaume]

@(s, e) max(flowtable.Flow(s:e)) is an anonymous function with two input arguments, s and e. As with normal function the names assigned to these inputs don't matter. You just have to use the same names in the body of the function. The anonymous function is more or less equivalent to the normal function:

function result = noname(s, e)
   result = max(flowtable.Flow(s:e));
end

except that the normal function, defined in its own m file, wouldn't be able to access flowtable.Flow (a nested function could).

I'm passing two arrays to arrayfun, pfstarts and pfends, hence the anonymous function needs two inputs. s and e receives each element of the respective array.

Overall, the arrayfun is exactly equivalent to the loop:

maxflow = zeros(size(pfstarts));
for iter = 1:numel(pfstarts)
    s = pfstarts(iter);
    e = pfends(iter);
   maxflow(iter) = max(flowtable.flow(s:e));
end

[Kasih Ditaningtyas Sari Pratiwi]

Awesome. I understand now why we use the function. You're a lifesaver! Many thanks Guillaume, for helping me in solving this matter and you always give me clear explanation about the matlab code you wrote.

[Kasih Ditaningtyas Sari Pratiwi]

Hi Guillaume. I'm sorry again to ask you. I forgot to find the location of the max flow data. I try to write the codes but still they did not work. Here are the codes:

peakflow = arrayfun(@(s, e) max(flowtable.Flow(s:e)), eventStarts, eventEnds);
peakflowlocs=flowtable.DateAndTime(peakflow);

After running the code, I got an error message: "Error using tabular/dotParenReference (line 108) Subscript indices must either be real positive integers or logicals." I wonder why I can not use it.

I try some example from the web by using find :

s=peakflow;
y=flowtable.Durchflusslm;
index=find(y==s);
X_point=flowtable.DateAndTime(index);

then now, it says "Matrix dimensions must agree." And If I use this below code, it can run, but the result are so many. It does not refer to the location of max flow itself. Probably it has unique value as well, so others with the same value are also considered:

c= ismember(flowtable.Durchflusslm,peakflow);
indexes=find(c);

Could you please help me to solve it? I am sorry for asking many things.

[Guillaume]

Well, the location of the maximum is the second output of max. This will give you the index relative to the start of the dry period. If you want the absolute index, add the relative indices to the event starts:

[peakflow, peakflowlocrel] = arrayfun(@(s, e) max(flowtable.Flow(s:e)), eventStarts, eventEnds);
peakflowlocabs = peakflowlocrel + eventStarts - 1;

[Kasih Ditaningtyas Sari Pratiwi]

Thank you so much Guillaume. You saved my life again! :) I use peakflowlocabs to then draw the max flow in original data and it works. But actually I do not understand why we have to add "eventStarts - 1" there in looking for "peakflowlocabs". Could you please explain it? Thank you again for your help