We’ve been using the code provided by @lmazzuca on Steinmetz dataset (we binned it to 0.005s) and we’re unable to find state sequences consistent across trials. We’ve run multiple simulations with several parameter configurations, on different brain areas, different number of trials and on different trial types. The only consistent dynamics we are able to recover seem to be based on longer time scales - not capturing trial related dynamics.

Thanks for sharing your results! Did you fit the same HMM to all simultaneously recorded neurons? In that case, I have a suggestion for your analysis. I recommend fitting a different HMM to each single area, instead of fitting the same model to all neurons from all areas, for several reasons. First, neural activity in some areas might not be well captured by a discrete state sequence - in our paper, we focused on secondary motor cortex, and found really consistent sequences, but in general that might not be the case for other areas… Second, different areas may have different dynamics, with transitions happening at different times, which you’d like to tease apart. Third, when you fit the same Poisson-HMM to all areas, your emissions are very high-dimensional (dim=number of neurons~1000), so the number of parameters is in the same ballpark as the number of observations (=bins per trial x number of trials~500 x 100), which is not the best regime for fitting generative models. But if you fit a single area you’ll be in much better shape.

Did you plot state sequences using MAP, say the Viterbi algorithm? If you did, I suggest instead making another plot, where you look at the posterior probabilities of the states in each trial, and check whether they are “bimodal” namely, zero most of the time, and then close to 1 for many consecutive bins where that state is detected with high confidence. You can then set a “detection threshold” at say 80% probability, and consider only intervals with confident detection. You’re going to be left with “uncertain” intervals, where no state goes above 80% prob. This gives you a measure of how well your model fit is doing, if you have a lot of “uncertain” intervals it means the model does not fit the data very well.

We are already fitting a different HMM to each single area (the results we attached follow this too). However, with all the areas (including secondary motor cortex, striatum, and others in the midbrain) we have fitted until now we have obtained results similar to the ones shown in the previous figure.

For the fitting we are using your matlab code, and we looked directly into the generated plots with the posterior probabilities for each trial. Then, we used the states (hmm_postfit.sequences) with higher confidence to make our plot (in which zero implies no state with confidence at all). For example, in the figure we showed, we have a lot of uncertainty in trials 40-100.

Considering that we were already following your suggestions (perhaps we are still misunderstanding parameters or any other concept), we were worried about the fact of the dynamics related to the task are not captured by the model. Also, due to the computational cost, exploration of different datasets becomes difficult.

We attach one example of the last simulation we run with the MRN area (only including correct trials). At least it seems we are capturing some pre and post-stimulus dynamics, but we still have a lot of uncertainty and variability. If you can see a clear problem that we are missing in the transition matrix (sparseness or whatever) or in any other thing, please let us know. Also, note that we have changed the color coding of the states to show uncertain states in white.

Thank you very much for your help and sorry for replying during the weekend,

Hey Miguel, thanks for the plots! You are already way ahead with the analyses, I think it’d be best to chat on Zoom if you’d like, so you can show me some raw data from MRN, I’m not familiar with its spike trains. Some quick comments. Maybe the model selection is overfitting (too many states), for example the first bunch of states look very similar from the epm. Also, sometimes HMMs don’t work, the underlying dynamics may be not attractor-like, or maybe if you tried a few other sessions it’d work better. Did you try to fit secondary motor cortex instead of MRN?
We could meet Monday or Wed at 10:30am PDT, feel free to send me an email with a Zoom link.

Thank you very much for your help. Monday is perfect for us (we will be working on the project at that time). I already sent you an email with the invitation to my zoom room. We will have some figures prepared to show you.

We are also trying to fit an HMM on the neural activity of different brain areas, using the Steinmetz dataset. We are using the ssm python library (not the matlab code).

After fitting a 4 hidden states HMM on the secondary motor area, we found that the states correspond to different movements of the wheel: one for no movement, one for right, one for left and the last one is active at the beginning of a movement in either direction (see the two trials below, 0s corresponds to stimulus onset).

Whereas in the primary visual area, we found that a 3 hidden states HMM result in one state (state 1 on the figures below) that is correlated with visual stimulus presented on the right side only.

Now we are thinking about some ideas to try to go a bit further. For example, fit a HMM on the prestimulus activity and then try to predict the performance of the mouse from the sequence of hidden states, and compare the results for different brain areas. Or to fit a HMM on the secondary motor area, for passive trials (only small wheel movements, mainly spontaneous activity), and then use it on active trials to try to recover the relationship between hidden states and wheel movement, similar to the one shown above if possible.

Thanks for sharing your data @cloosnathan, it looks very promising. It seems you can clearly predict when the animal is going to start engaging in movement by the transition into the red state. It would be interesting to see what happens when you compare states encoding for the same movements during the passive and during the goal-directed conditions. So, you can either decode the state probabilities of passive trials using the model trained on active, or viceversa. To check that the states encode what you claim, please align each state onset to the behavioral event in each trial, then plot the distribution of intervals across trials between state onset and event. The peak of the distribution for the eg red state should precede the movement onset - and so on for the other states. Please let me know if you’d like to chat on Zoom, eg Wed at 11am PDT works for me.

Thank you for your answer @lmazzuca. The figure below shows histograms of the intervals between the first state onset and the first movement onset in each active trial, for the 4 states HMM fitted on the secondary motor area. So we can observe that the red state is indeed active before wheel movement.

We tried to take the HMM trained on active trials, and apply it on passive trials. The nice relation between hidden state and wheel movement is lost, instead they are fast transitions (compared to active trials) between the green and the red states. So the red state can be active even if no movement follows.