hnn_core.dipole.Dipole

class hnn_core.dipole.Dipole(times, data, nave=1)[source]

Dipole class.

An instance of the Dipole-class contains the simulated dipole moment timecourses for L2 and L5 pyramidal cells, as well as their aggregate ('agg'). The units of the dipole moment are in nAm (1e-9 Ampere-meters).

Parameters:
timesarray (n_times,)

The time vector (in ms)

dataarray, shape (n_times x n_layers)

The data. The first column represents ‘agg’ (the total diple), the second ‘L2’ layer and the last one ‘L5’ layer. For experimental data, it can contain only one column.

naveint

Number of trials that were averaged to produce this Dipole. Defaults to 1

Attributes:
timesarray-like

The time vector (in ms)

sfreqfloat

The sampling frequency (in Hz)

datadict of array

Dipole moment timecourse arrays with keys ‘agg’, ‘L2’ and ‘L5’

naveint

Number of trials that were averaged to produce this Dipole

scale_appliedint or float

The total factor by which the dipole has been scaled (using scale()).

Methods

copy()

Return a copy of the Dipole instance

plot([tmin, tmax, layer, decim, ax, color, show])

Simple layer-specific plot function.

plot_psd([fmin, fmax, tmin, tmax, layer, ...])

Plot power spectral density (PSD) of dipole time course

plot_tfr_morlet(freqs[, n_cycles, tmin, ...])

Plot Morlet time-frequency representation of dipole time course

savgol_filter(h_freq)

Smooth the dipole waveform using Savitzky-Golay filtering

scale(factor)

Scale (multiply) the dipole moment by a fixed factor

smooth(window_len)

Smooth the dipole waveform using Hamming-windowed convolution

write(fname[, overwrite])

Write dipole values to a txt or hdf5 file.

copy()[source]

Return a copy of the Dipole instance

Returns:
dpl_copyinstance of Dipole

A copy of the Dipole instance.

plot(tmin=None, tmax=None, layer='agg', decim=None, ax=None, color='k', show=True)[source]

Simple layer-specific plot function.

Parameters:
layerstr

The layer to plot. Can be one of ‘agg’, ‘L2’, and ‘L5’

decimateint

Factor by which to decimate the raw dipole traces (optional)

axinstance of matplotlib figure | None

The matplotlib axis

colortuple of float

RGBA value to use for plotting. By default, ‘k’ (black)

showbool

If True, show the figure

Returns:
figinstance of plt.fig

The matplotlib figure handle.

plot_psd(fmin=0, fmax=None, tmin=None, tmax=None, layer='agg', color=None, label=None, ax=None, show=True)[source]

Plot power spectral density (PSD) of dipole time course

Applies periodogram from SciPy with window='hamming'. Note that no spectral averaging is applied across time, as most hnn_core simulations are short-duration. However, passing a list of Dipole instances will plot their average (Hamming-windowed) power, which resembles the Welch-method applied over time.

Parameters:
dplinstance of Dipole | list of Dipole instances

The Dipole object.

fminfloat

Minimum frequency to plot (in Hz). Default: 0 Hz

fmaxfloat

Maximum frequency to plot (in Hz). Default: None (plot up to Nyquist)

tminfloat or None

Start time of data to include (in ms). If None, use entire simulation.

tmaxfloat or None

End time of data to include (in ms). If None, use entire simulation.

layerstr, default ‘agg’

The layer to plot. Can be one of ‘agg’, ‘L2’, and ‘L5’

colorstr | tuple | None

The line color of PSD

labelstr | None

Line label for PSD

axinstance of matplotlib figure | None

The matplotlib axis.

showbool

If True, show the figure

Returns:
figinstance of matplotlib Figure

The matplotlib figure handle.

plot_tfr_morlet(freqs, n_cycles=7.0, tmin=None, tmax=None, layer='agg', decim=None, padding='zeros', ax=None, colormap='inferno', colorbar=True, colorbar_inside=False, show=True)[source]

Plot Morlet time-frequency representation of dipole time course

NB: Calls tfr_array_morlet, so mne must be installed.

Parameters:
dplinstance of Dipole | list of Dipole instances

The Dipole object. If a list of dipoles is given, the power is calculated separately for each trial, then averaged.

freqsarray

Frequency range of interest.

n_cyclesfloat or array of float, default 7.0

Number of cycles. Fixed number or one per frequency.

tminfloat or None

Start time of plot in milliseconds. If None, plot entire simulation.

tmaxfloat or None

End time of plot in milliseconds. If None, plot entire simulation.

layerstr, default ‘agg’

The layer to plot. Can be one of ‘agg’, ‘L2’, and ‘L5’

decimint or list of int or None (default)

Optional (integer) factor by which to decimate the raw dipole traces. The SciPy function decimate() is used, which recommends values <13. To achieve higher decimation factors, a list of ints can be provided. These are applied successively.

paddingstr or None

Optional padding of the dipole time course beyond the plotting limits. Possible values are: ‘zeros’ for padding with 0’s (default), ‘mirror’ for mirror-image padding.

axinstance of matplotlib figure | None

The matplotlib axis

colormapstr

The name of a matplotlib colormap, e.g., ‘viridis’. Default: ‘inferno’

colorbarbool

If True (default), adjust figure to include colorbar.

colorbar_inside: bool, default False

Put the color inside the heatmap if True.

showbool

If True, show the figure

Returns:
figinstance of matplotlib Figure

The matplotlib figure handle.

savgol_filter(h_freq)[source]

Smooth the dipole waveform using Savitzky-Golay filtering

Note that this method operates in-place, i.e., it will alter the data. If you prefer a filtered copy, consider using the copy()-method. The high-frequency cutoff value of a Savitzky-Golay filter is approximate; see the SciPy reference: savgol_filter().

Parameters:
h_freqfloat or None

Approximate high cutoff frequency in Hz. Note that this is not an exact cutoff, since Savitzky-Golay filtering is done using polynomial fits instead of FIR/IIR filtering. This parameter is thus used to determine the length of the window over which a 5th-order polynomial smoothing is applied.

Returns:
dpl_copyinstance of Dipole

A copy of the modified Dipole instance.

scale(factor)[source]

Scale (multiply) the dipole moment by a fixed factor

The attribute Dipole.scale_applied is updated to reflect factors applied and displayed in plots.

Parameters:
factorint

Scaling factor, applied to the data in-place.

smooth(window_len)[source]

Smooth the dipole waveform using Hamming-windowed convolution

Note that this method operates in-place, i.e., it will alter the data. If you prefer a filtered copy, consider using the copy()-method.

Parameters:
window_lenfloat

The length (in ms) of a hamming window to convolve the data with.

Returns:
dpl_copyinstance of Dipole

A copy of the modified Dipole instance.

write(fname, overwrite=True)[source]

Write dipole values to a txt or hdf5 file.

Parameters:
fnamestr | Path object

Full path to the output file (.txt or .hdf5)