stimulus trigger times from CED files

rebecca-mease commented 3 years ago

Hi all, Rebecca here. I am trying to make some of our home-grown analysis interact with NWB files in the best way possible. Would it be possible to include TTL pulse info (stimulus information) from the .smrx files? Current there is nothing there as far as I can tell.

CodyCBakerPhD commented 3 years ago

@rivkahrivkah I'll have to re-download the example data you originally shared to fully debug (eta of a few days), but could you check and show what is written in the 'processing' tab? The mechanical and laser stimulus information is being converted by the CEDStimulusInterface and stored as a table of TimeIntervals (for start/stop times of the TTL pulses); this should be writing to a module called 'Stimulus'.

rebecca-mease commented 3 years ago

Hi, thanks for the rapid reply--here's what I get using util.nwbTree. Is this what you mean? apologies, it just the last few days that I have finally had time to really try to interact with the NWB files myself!

rebecca-mease commented 3 years ago

I'll doublecheck w/ Ross that you have access to a file that includes 1) mechanical stim TTL pulse, 2) raw recorded mechanical stimulus (actual stimulus is not 1:1 with command TTL pulse), and 3) optogenetic stimuli.

CodyCBakerPhD commented 3 years ago

For exploring the contents of an NWBFile I typically use the following lines

from pynwb import NWBHDF5IO

io= NWBHDF5IO(path="name_of_file.nwb", mode="r")
nwbfile = io.read()

# Do the things to interact with nwbfile, such as
print(nwbfile.processing)
...

# Always remember to close when done to prevent strange errors elsewhere
io.close()

and starting from a simple print(nwbfile), you can see all the attributes that you can access such as nwbfile.acquisition, etc. Direct data access sometimes requires further indexing commands such as nwbfile.acquisition["name_of_electrical_series"].data[slice_in_time, slice_in_channels].

I'll doublecheck w/ Ross that you have access to a file that includes 1) mechanical stim TTL pulse, 2) raw recorded mechanical stimulus (actual stimulus is not 1:1 with command TTL pulse), and 3) optogenetic stimuli.

I do recall some data we prototyped on ("pt1 15 + mech.smrx" of M365 or "M113_C4.smrx") had channels for mechanical and laser TTL pulses, which we decoded and stored. It also had a TimeSeries of mechanical stimulus pressure.

Is this NWBFile the result of the SpikeInterface pipeline only? If so, it did not include the stimulus or raw acquisition at that step, but only the LFP and spike sorted units. The convert_ced.py script is what writes those two outside the SpikeInterface pipeline; note that to append only the stimulus data to your NWBFile, set overwrite=False and comment out the CEDRecording lines from the source_data and conversion_options.

Let me know if you have any trouble running those steps; the API may have changed since the last time we debugged it.

rebecca-mease commented 3 years ago

Thanks Cody, I think this answers my question for now; we will try it out.

ross-folkard commented 3 years ago

Hi all, Rebecca here. I am trying to make some of our home-grown analysis interact with NWB files in the best way possible. Would it be possible to include TTL pulse info (stimulus information) from the .smrx files? Current there is nothing there as far as I can tell.

Hey @CodyCBakerPhD,

Thanks for making this quick script to get the TTL trigger times. It looks like the data I'm giving the script is suitable and should run (everything was valid apparently), but I'm getting an invalid argument regarding timestamps and timezone. Would you mind taking a quick look to see if I've done something wrong here?

Many thanks!

Ross

CodyCBakerPhD commented 3 years ago

@rf13734 Can you show me how you're setting the session_start_time there? That ought to be the only datetime object that could cause that unfortunate upstream error.

ross-folkard commented 3 years ago

Hi @CodyCBakerPhD ,

I add my session start time here:

CodyCBakerPhD commented 3 years ago

@rf13734 To clarify, these screenshots indicate two different pipelines; one for adding the sorting_phy and recording_lfp objects, and the other for adding the entire raw recording (can take a while) with behavioral stuff (including mechnical/laser TTL, etc).

When using both of these pipelines in tandem to write to a single file, one pipeline or the other should have overwrite=False so that it appends to the file. My best guess to what is happening above is when you try running the Recording/Stimulus pipeline, without a session_start_time specified and overwrite=True, it tries to create a new NWBFile, but doesn't have the critical required information for the session_start_time. It's unfortunate this terminates in a classic messy hdmf error at the final io.write(), I'll have to add some checks for that at some level...

Two ways to fix this should be (a) set overwrite=False in the recording/stimulus pipeline that was originally causing the issue, no start time should be needed in that case, or (b) set the session start time at the same level of metadata as the session_description in that same pipeline, and write that first. Remember, if you don't want to wait on writing the entire recording data, you can comment or otherwise remove that field from the source_data and it should just write the stimulus.

Also, if you are making the NWBFile first from the spikinterface pipeline, I would recommend copying/pasting that file to keep a backup in case you get more failures during io.write() in the other pipeline, which has the potential of irreparably corrupting the file you're intending to append.

ross-folkard commented 3 years ago

@CodyCBakerPhD Thanks very much for the detailed breakdown! I think I might have accidentally overwritten said file, so will try with a different (backed up) file and will let you know.

ross-folkard commented 3 years ago

Hey @CodyCBakerPhD , thanks for the advice, I'm no longer getting the same error now that I manually added the session_start_time.

I am however getting another error message suggesting that the nwb file may not have saved correctly in the first place perhaps.

There is an nwb file in here:

...which was generated by a separate pipeline.

Any ideas on what I've missed? No rush, I'm going on holiday for a week so probably won't progress this much till I get back.

Once again, many thanks!!

CodyCBakerPhD commented 3 years ago

Sadly, that exact error is what you see when a .nwb file has been created, but is somehow corrupted internally and so cannot be read.

What I would maybe suggest doing first, within each pipeline, is trying to run a stub_test for all data types, and writing both to separate individual files. If the process completes without error, then it should be safe to attempt the combination to a single file as you were just trying to do.

ross-folkard commented 3 years ago

Hi @CodyCBakerPhD , Just to update, appending the file with the TTL data deems to have worked i.e. the file successfully saved.

However I'm struggling to find any sign of them in the stimulus view of the widgets, though I have both laser and mechanical TTL data from Spike 2.

I also could not seem to find the MechStim data from the pressure sensor.

Any help on this would be amazing, looking forward to hearing from you.

Ross

CodyCBakerPhD commented 3 years ago

Hi @CodyCBakerPhD , Just to update, appending the file with the TTL data deems to have worked i.e. the file successfully saved.

Great!

However I'm struggling to find any sign of them in the stimulus view of the widgets, though I have both laser and mechanical TTL data from Spike 2.

Not sure I understand - from your screenshots the IntervalSeries seem to be showing up under LasterStimulus and MechanicalStimulus. Did you want an actual plot of those vs. time to see side-by-side with the other TimeSeries?

I also could not seem to find the MechStim data from the pressure sensor.

From the interface, this looks like it should be under NWBFile -> Stimulus -> MechanicalPressure and should be a TimeSeries

ross-folkard commented 3 years ago

Hi @CodyCBakerPhD, thanks for the quick response!

From the interface, this looks like it should be under NWBFile -> Stimulus -> MechanicalPressure and should be a TimeSeries

So should the MechStim time series be here?

CodyCBakerPhD commented 3 years ago

Hmm, so the fact that the tab is present means it must be, but from the sliders it appears something got transposed in terms of channels vs. time... Assuming this file you're working with is small could you send it my way so I can inspect the inner elements?

ross-folkard commented 3 years ago

@CodyCBakerPhD Sure, have you a preferred email address?

CodyCBakerPhD commented 3 years ago

Let's try cody.baker@catalystneuro.com to start.

If the file is too large and requires a google account for a drive upload, cody.c.baker.phd@gmail.com otherwise.

ross-folkard commented 3 years ago

@CodyCBakerPhD Okay, sent it over, thank you very much!

ross-folkard commented 3 years ago

@CodyCBakerPhD Interestingly, the file that supposedly has the TTL pulses in won't open on Matlab, but the one with just the spikes in will:

ross-folkard commented 3 years ago

@CodyCBakerPhD , apologies, I think I generated a file minus the spikes this time. However, when I generate a file with what appear to be both (with the siZe of the TTL file + the Spikes file (1.6gb and 37MB respectively) I still get the same error trying to read the file in Matlab:

I don't get this error when reading the file that only contains spike time information.

CodyCBakerPhD commented 3 years ago

@rf13734 Thanks for the data share - I have confirmed that the mechanical stimulus was indeed transposed by mistake, which is why the widget didn't work on it. I've fixed that on master if you update the code, it should now correctly index time as the first axis.

As for the file not loading through the Matlab reader, I'll have to refer you to MatNWB as I don't have much experience with that. I can confirm on my end that the file you shared loaded without issue in PyNWB and HDFView, so the file isn't corrupted or anything like that.

ross-folkard commented 3 years ago

@CodyCBakerPhD

Beautiful!

Thanks so much!

bendichter commented 3 years ago

@rf13734 Could you try reinstalling MatNWB? There have been a few recent changes that might fix this

ross-folkard commented 3 years ago

@rf13734 Could you try reinstalling MatNWB? There have been a few recent changes that might fix this

Hey @bendichter. I just tried, and go the same error.

rebecca-mease commented 2 years ago

Hi, I have a question about the optogenetic stimulus--is there a particular reason that it is not incorporated as an OptogeneticSeries? Ultimately we would like to include laser power etc. and all of the OptogeneticStimulusSite info. Thanks for any input, we are trying to coordinate generation of the NWB files with the SI pipeline and proper entry of the stimulation metadata.

CodyCBakerPhD commented 2 years ago

I don't recall seeing that specific type of info in the CED file we prototyped on, so without the extra metadata context it would belong in acquisition as a more basic type.

If you are willing to add that metadata at this point, though, then yes it would be quite possible for you to extend it to the OptogeneticSeries, which can be easily written from these commands in the ogen modules: https://nwb-schema.readthedocs.io/en/latest/format.html?highlight=optogenetic#optogeneticseries

bendichter commented 2 years ago

rebecca-mease commented 2 years ago

Thank you both for the information. I think we are getting close to having the data organized in a way which will let us do the analysis we would like. Our (hopefully!) final bottleneck is making sure that our various stimulus conditions are stored in the correct way in NWB files.

The goal is to be able to separate different types of trials for viewing in NWBWidgets and programmatic access to triggered spike time and LFP data for later analysis beyond PSTHs, rasters, etc.
Different trials are inferred from combinations of optogenetic and sensory stimuli, which we currently have implemented in Matlab for non-NWB data--example "Conditions" structure here: https://heibox.uni-heidelberg.de/f/0f7de498c64b4b2784d2/?dl=1
(this example just has different frequencies and powers of laser stimulation, others might have additional trials with mechanical stimulation, etc.)
What is the best place to store this sort of trial information? I'd like to make sure that this trial-based data goes into the right NWB container so that we can best leverage NWB-compatible tools, but it’s not clear to me where it should go.
Could you point us in the direction of any examples on Github with different types of trials? Thank you in advance for any advice on how to best organize our data!

CodyCBakerPhD commented 2 years ago

@rivkahrivkah Adding trial structure is a great way to encapsulate those details, and pynwb has easy to use tools for adding this information to the file! Asssuming you have an NWBFile object ("nwbfile") and the frequencies and powers loaded as some kind of iterables ("all_trial_frequencies/all_trial_powers")

# Add custom columns, specifying scientific units in description
# See https://pynwb.readthedocs.io/en/stable/pynwb.file.html#pynwb.file.NWBFile.add_trial_column
nwbfile.add_trial_column(name=frequency, description="The frequency of the laser in Hz.") 
nwbfile.add_trial_column(name=power, description="The average power of the laser in watts.")

for start_time, stop_time, frequency_per_trial, power_per_trial in zip(start_times, stop_times, all_trial_frequencies, all_trial_powers):
    # Add each trial iteratively (easiest way, imo)
    # See https://pynwb.readthedocs.io/en/stable/pynwb.file.html#pynwb.file.NWBFile.add_trial
    nwbfile.add_trial(start_time=start_time, stop_time=stop_time, frequency=frequency_per_trial, power=power_per_trial)

Though, a couple notes here. (a) For best performance in the NWBWidgets, follow these Best Practices (DynamicTable section) for following naming conventions so our algorithms can auto-detect things like alignment times, (b) it is critical that all references to timing throughout the NWBFile (here, start_time and stop_time of each trial) are made with respect to the session_start_time of the file, which is typically also the start time of the recording.

Most of our examples on GitHub are a tad more programatic in their extraction from .csv and .mat files than the demo code above, but here is a decently simple one illustrating the process in action in one of our behavior interfaces from another project (which I'd recommend making and integrating with your existing NWBConverter/DataInterface setup, let me know if you have questions about that): https://github.com/catalystneuro/buzsaki-lab-to-nwb/blob/add_watson/buzsaki_lab_to_nwb/watson_code/watsonbehaviordatainterface.py#L80-L95

catalystneuro / mease-lab-to-nwb

stimulus trigger times from CED files #46