VBN trials change time disagrees with stim table change time

[corbennett]

Describe the bug It appears that the change time reported in the trials df in the VBN nwb files precedes the stim table change time by ~8 ms. I suspect this is because the 8 ms offset described here wasn't applied.

To Reproduce Load an nwb file into 'session'. Then:

st = session.stimulus_presentations
trials = session.trials

st_change_times= st[(st.is_change==True)&(st.active)]['start_time']
trials_change_times = trials[trials.stimulus_change==True].change_time
print(np.mean(st_changes.values - trials_change_times))

which gives ' 0.0088...'

Expected behavior These values should be identical.

[danielsf]

@corbennett I ran your test on all of the NWB files. The vast majority of them show disagreement at the 8 ms level as you describe. However, there are some that show disagreement at the 0.02 second level (which, eerily, is the value we used for monitor_delay on the NP.0, NP.1 rigs)

This table shows the NWB path along with the mean +/- stddev of the difference between stimulus change times and trials change times as you calculated above.

ecephys_session_1084427055.nwb -- 2.50e-02 +/ 1.59e-04
ecephys_session_1086200042.nwb -- 2.50e-02 +/ 1.65e-04
ecephys_session_1093642839.nwb -- 2.49e-02 +/ 1.86e-04
ecephys_session_1087723305.nwb -- 2.50e-02 +/ 1.49e-04
ecephys_session_1052533639.nwb -- 1.88e-02 +/ 7.98e-03
ecephys_session_1086410738.nwb -- 2.50e-02 +/ 1.58e-04
ecephys_session_1044597824.nwb -- 2.42e-02 +/ 1.82e-04
ecephys_session_1044389060.nwb -- 2.42e-02 +/ 1.66e-04
ecephys_session_1087993643.nwb -- 2.50e-02 +/ 1.58e-04
ecephys_session_1093864136.nwb -- 2.49e-02 +/ 1.68e-04

As we discussed yesterday, while the stimulus table does its own careful calculation of monitor_delay, the trials table just applies a pre-computed scalar value of 20 ms as monitor delay. If not for the above cases, I was going to suggest we just adjust that value by 8ms. This looks like it is going to be more tricky to solve.

[danielsf]

One possible "all I have is this hammer" solution is to drop the change_time column from the trials table. On some level, we are suffering from the fact that we are storing the same data (the time at which the stimulus changed) in two places: the stimulus_presentations table and the trials table. Given that the trials table has a change_frame column and the stimulus presentations column has both a start_frame and a start_time column, we can use the former to lookup the timestamp from the latter using something like this.

>>> change_frames = trials[trials.stimulus_change].change_frame
>>> stim_frames = stim[stim.is_change & stim.active].start_frame
>>> stim_times = stim[stim.is_change & stim.active].start_time
>>> frame_to_time = {f:t for f, t in zip(stim_frames, stim_times)}
>>> change_times = np.array([frame_to_time[f] for f in change_frames])
>>> np.array_equal(change_times, stim[stim.is_change & stim.active].start_time.values)
True
>>> 

(There are probably easier ways, but I am not as comfortable with pandas as most people; I prefer carrying around raw numpy arrays of data).

If we trust that the the monitor delay was handled correctly in the stimulus table and are dubious of what is being done in the trials table, I think this is the easiest solution to bring the two datastreams into alignment. It will require our users to do some extra work to get change times from the trials table, but it minimizes the risk that we are contradicting ourselves in our data.

The alternative solution is to port over the monitor delay calculation code from ecephys_etl_pipelines into the SDK and connect it to the part of the SDK that creates the trials table. That is probably a few days of work, since the trials table is still in a very difficult to refactor state.

[danielsf]

I guess an even more brute-force solution would be to calculate the trials and stimulus_presentations as they are, and then add a postprocessing step that uses the lookup table I propose above to patch the change_time column in the trials table. Then, no code needs to be ported over into the SDK and we are again guaranteed that the two tables agree.

That is a very ugly solution and my team's eyebrows will go through the roof, but that's my problem, I guess.

[danielsf]

Maybe this won't work. I see that it is possible for trials.change_frame to be non-NULL even when trials.stimulus_change == False (I guess these are "shame changes"?). In that case, there is no frame-to-time mapping that we can look up from the stimulus_presentations table.

[corbennett]

I chatted with @danielsf about this yesterday, but just to document:

@corbennett I ran your test on all of the NWB files. The vast majority of them show disagreement at the 8 ms level as you describe. However, there are some that show disagreement at the 0.02 second level (which, eerily, is the value we used for monitor_delay on the NP.0, NP.1 rigs)

This table shows the NWB path along with the mean +/- stddev of the difference between stimulus change times and trials change times as you calculated above.

ecephys_session_1084427055.nwb -- 2.50e-02 +/ 1.59e-04
ecephys_session_1086200042.nwb -- 2.50e-02 +/ 1.65e-04
ecephys_session_1093642839.nwb -- 2.49e-02 +/ 1.86e-04
ecephys_session_1087723305.nwb -- 2.50e-02 +/ 1.49e-04
ecephys_session_1052533639.nwb -- 1.88e-02 +/ 7.98e-03
ecephys_session_1086410738.nwb -- 2.50e-02 +/ 1.58e-04
ecephys_session_1044597824.nwb -- 2.42e-02 +/ 1.82e-04
ecephys_session_1044389060.nwb -- 2.42e-02 +/ 1.66e-04
ecephys_session_1087993643.nwb -- 2.50e-02 +/ 1.58e-04
ecephys_session_1093864136.nwb -- 2.49e-02 +/ 1.68e-04

As we discussed yesterday, while the stimulus table does its own careful calculation of monitor_delay, the trials table just applies a pre-computed scalar value of 20 ms as monitor delay. If not for the above cases, I was going to suggest we just adjust that value by 8ms. This looks like it is going to be more tricky to solve.

I've confirmed that all of these sessions are ones for which the monitor lag was measured to be 37 ms, which makes these numbers make sense. Session 1052533639 had a variable monitor lag that bounced between 20 and 37, explaining why it has about half the error (and is more variable).

Maybe this won't work. I see that it is possible for trials.change_frame to be non-NULL even when trials.stimulus_change == False (I guess these are "shame changes"?). In that case, there is no frame-to-time mapping that we can look up from the stimulus_presentations table.

Yes. That's correct. The trials table gives 'sham_changes' a change frame.

I think the correct strategy is close to what you suggested above. 1) We report an 'uncorrected' change time in the trials table. This change time is when the task computer thought the change occurred and doesn't take display lag into account. We will call this 'change_time_no_display_delay' 2) To get actual change times, the user will use the stimulus presentations table. So something like the following would add the display corrected change times to the trials table:

from functools import partial
def get_change_time_from_stim_table(row, table):

    change_frame = row['change_frame']
    if np.isnan(change_frame):
        return np.nan

    change_time = table[table.start_frame==change_frame]\
                    ['start_time'].values[0]

    return change_time

get_change_times = partial(get_change_time_from_stim_table, table=stimulus_presentations)
change_times = trials.apply(get_change_times, axis=1)
trials['change_time_with_display_delay'] = change_times

[danielsf]

addressed by #2469