Bugs in Example_run and nan values in linprog

[fvalerio]

Hi Augustine,

First of all, thank you for writing this cool pipeline. It's going to be invaluably helpful to a lot of people.

I'm trying to run your example script and I get a few consecutive errors. The first one is very easily trackable - In line 36 & 37 of Example_run.m, the code should read
input.data_path1 = ['Example\D',num2str(id)]; % frist dataset, NEED CHANGE input.data_path2 = ['Example\D',num2str(id+1)]; % second dataset, NEED CHANGE since the days folders are inside the example folder.

The second one is in EMD_unit_match.m. Here, line 70, f1 and f2 variables are indexed and I think they shouldn't. So I simply corrected as [C,~] = gdm_nt(f1, f2, mw1, mw2, chan_pos, dim_mask, l2_weight, xStep, zStep, @weighted_gdf_nt);

In example_Run.m, chain_summary and chain_stats functions are run with a missing input, which should be the output folder.

Finally, plot_waveforem.m is called but the plot_waveform function itself is commented out and the file is run as a script. This produces an error. Easy fix.

After these fixes, the code runs and produces an output that looks right on your example data.

Now, I'm trying to run the same code on my own data. I've created a D1-D3 folders in a way analogous to your example code. In it, I have created all the files that should be needed to run the Example_run script. The code runs indeed for a while and it generates the EMD_input folder with input1.mat in it, but eventually it spits out an error which you can find below.

Could you please help with this? Happy to provide more info as needed. Thanks a lot, Valerio

-------ERROR---------

Error using linprog (line 369) LINPROG stopped because some coefficients in the objective or constraint matrices are NaN.

Error in emd_nt (line 79) [x, fval] = linprog(f, A, b, Aeq, beq, lb);

Error in EMD_unit_match (line 62) [x, fval, L2] = emd_nt(f1, f2, w1, w2, mw1, mw2, chan_pos, dim_mask, l2_weight, xStep, zStep, @weighted_gdf_nt);

Error in NT_main (line 16) output = EMD_unit_match(input,output,'pre');

Error in Example_run (line 46) NT_main(input,chan_pos,mwf1,mwf2);

[jenniferColonell]

Hi Valerio, Sorry you ran into these troubles. We found these bugs (and some other issues) and posted a fix a few days ago, but after you cloned the repo, I think. So, as a first step, please make a fresh clone of the repo, and try re-running your data. You'll want to start from the new Example script.

However, I'm not sure if that will fix the issue with your data. This error occurs when the LINPROG module can't converge. This is more common when trying to match very large sets of units. How many units are in your datasets?

[fvalerio]

Hi Jennifer, thank you for your reply

Yeah, Updated the repo and everything runs fine! Thank you.

As per my data: I am trying to match about 600 units across 3 days. But... Is it just a convergence issue? When I run linprog, f contains nan values for me. I tracked the problem down the rabbit hole and ultimately, the issue seems to be with calcCenteredCorrL2. Here centL2 returns nan values because cl2 is a nan itself. This is because, wf1 and wf2 are 0s vectors. Now... I'm not quite sure what's going on here, but what I can see, is that despite mw1 and mw2 look legit (i.e., they posses a waveform detected on some specific channel), the channel selected for comparisons, are "empty channels" (i.e., channels in which the units produce no signal).

Let me know if I can provide further details or if I should just give up. Thanks :):)

[jenniferColonell]

Hi Valerio, If you don't mind giving me your input data, I can take a quick look to see what's going wrong. Naturally the waveform similarity function shouldn't be running the comparison on the zero valued vectors, but I don't know why that's happening with your data (yet). Here's a sharing link to a dropbox folder: https://www.dropbox.com/scl/fo/grsvhx0gd9z6pvdx5b60i/h?rlkey=t2zn6ss3q2jin48wzbrxo33i3&dl=0

[fvalerio]

Hi, thank you for looking into this. You can find the data here

You might notice that 'metrics.csv' is missing, but the code will break before. Normally, I use spike_cluster.npy and spike_times.npy as a substitute for metrics.csv... but that doesn't matter for this purpose. I wanted to clarify just in case

Also I should say.. calcCenteredCorrL2.m has 2 outputs: centCorr and centL2. centCorr always returns a 0 because cv is always 0. Is this supposed to be?

Thanks, Valerio

[jenniferColonell]

Hi Valerio, Sorry it took me a while to get back to you. In answer to your first question: yes, calcCenteredCorrL2.m always returns zero for centCorr because we're not using it -- and the extra calculations take time. From the data, I'm guessing your ksproc_mean_waveforms are actually the templates out of Kilosort. We haven't tested the algorithm using the templates. The actual mean waveforms, calculated from either the filtered data or the drift corrected data are preferred, because the amplitudes are more accurate. Using templates presents two problems: (1) that we only include site with signal > 60 uV for the position fitting. You can make the templates look like that by scaling them (300 is a reasonable value). (2) All sites outside the template radius are set identically to zero. When the peak of the template isn't centered in the set of template sites, this can lead to running the L2 comparison on two sites with exactly zero signal, which the original code had not allowed for. I've made that fix in my fork of the repo (https://github.com/jenniferColonell/Neuron_Tracking). We'll get this repo caught up soon. I encourage you to use real mean waveforms (e.g. using the C_Waves too, but there are others). However, if you want to run with the templates, scale them first by multiplying by 300. The only metric that is getting used in the posted code is the firing rate, which you can calculate from spike_times.npy and spike_clusters.npy. Let us know if you have any other questions!