erfanzabeh / WaveMonk

Detection and Quantification of Traveling Wave in Monkey Utah array recordings
5 stars 0 forks source link

Bugs in Periodic2ContinuesPhase2.m #3

Open OccassionalBioinformatician opened 1 year ago

OccassionalBioinformatician commented 1 year ago
[idx,c] = kmeans(clustervector,ceil(range(clustervector)));

It is very plausible that the number of rows in clustervector — which is equal to the number of electrodes (say, x) — exceeds the ceiling-ed range!

Let's reformulate the code for a single electrode:

% phase_interval is a t*1 array of phase angles for one electrode for one epoch.

timesample = 1:size(phase_interval,1)
K = zeros(size(phase_interval));

[~, loc] = findpeaks(phase_interval);

for k = 1:numel(loc)
    K = K + (timesample - loc(k) > 0)';
end
if phase_interval_2(1)>5
    K = K-1;
end

Phase_Out = phase_interval+ K*2*pi;
Phase_Out_Proportion = Phase_Out./(2*pi);

The number of peaks detected for an electrode can obviously exceed x (an independent variable) if the time span is sufficiently large! Accordingly, all time points roughly exceeding (roughly since the original phase angle can be negative and slightly offset the gain) the point when the xth peak was observed will return values greater than x in Phase_Out_Proportion. So, let's say that my (2,4) electrode returns a value exceeding x for the first time for the 3461st second.

Now, when Periodic2ContinuesPhase2 runs the loop for trailtime = 3461 (and onward), ceil(range(clustervector)) returns a value exceeding x. The clustering fails. So, it is perhaps wise to (1) detect peaks setting NPeaks = x, (2) calculate the minimum time t_min taken by an electrode (of all the electrodes) to hit the limit across epochs, and (3) suggest trimming down epoch-size to it for a meaningful analysis.