ShanechiLab / PyPSID

Preferential Subspace Identification Algorithm
Other
48 stars 6 forks source link

(I)PSID: (Input) Preferential subspace identification
[Python implementation]

For MATLAB implementation see http://github.com/ShanechiLab/PSID

Given signals y_t (e.g. neural signals) and z_t (e.g behavior), PSID learns a dynamic model for y_t while prioritizing the dynamics that are relevant to z_t.

IPSID is an extension of PSID that also supports taking a third signal u_t (e.g., task instructions) that is simultaneously measured with y_t. In the learned dynamical model, u_t plays the role of input to the latent states.

Publications

PSID

For the derivation of PSID and results in real neural data see the paper below.

Omid G. Sani, Hamidreza Abbaspourazad, Yan T. Wong, Bijan Pesaran, Maryam M. Shanechi. Modeling behaviorally relevant neural dynamics enabled by preferential subspace identification. Nature Neuroscience, 24, 140–149 (2021). https://doi.org/10.1038/s41593-020-00733-0

View-only full-text link: https://rdcu.be/b993t

Original preprint: https://doi.org/10.1101/808154

You can also find a summary of the paper in the following Twitter thread: https://twitter.com/MaryamShanechi/status/1325835609345122304

IPSID

For the derivation of IPSID and results in real neural data see the paper below.

Parsa Vahidi, Omid G. Sani, Maryam M. Shanechi. Modeling and dissociation of intrinsic and input-driven neural population dynamics underlying behavior. PNAS (2024). https://doi.org/10.1073/pnas.2212887121

Usage guide

Installation

Download the source code from the GitHub repository, or install PSID in your Python environment using pip, by running:

pip install PSID --upgrade

You can find the usage license in LICENSE.md.

Initialization

Import the PSID module.

import PSID

Main learning function

The main functions for the Python implementation are the follwing:

A complete usage guide is available in as comments in each function. The following shows example use cases:

idSys = PSID.PSID(y, z, nx, n1, i)
# Or, if modeling effect of input u is also of interest
idSys = PSID.IPSID(y, z, u, nx, n1, i)

Inputs:

Output:

Extracting latent states using learned model

Once a model is learned using (I)PSID, you can apply the model to new data (i.e. run the associated Kalman filter) as follows:

zPred, yPred, xPred = idSys.predict(y)
# Or, for IPSID:
zPred, yPred, xPred = idSys.predict(y, u)

Input:

Outputs:

Required preprocessing

Choosing the hyperparameters

How to pick the state dimensions nx and n1?

nx determines the total dimension of the latent state and n1 determines how many of those dimensions will be prioritizing the inclusion of behaviorally relevant neural dynamics (i.e. will be extracted using stage 1 of (I)PSID). So the values that you would select for these hyperparameters depend on the goal of modeling and on the data. Some examples use cases are:

If you want to perform dimension reduction, nx will be your desired target dimension. For example, to reduce dimension to 2 to plot low-dimensional visualizations of neural activity, you would use nx=2. Now if you want to reduce dimension while preserving as much behaviorally relevant neural dynamics as possible, you would use n1=nx. If you want to find the best fit to data overall, you can perform a grid search over values of nx and n1 and pick the value that achieves the best performance metric in the training data. For example, you could pick the nx and n1 pair that achieves the best cross-validated behavior decoding in an inner-cross-validation within the training data.

How to pick the horizon i?

The horizon i does not affect the model structure and only affects the intermediate linear algebra operations that (I)PSID performs during the learning of the model. Nevertheless, different values of i may have different model learning performance. i needs to be at least 2, but also also determines the maximum n1 and nx that can be used per:

n1 <= nz * i
nx <= ny * i

So if you have a low dimensional y_k or z_k (small ny or nz), you typically would choose larger values for i, and vice versa. It is also possible to select the best performing i via an inner cross-validation approach similar to nx and n1 above. Overall, since i affects the learning performance, it is important for reproducibility that the i that was used is reported.

For more information, see the notebook(s) referenced in the next section.

Usage examples

PSID

Example code for running PSID is provided in source/example/PSID_example.py This script performs PSID model identification and visualizes the learned eigenvalues similar to in Supplementary Fig 1 in (Sani et al, 2021).

The following notebook also contains some examples along with more descriptions: source/example/PSID_tutorial.ipynb

IPSID

Example code for running IPSID is provided in source/example/IPSID_example.py This script performs IPSID model identification and visualizes the learned eigenvalues similar to in Fig. 2A in (Vahidi, Sani, et al, 2024).

The following notebook also contains some examples along with more descriptions: source/example/IPSID_tutorial.ipynb

Change Log

You can see the change log in in ChangeLog.md

Licence

Copyright (c) 2020 University of Southern California
See full notice in LICENSE.md
Omid G. Sani, Parsa Vahidi and Maryam M. Shanechi
Shanechi Lab, University of Southern California