nstat.Analysis

class nstat.Analysis[source]

Bases: object

Collection of static methods for GLM analysis of point-process data.

Every public method is a @staticmethod; the class acts as a pure namespace that mirrors the Matlab @Analysis class. Two naming conventions coexist:

  • PEP 8 (snake_case): run_analysis_for_neuron, run_analysis_for_all_neurons

  • Matlab-facing (camelCase): RunAnalysisForNeuron, RunAnalysisForAllNeurons

See also

Trial, ConfigCollection, SpikeTrainCollection, History

__init__()
colors = ['b', 'g', 'r', 'c', 'm', 'y', 'k']
static psth(spike_trains, bin_edges)[source]

Compute peri-stimulus time histogram. Delegates to module-level psth().

Parameters:
Return type:

tuple[ndarray, ndarray]

static run_analysis_for_neuron(trial, neuron_index, config_collection, *, algorithm='GLM', l2=0.0, max_iter=120)[source]

Run GLM analysis for one neuron across all configurations.

Iterates over the configurations in config_collection, fits a GLM for each, computes KS diagnostics, and returns a single FitResult that aggregates all fits.

Parameters:

  • trial (Trial) – Trial object containing spike trains and covariates.

  • neuron_index (int) – Zero-based neuron index.

  • config_collection (ConfigCollection) – Configurations describing the fits to perform (covariates, history, ensemble history).

  • algorithm ({‘GLM’, ‘BNLRCG’}, default 'GLM') – Regression algorithm.

  • l2 (float, default 0.0) – L-2 regularisation strength.

  • max_iter (int, default 120) – Maximum iterations for the GLM solver.

Returns:

Fit result with KS statistics already populated.

Return type:

FitResult

static run_analysis_for_all_neurons(trial, config_collection, *, algorithm='GLM', l2=0.0, max_iter=120)[source]

Run GLM analysis for every unmasked neuron in the trial.

Calls run_analysis_for_neuron() for each neuron in the trial’s spike-train collection.

Parameters:

  • trial (Trial) – Trial to analyse.

  • config_collection (ConfigCollection) – Configurations describing the fits to perform.

  • algorithm ({‘GLM’, ‘BNLRCG’}, default 'GLM') – Regression algorithm.

  • l2 (float, default 0.0) – L-2 regularisation strength.

  • max_iter (int, default 120) – Maximum iterations for the GLM solver.

Returns:

One FitResult per neuron.

Return type:

list of FitResult

static RunAnalysisForNeuron(tObj, neuronNumber, configColl, makePlot=1, Algorithm='GLM', DTCorrection=1, batchMode=0)[source]

Matlab-facing wrapper for run_analysis_for_neuron().

Parameters:

  • tObj (Trial) – Trial to analyse.

  • neuronNumber (int or str or sequence) – Matlab-style 0-based neuron index, name, or vector of indices. If more than one neuron is specified the return value is a list.

  • configColl (ConfigCollection) – Configurations describing the fits.

  • makePlot (int, default 1) – If 1, plot a summary for the neuron.

  • Algorithm ({‘GLM’, ‘BNLRCG’}, default 'GLM') – Regression algorithm.

  • DTCorrection (int, default 1) – Discrete-time KS correction flag (kept for API parity; unused).

  • batchMode (int, default 0) – Batch-mode flag (kept for API parity; unused).

Returns:

Single result when one neuron is specified, list otherwise.

Return type:

FitResult or list of FitResult

static RunAnalysisForAllNeurons(tObj, configs, makePlot=1, Algorithm='GLM', DTCorrection=1, batchMode=0)[source]

Matlab-facing wrapper for run_analysis_for_all_neurons().

Runs the fits specified by configs on every unmasked neuron in the trial.

Parameters:

  • tObj (Trial) – Trial to analyse.

  • configs (ConfigCollection) – Configurations describing the fits.

  • makePlot (int, default 1) – If 1, generate a summary plot for each neuron.

  • Algorithm ({‘GLM’, ‘BNLRCG’}, default 'GLM') – Regression algorithm.

  • DTCorrection (int, default 1) – Discrete-time KS correction flag (unused).

  • batchMode (int, default 0) – Batch-mode flag (unused).

Returns:

Single result when the trial has one neuron, list otherwise.

Return type:

FitResult or list of FitResult

static GLMFit(tObj, neuronNumber, lambdaIndex, Algorithm='GLM', *, l2=0.0, max_iter=120)[source]

Fit a point-process GLM for a single neuron from a Trial.

Extracts the design matrix X from the current covariate masks, history, and ensemble history in the Trial, and the observation vector Y, then performs the GLM regression.

Parameters:

  • tObj (Trial) – Trial containing spike trains and covariates.

  • neuronNumber (int or str or sequence) – Matlab-style 0-based neuron index, name, or sequence thereof.

  • lambdaIndex (int) – Configuration index used for labelling the returned λ.

  • Algorithm ({‘GLM’, ‘BNLRCG’}, default 'GLM') – 'GLM' — standard Poisson GLM regression. 'BNLRCG' — truncated, L-2 regularised binomial logistic regression (requires binary spike representation).

  • l2 (float, default 0.0) – L-2 regularisation strength.

  • max_iter (int, default 120) – Maximum IRLS / CG iterations.

Returns:

  • lambda_sig (Covariate) – Conditional intensity function evaluated on the design-matrix time grid.

  • b (ndarray) – GLM regression coefficients.

  • dev (float) – Deviance of the fit.

  • stats (dict) – Fit statistics (standard errors, convergence info, covariance matrix).

  • AIC (float) – Akaike information criterion.

  • BIC (float) – Bayesian information criterion.

  • logLL (float) – Log-likelihood evaluated with the fit parameters.

  • distribution (str) – 'poisson' or 'binomial'.

Notes

A rank-deficient design in the unpenalized (l2 = 0) poisson 'GLM' fit is handled as MATLAB glmfit handles it: a column-pivoted QR finds the rank, the fit uses the independent columns, and each dependent column gets coefficient 0 and standard error 0. A design (or spike vector) with a NaN entry has that row dropped before fitting – b, dev and stats["se"] come from the NaN-row-dropped fit, as MATLAB glmfit’s statremovenan does – but lambda_signal (and so AIC, BIC and logLL, which read it) is evaluated over the original full design, as MATLAB’s Analysis.GLMFit does (data = exp(X*b)); a NaN row there still gives NaN in lambda_signal at that row, though logLL’s own eps floor (MATLAB’s NaN-ignoring max) keeps that row’s contribution to logLL finite. Only the poisson 'GLM' path mirrors statremovenan; see the binomial note below.

The binomial 'BNLRCG' fit is NOT a MATLAB mirror on a rank-deficient design, and deliberately so: MATLAB’s bnlrCG has no rank handling at all, and its standard errors there come out complex (an eigenvalue-clipped inverse) – itself a MATLAB defect, not something to reproduce. Rather than leave the previous Python-only behavior (SEs from the singular inv(X'WX) with negative variances clipped to 0, letting different – and larger, up to ~2.8e3 in the EM GLM M-step – coefficients pass the M-step’s se < 100 filter), a rank-deficient binomial design now gets the SAME column-pivoted-QR handling as the poisson path above: dependent columns get b = 0, se = 0. On full-rank designs the binomial fit reaches the MLE, which MATLAB’s truncated bnlrCG approaches (about 1e-4 apart on the EM GLM M-step gold). See parity/matlab_defects.yml (em-binomial-glm-bnlrcg, binomial-rank-deficiency-improvement).

static plotInvGausTrans(fitResults, makePlot=0)[source]

Compute and optionally plot the inverse-Gaussian transform ACF.

Parameters:

  • fitResults (FitResult) – Fit result to compute the transform for.

  • makePlot (int, default 0) – If 1, generate the ACF plot.

Returns:

Plot handles (empty list when makePlot is 0).

Return type:

list

static plotFitResidual(fitResults, windowSize=0.01, makePlot=1)[source]

Compute and plot the point-process residual.

Parameters:

  • fitResults (FitResult) – Fit result to compute the residual for.

  • windowSize (float, default 0.01) – Integration window size (seconds).

  • makePlot (int, default 1) – If 1, generate the residual plot.

Returns:

Plot handles (empty list when makePlot is 0).

Return type:

list

static KSPlot(fitResults, DTCorrection=1, makePlot=1)[source]

Compute KS statistics and optionally generate the KS plot.

Parameters:

  • fitResults (FitResult) – Fit result to compute KS statistics for.

  • DTCorrection (int, default 1) – Discrete-time correction flag.

  • makePlot (int, default 1) – If 1, generate the KS plot.

Returns:

Plot handles (empty list when makePlot is 0).

Return type:

list

static plotSeqCorr(fitResults)[source]

Plot the sequential correlation of rescaled ISIs (z_j vs z_{j-1}).

Parameters:

fitResults (FitResult) – Fit result (inverse-Gaussian transform is computed if needed).

Returns:

Plot handles.

Return type:

list

static plotCoeffs(fitResults)[source]

Plot regression coefficients for all fits in fitResults.

Parameters:

fitResults (FitResult) – Fit result whose coefficients to plot.

Returns:

Plot handles.

Return type:

list

static computeInvGausTrans(Z)[source]

Compute the inverse-Gaussian transform of rescaled spike times.

Transforms rescaled spike times Z to uniform(0, 1) via U = 1 − exp(−Z), then applies the inverse-Gaussian (probit) transform X = Φ⁻¹(U). The autocorrelation of X is used to test for independence of the rescaled ISIs (a condition for the time-rescaling theorem).

Parameters:

Z (array_like) – Rescaled spike times (exponential rate-1 under H₀).

Returns:

  • X (ndarray) – Inverse-Gaussian transformed values.

  • rhoSig (SignalObj) – Autocorrelation function of X.

  • confBoundSig (SignalObj) – ±1.96 / √N confidence bounds for zero autocorrelation.

static computeKSStats(nspikeObj, lambdaInput, DTCorrection=1, *, random_values=None)[source]

Compute KS goodness-of-fit statistics via the time-rescaling theorem.

Given a neural spike train and a candidate conditional intensity function, computes the rescaled ISIs and the KS plot data.

Parameters:

  • nspikeObj (nspikeTrain or SpikeTrainCollection or sequence) – Neural spike train(s).

  • lambdaInput (Covariate) – Candidate conditional intensity function.

  • DTCorrection (int, default 1) – If 1, apply discrete-time correction to KS plot.

  • random_values (array_like, optional) – Pre-drawn uniform random values for reproducibility.

Returns:

  • Z (ndarray) – Rescaled spike times.

  • U (ndarray) – Z transformed to uniform(0, 1).

  • xAxis (ndarray) – x-axis of the KS plot.

  • KSSorted (ndarray) – Sorted rescaled times (y-axis of KS plot).

  • ks_stat (ndarray) – KS statistic — maximum deviation from the 45° line for each conditional intensity function.

static ksdiscrete(pk, st, spikeflag, *, random_values=None)[source]

Thin static wrapper around the file-local ksdiscrete helper.

Mirrors the MATLAB Analysis.ksdiscrete static method (see Analysis.m), which exposes the file-local discrete-time time-rescaling helper for unit testing. Production code calls the same helper via computeKSStats(); both routes share the same implementation.

Parameters:

  • pk (array_like) – Discrete-time conditional probabilities (each value in [0, 1]).

  • st (array_like) – Spike train (binary vector) or spike-index vector, selected by spikeflag.

  • spikeflag ({'spiketrain', 'spikeind'}) – Interpret st as a binary spike train or as spike indices.

  • random_values (array_like, optional) – Uniform random draws used by the discrete-time rescaling. When None the routine draws its own with numpy.random.random_sample().

Returns:

rst, rstsort, xks, cb, rstoldsort – Same five outputs as the MATLAB function: unsorted rescaled times, sorted rescaled times, x-axis values for the KS plot, the 95% confidence-bound scalar, and sorted rescaled times without the discrete-time correction.

Return type:

tuple of numpy.ndarray

References

Haslinger, Pipa & Brown (2010), discrete-time time-rescaling correction.

static computeFitResidual(nspikeObj, lambdaInput, windowSize=0.01)[source]

Compute the point-process residual.

Defined as the difference between the observed spike count and the integral of the candidate conditional intensity function in each time window, following Truccolo et al. (2005).

Parameters:

  • nspikeObj (nspikeTrain or SpikeTrainCollection or sequence) – Neural spike train(s).

  • lambdaInput (Covariate) – Candidate conditional intensity function.

  • windowSize (float, default 0.01) – Size of the integration window (seconds).

Returns:

Point-process residual M(t_k).

Return type:

Covariate

References

Truccolo, W., Eden, U. T., Fellows, M. R., Donoghue, J. P., & Brown, E. N. (2005). J Neurophysiol, 93(2), 1074–1089.

static compHistEnsCoeffForAll(tObj, history, makePlot=1)[source]

Compute ensemble-history coefficients for all unmasked neurons.

Calls compHistEnsCoeff() for each neuron that is not masked.

Parameters:

  • tObj (Trial) – Trial to analyse.

  • history (History) – History object defining the window structure.

  • makePlot (int, default 1) – Summary plot flag.

Returns:

  • fit_results (list of FitResult) – One fit result per neuron.

  • ensemble_cov (CovariateCollection or None) – Ensemble covariates from the last neuron.

  • config_collections (list of ConfigCollection) – Configuration collections used per neuron.

static compHistEnsCoeff(tObj, history, neuronNum=None, neighbors=None, ensembleCov=None, makePlot=1)[source]

Compute ensemble-history coefficients for one neuron.

Builds a covariate collection from the spiking history of neighbouring neurons and fits a GLM with ensemble history as the design matrix.

Parameters:

  • tObj (Trial) – Trial containing spike trains and covariates.

  • history (History) – History object defining the window structure.

  • neuronNum (int or None) – Matlab-style 0-based neuron index. Defaults to the first unmasked neuron.

  • neighbors (array_like or None) – Indices of neighbouring neurons. Defaults to tObj.getNeuronNeighbors(neuronNum).

  • ensembleCov (CovariateCollection or None) – Pre-computed ensemble covariates. If None, computed automatically.

  • makePlot (int, default 1) – Summary plot flag.

Returns:

  • fitResults (FitResult) – Fit result for the ensemble-history model.

  • ensembleCov (CovariateCollection) – Ensemble covariates used in the fit.

  • tcc (ConfigCollection) – Configuration collection used.

static computeGrangerCausalityMatrix(tObj, Algorithm='GLM', confidenceInterval=0.95, batchMode=0)[source]

Compute the Granger-causality matrix for the neural ensemble.

For every pair of neurons, fits a baseline model (full ensemble history) and a reduced model (one neighbour excluded), then computes the log-likelihood ratio. Statistical significance is corrected for multiple comparisons with Benjamini–Hochberg FDR.

Parameters:

  • tObj (Trial) – Trial with ensemble history configured.

  • Algorithm ({‘GLM’, ‘BNLRCG’}, default 'GLM') – Regression algorithm.

  • confidenceInterval (float, default 0.95) – Confidence level for the significance test.

  • batchMode (int) – Unused (Matlab API parity).

Returns:

  • fitResults (list of list of FitResult) – fitResults[i][j] is the fit result for the test of neighbour j → neuron i.

  • gammaMat (ndarray, shape (N, N)) – Log-likelihood ratio Γ matrix.

  • phiMat (ndarray, shape (N, N)) – Signed Γ matrix (sign from sum of excluded coefficients).

  • devianceMat (ndarray, shape (N, N)) – Deviance (−2Γ) matrix.

  • sigMat (ndarray, shape (N, N)) – Binary significance matrix after FDR correction.

static computeHistLag(tObj, neuronNum=None, windowTimes=None, CovLabels=None, Algorithm='GLM', batchMode=0, sampleRate=None, makePlot=1, histMinTimes=None, histMaxTimes=None)[source]

Sweep self-history window orders for a single neuron.

Fits a sequence of GLMs with increasing numbers of history windows (no extrinsic covariates, no ensemble history) and returns the fit results for model selection via AIC / BIC / KS.

Parameters:

  • tObj (Trial) – Trial to analyse.

  • neuronNum (int or None) – Matlab-style 0-based neuron index. If None, uses the first unmasked neuron.

  • windowTimes (array_like) – Vector of window boundary times. len(windowTimes) - 1 configurations are created with increasing history order.

  • CovLabels (list of str or None) – Covariate labels to include in each configuration.

  • Algorithm ({‘GLM’, ‘BNLRCG’}, default 'GLM') – Regression algorithm.

  • batchMode (int) – Unused (Matlab API parity).

  • sampleRate (float or None) – Sample rate override; defaults to tObj.sampleRate.

  • makePlot (int, default 1) – If 1, generate a summary plot.

  • histMinTimes (float or None) – Optional time bounds passed to the History object.

  • histMaxTimes (float or None) – Optional time bounds passed to the History object.

Returns:

  • fitResults (FitResult) – Fit result containing all history-order configurations.

  • tcc (ConfigCollection) – The generated configuration collection.

static computeHistLagForAll(tObj, windowTimes, CovLabels=None, Algorithm='GLM', batchMode=0, sampleRate=None, makePlot=1, histMinTimes=None, histMaxTimes=None)[source]

Sweep self-history window orders for all unmasked neurons.

Calls computeHistLag() for each unmasked neuron in the trial.

Parameters:

  • tObj (Trial) – Trial to analyse.

  • windowTimes (array_like) – Vector of window boundary times.

  • CovLabels (list of str or None) – Covariate labels for each configuration.

  • Algorithm ({‘GLM’, ‘BNLRCG’}, default 'GLM') – Regression algorithm.

  • batchMode (int) – Unused (Matlab API parity).

  • sampleRate (float or None) – Sample rate override.

  • makePlot (int, default 1) – Summary plot flag.

  • histMinTimes (float or None) – Optional time bounds for the History object.

  • histMaxTimes (float or None) – Optional time bounds for the History object.

Returns:

One fit result per unmasked neuron.

Return type:

list of FitResult

static computeNeighbors(tObj, neuronNum=None, sampleRate=None, windowTimes=None, makePlot=1)[source]

Sweep ensemble-history orders for one neuron (no self-history).

Fits a sequence of GLMs with increasing ensemble-history window orders but no self-history and no extrinsic covariates, for model selection of the ensemble effect.

Parameters:

  • tObj (Trial) – Trial to analyse.

  • neuronNum (int or None) – Matlab-style 0-based neuron index.

  • sampleRate (float or None) – Sample rate override.

  • windowTimes (array_like) – Vector of window boundary times.

  • makePlot (int, default 1) – Summary plot flag.

Returns:

  • fitResults (FitResult) – Fit result with all configurations.

  • tcc (ConfigCollection) – Generated configuration collection.

static spikeTrigAvg(tObj, neuronNum, windowSize)[source]

Compute the spike-triggered average of all covariates.

Each covariate dimension is sampled at every spike time of the specified neuron ± windowSize / 2. The returned collection contains one covariate per original dimension, where each column corresponds to a single spike. Use plotVariability on the returned signals to visualise the average and spread.

Parameters:

  • tObj (Trial) – Trial containing spike trains and covariates.

  • neuronNum (int) – Matlab-style 0-based neuron index.

  • windowSize (float) – Total window length (seconds) centred on each spike.

Returns:

Collection of spike-triggered covariate samples.

Return type:

CovariateCollection