pynapple.process.tuning_curves#
Functions to compute n-dimensional tuning curves.
Functions
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
Computes mutual information from n-dimensional tuning curves. |
|
Compute mean response per epoch, given a dictionary of epochs. |
|
Computes n-dimensional tuning curves relative to n features. |
- pynapple.process.tuning_curves.compute_1d_mutual_info(tc, feature, ep=None, minmax=None, bitssec=False)[source]#
Deprecated since version 0.9.2: compute_1d_mutual_info will be removed in Pynapple 1.0.0, it is replaced by compute_mutual_information because the latter works for N dimensions.
- pynapple.process.tuning_curves.compute_1d_tuning_curves(group, feature, nb_bins, ep=None, minmax=None)[source]#
Deprecated since version 0.9.2: compute_1d_tuning_curves will be removed in Pynapple 1.0.0, it is replaced by compute_tuning_curves because the latter works for N dimensions.
- pynapple.process.tuning_curves.compute_1d_tuning_curves_continuous(tsdframe, feature, nb_bins, ep=None, minmax=None)[source]#
Deprecated since version 0.9.2: compute_1d_tuning_curves will be removed in Pynapple 1.0.0, it is replaced by compute_tuning_curves because the latter works for N dimensions and continuous data.
- pynapple.process.tuning_curves.compute_2d_mutual_info(dict_tc, features, ep=None, minmax=None, bitssec=False)[source]#
Deprecated since version 0.9.2: compute_2d_mutual_info will be removed in Pynapple 1.0.0, it is replaced by compute_mutual_information because the latter works for N dimensions.
- pynapple.process.tuning_curves.compute_2d_tuning_curves(group, features, nb_bins, ep=None, minmax=None)[source]#
Deprecated since version 0.9.2: compute_1d_tuning_curves will be removed in Pynapple 1.0.0, it is replaced by compute_tuning_curves because the latter works for N dimensions.
- pynapple.process.tuning_curves.compute_2d_tuning_curves_continuous(tsdframe, features, nb_bins, ep=None, minmax=None)[source]#
Deprecated since version 0.9.2: compute_1d_tuning_curves will be removed in Pynapple 1.0.0, it is replaced by compute_tuning_curves because the latter works for N dimensions and continuous data.
- pynapple.process.tuning_curves.compute_discrete_tuning_curves(group, dict_ep)[source]#
Deprecated since version 0.9.2: compute_discrete_tuning_curves will be removed in Pynapple 1.0.0, it is replaced by compute_response_per_epoch.
- pynapple.process.tuning_curves.compute_mutual_information(tuning_curves)[source]#
Computes mutual information from n-dimensional tuning curves.
This function implements Skaggs et al.’s [1] metric to quantify the information content of a neuron’s firing with respect to a variable (e.g., position), based on its tuning curve.
The mutual information in bits per second is given by:
\[I_\text{bits/s} = \sum_x P(x) \lambda(x) \log_2 \left( \frac{\lambda(x)}{\bar{\lambda}} \right)\]where:
\(P(x)\) is the probability of being in bin \(x\) (occupancy),
\(\lambda(x)\) is the firing rate of the neuron in bin \(x\),
\(\bar{\lambda}\) is the overall mean firing rate.
The information per spike is computed by dividing the result by the mean firing rate:
\[I_\text{bits/spike} = \frac{I}{\bar{\lambda}}\]References
- Parameters:
tuning_curves (xarray.DataArray) – Tuning curves as computed by
compute_tuning_curves().- Returns:
A table containing the spatial information per unit, in both bits/sec and bits/spike.
- Return type:
Examples
We can compute the mutual information between a variable and a set of neurons’ firing from the tuning curves:
>>> import pynapple as nap >>> import numpy as np; np.random.seed(42) >>> epoch = nap.IntervalSet([0, 100]) >>> t = np.arange(0, 100, 0.01) >>> feature = nap.Tsd(t=t, d=np.clip(t*0.01 + np.random.normal(0, 0.02, len(t)), 0, 1), time_support=epoch) >>> group = nap.TsGroup({ ... 1: nap.Ts(t[(feature.values >= 0.2) & (feature.values < 0.3)]), ... 2: nap.Ts(t[(feature.values >= 0.7) & (feature.values < 0.8)]) ... }, time_support=epoch) >>> tcs = nap.compute_tuning_curves(group, feature, bins=10) >>> tcs <xarray.DataArray (unit: 2, 0: 10)> Size: 160B array([[ 0., 0., 100., 0., 0., 0., 0., 0., 0., 0.], [ 0., 0., 0., 0., 0., 0., 0., 100., 0., 0.]]) Coordinates: * unit (unit) int64 16B 1 2 * 0 (0) float64 80B 0.05 0.15 0.25 0.35 0.45 0.55 0.65 0.75 0.85 0.95 Attributes: occupancy: [ 985. 1009. 1014. 996. 993. 1008. 991. 1008. 999. 997.] bin_edges: [array([0. , 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1. ])] >>> MI = nap.compute_mutual_information(tcs) >>> MI bits/sec bits/spike 1 33.480966 3.301870 2 33.369159 3.310432
- pynapple.process.tuning_curves.compute_response_per_epoch(data, epochs_dict, return_pandas=False)[source]#
Compute mean response per epoch, given a dictionary of epochs.
- Parameters:
- Returns:
A tensor containing the tuning curves with labeled epochs.
- Return type:
Examples
This function is typically used for a set of discrete stimuli being presented for multiple epochs. The stimulus epochs can overlap, though note that epochs within an IntervalSet can not overlap.
>>> import pynapple as nap >>> import numpy as np; np.random.seed(42) >>> epochs_dict = { ... "stim0": nap.IntervalSet(start=0, end=30), ... "stim1":nap.IntervalSet(start=60, end=90) ... } >>> group = nap.TsGroup({ ... 1: nap.Ts(np.arange(0, 100, 0.1)), ... 2: nap.Ts(np.arange(0, 100, 0.2)) ... }) >>> tcs = nap.compute_response_per_epoch(group, epochs_dict) >>> tcs <xarray.DataArray (unit: 2, epochs: 2)> Size: 32B array([[10.03333333, 10.03333333], [ 5.03333333, 5.03333333]]) Coordinates: * unit (unit) int64 16B 1 2 * epochs (epochs) <U5 40B 'stim0' 'stim1'
You can also pass a TsdFrame (e.g. calcium imaging data), in that case the response is computed:
>>> frame = nap.TsdFrame(d=np.random.rand(2000, 3), t=np.arange(0, 100, 0.05)) >>> tcs = nap.compute_response_per_epoch(frame, epochs_dict) >>> tcs <xarray.DataArray (unit: 3, epochs: 2)> Size: 48B array([[0.50946668, 0.50897635], [0.48343249, 0.48191892], [0.50063158, 0.48748094]]) Coordinates: * unit (unit) int64 24B 0 1 2 * epochs (epochs) <U5 40B 'stim0' 'stim1'
- pynapple.process.tuning_curves.compute_tuning_curves(data, features, bins=10, range=None, epochs=None, fs=None, feature_names=None, return_pandas=False, return_counts=False)[source]#
Computes n-dimensional tuning curves relative to n features.
- Parameters:
data (TsGroup, TsdFrame, Ts, Tsd) – The data for which the tuning curves will be computed. This usually corresponds to the activity of the neurons, either as spike times (TsGroup or Ts) or continuous values (TsdFrame or Tsd).
features (Tsd, TsdFrame) – The features (i.e. one column per feature). This usually corresponds to behavioral variables such as position, head direction, speed, etc.
bins (sequence or int) –
The bin specification:
A sequence of arrays describing the monotonically increasing bin edges along each dimension.
The number of bins for each dimension (nx, ny, … =bins)
The number of bins for all dimensions (nx=ny=…=bins).
range (sequence, optional) – A sequence of entries per feature, each an optional (lower, upper) tuple giving the outer bin edges to be used if the edges are not given explicitly in bins. An entry of None in the sequence results in the minimum and maximum values being used for the corresponding dimension. The default, None, is equivalent to passing a tuple of D None values.
epochs (IntervalSet, optional) – The epochs on which tuning curves are computed. If None, the epochs are the time support of the features.
fs (float, optional) – The exact sampling frequency of the features used to normalise the tuning curves. Unit should match that of the features. If not passed, it is estimated.
feature_names (list, optional) – A list of feature names. If not passed, the column names in features are used.
return_pandas (bool, optional) – If True, the function returns a pandas.DataFrame instead of an xarray.DataArray. Note that this will not work if the features are not 1D and that occupancy and bin edges will not be stored as attributes.
return_counts (bool, optional) – If True, does not divide the spike counts by occupancy, but returns the counts directly. The occupancy is stored in the xarray attributes, so the division can be performed after any particular processing steps. If the input is a TsdFrame, this does not do anything.
- Returns:
A tensor containing the tuning curves with labeled bin centres. The bin edges and occupancy are stored as attributes.
- Return type:
Examples
In the simplest case, we can pass a group of spikes per neuron and a single feature:
>>> import pynapple as nap >>> import numpy as np; np.random.seed(42) >>> group = nap.TsGroup({ ... 1: nap.Ts(np.arange(0, 100, 0.1)), ... 2: nap.Ts(np.arange(0, 100, 0.2)) ... }) >>> feature = nap.Tsd(d=np.arange(0, 100, 0.1) % 1, t=np.arange(0, 100, 0.1)) >>> tcs = nap.compute_tuning_curves(group, feature, bins=10) >>> tcs <xarray.DataArray (unit: 2, 0: 10)> Size: 160B array([[10., 10., 10., 10., 10., 10., 10., 10., 10., 10.], [10., 0., 10., 0., 10., 0., 10., 0., 10., 0.]]) Coordinates: * unit (unit) int64 16B 1 2 * 0 (0) float64 80B 0.045 0.135 0.225 0.315 ... 0.585 0.675 0.765 0.855 Attributes: occupancy: [100. 100. 100. 100. 100. 100. 100. 100. 100. 100.] bin_edges: [array([0. , 0.09, 0.18, 0.27, 0.36, 0.45, 0.54, 0.63, 0.72,...
The function can also take multiple features, in which case it computes n-dimensional tuning curves. We can specify the number of bins for each feature:
>>> features = nap.TsdFrame( ... d=np.stack( ... [ ... np.arange(0, 100, 0.1) % 1, ... np.arange(0, 100, 0.1) % 2 ... ], ... axis=1 ... ), ... t=np.arange(0, 100, 0.1) ... ) >>> tcs = nap.compute_tuning_curves(group, features, bins=[5, 3]) >>> tcs <xarray.DataArray (unit: 2, 0: 5, 1: 3)> Size: 240B array([[[10., 10., nan], [10., 10., 10.], [10., nan, 10.], [10., 10., 10.], [nan, 10., 10.]], ... [[ 5., 5., nan], [ 5., 10., 0.], [ 5., nan, 5.], [10., 0., 5.], [nan, 5., 5.]]]) Coordinates: * unit (unit) int64 16B 1 2 * 0 (0) float64 40B 0.09 0.27 0.45 0.63 0.81 * 1 (1) float64 24B 0.3167 0.95 1.583 Attributes: occupancy: [[100. 100. nan]\n [100. 50. 50.]\n [100. nan 100.]\n [ 5... bin_edges: [array([0. , 0.18, 0.36, 0.54, 0.72, 0.9 ]), array([0. ...
Or even specify the bin edges directly:
>>> tcs = nap.compute_tuning_curves( ... group, ... features, ... bins=[np.linspace(0, 1, 5), np.linspace(0, 2, 3)] ... ) >>> tcs <xarray.DataArray (unit: 2, 0: 4, 1: 2)> Size: 128B array([[[10. , 10. ], [10. , 10. ], [10. , 10. ], [10. , 10. ]], ... [[ 6.66666667, 6.66666667], [ 5. , 5. ], [ 3.33333333, 3.33333333], [ 5. , 5. ]]]) Coordinates: * unit (unit) int64 16B 1 2 * 0 (0) float64 32B 0.125 0.375 0.625 0.875 * 1 (1) float64 16B 0.5 1.5 Attributes: occupancy: [[150. 150.]\n [100. 100.]\n [150. 150.]\n [100. 100.]] bin_edges: [array([0. , 0.25, 0.5 , 0.75, 1. ]), array([0., 1., 2.])]
In all of these cases, it is also possible to pass continuous values instead of spikes (e.g. calcium imaging data), in that case the mean response is computed:
>>> frame = nap.TsdFrame(d=np.random.rand(2000, 3), t=np.arange(0, 100, 0.05)) >>> tcs = nap.compute_tuning_curves(frame, feature, bins=10) >>> tcs <xarray.DataArray (unit: 3, 0: 10)> Size: 240B array([[0.49147343, 0.50190395, 0.50971339, 0.50128013, 0.54332711, 0.49712328, 0.49594611, 0.5110517 , 0.52247351, 0.52057658], [0.51132036, 0.46410557, 0.47732505, 0.49830908, 0.53523019, 0.53099429, 0.48668499, 0.44198555, 0.49222208, 0.47453398], [0.46591801, 0.50662914, 0.46875882, 0.48734997, 0.51836574, 0.50722266, 0.48943577, 0.49730095, 0.47944075, 0.48623693]]) Coordinates: * unit (unit) int64 24B 0 1 2 * 0 (0) float64 80B 0.045 0.135 0.225 0.315 ... 0.585 0.675 0.765 0.855 Attributes: occupancy: [100. 100. 100. 100. 100. 100. 100. 100. 100. 100.] bin_edges: [array([0. , 0.09, 0.18, 0.27, 0.36, 0.45, 0.54, 0.63, 0.72,...