API Reference

multitaper_spectrogram(data, fs, time_step, window_length=None, NW=4.0, n_tapers=None, freq_range=None, weight_type='unity', detrend='constant', nfft=None, db_scale=True, p_ref=2e-05, boundary_pad=False)

Compute the multitaper PSD of the input data.

Parameters:

Name	Type	Description	Default
data	NDArray	(n_samples,) Input data	required
fs	float	Sampling frequency	required
time_step	float	Time step between frames in seconds	required
window_length	float	Window length in seconds. If None, will be set to the same as time_step. Defaults to None.	None
NW	float	NW value, see notes for details. Defaults to 4.0.	4.0
n_tapers	Optional[int]	The max number of tapers, if None, will be set to NW*2-1. Defaults to None.	None
freq_range	Optional[list]	The desired frequency range. If None, will be set to [0, fs/2]. Defaults to None.	None
weight_type	Literal['unity', 'eig']	The type of weights among tapers. Defaults to "unity".	'unity'
detrend	Literal['constant', 'linear', 'off']	Whether and how to detrend the signal. Defaults to "constant".	'constant'
nfft	Optional[int]	The number of FFT points. If None, will be set to the smallest power of 2 that is larger than the window length. Defaults to None.	None
db_scale	bool	Whether convert the result to db scale, i.e. 10log10(psd/p_ref**2). Defaults to True.	True
p_ref	float	If db_scale is True, the p_ref value used in the dB conversion. Defaults to 2e-5.	2e-05
boundary_pad	bool	Whether to pad the data with zeros at the beginning and end. This is useful when the data is not evenly divisible by the window length and time step. By default False. If True, the data will be padded with zeros at the beginning and end, so that the first frame is centered on the first sample of data, and all samples are included in (at least) one frame. If False, the first frame is centered at window_length/2 seconds after the first sample, and samples after n_frames*time_step+window_length seconds are ignored.	False

Notes

The value of 2W is the regularization bandwidth. Typically, we choose W to be a small multiple of the fundamental frequency 1/(Ndt) (where N is the number of samples in the data), i.e. W=i/(Ndt). The value of the parameter NW here is in fact the value of i (when dt is seen as 1). There's a trade-off between frequency resolution and variance reduction: A larger NW will reduce the variance of the PSD estimate, but also reduce the frequency resolution.

Returns:

Name	Type	Description
freqs	NDArray	(n_freqs,) Frequency points of the spectrogram
times	NDArray	(n_frames,) Time points of each frame
psd	NDArray	(n_freqs,n_frames) PSD spectrogram

Examples:

>>> times,freqs,psd = multitaper_spectrogram(data,fs,time_step=0.001,window_length=0.005,NW=4)

Source code in src/pymultitaper/spectral.py

def multitaper_spectrogram(data:NDArray,fs:float,time_step:float,window_length:Optional[float]=None,NW:float=4.0,n_tapers:Optional[int]=None,freq_range:Optional[list]=None,weight_type:Literal["unity","eig"]="unity",detrend:Literal["constant","linear","off"]="constant",nfft:Optional[int]=None,db_scale:bool=True,p_ref:float=2e-5,boundary_pad:bool=False)-> Tuple[NDArray,NDArray,NDArray]:
    """
    Compute the multitaper PSD of the input data.

    Args:
        data (NDArray): (n_samples,) Input data
        fs (float): Sampling frequency
        time_step (float): Time step between frames in seconds
        window_length (float, optional): Window length in seconds. If `None`, will be set to the same as `time_step`. Defaults to None.
        NW (float, optional): NW value, see notes for details. Defaults to 4.0.
        n_tapers (Optional[int], optional): The max number of tapers, if `None`, will be set to NW*2-1. Defaults to None.
        freq_range (Optional[list], optional): The desired frequency range. If `None`, will be set to [0, fs/2]. Defaults to None.
        weight_type (Literal["unity","eig"], optional): The type of weights among tapers. Defaults to "unity".
        detrend (Literal["constant","linear","off"], optional): Whether and how to detrend the signal. Defaults to "constant".
        nfft (Optional[int], optional): The number of FFT points. If `None`, will be set to the smallest power of 2 that is larger than the window length. Defaults to None.
        db_scale (bool, optional): Whether convert the result to db scale, i.e. 10log10(psd/p_ref**2). Defaults to True.
        p_ref (float, optional): If `db_scale` is `True`, the `p_ref` value used in the dB conversion. Defaults to 2e-5.
        boundary_pad (bool, optional): Whether to pad the data with zeros at the beginning and end. This is useful when the data is not evenly divisible by the window length and time step. By default `False`.

            - If `True`, the data will be padded with zeros at the beginning and end, so that the first frame is centered on the first sample of data, and all samples are included in (at least) one frame.
            - If `False`, the first frame is centered at `window_length/2` seconds after the first sample, and samples after `n_frames*time_step+window_length` seconds are ignored.

    Notes:
        The value of 2W is the regularization bandwidth. Typically, we choose W to be a small multiple of the fundamental frequency 1/(N*dt) (where N is the number of samples in the data), i.e. W=i/(N*dt). The value of the parameter `NW` here is in fact the value of i (when dt is seen as 1). There's a trade-off between frequency resolution and variance reduction: A larger `NW` will reduce the variance of the PSD estimate, but also reduce the frequency resolution. 

    Returns:
        freqs (NDArray): (n_freqs,) Frequency points of the spectrogram
        times (NDArray): (n_frames,) Time points of each frame
        psd (NDArray): (n_freqs,n_frames) PSD spectrogram

    Examples:
        >>> times,freqs,psd = multitaper_spectrogram(data,fs,time_step=0.001,window_length=0.005,NW=4)
    """
    # (nfft,n_frames)
    if n_tapers is None:
        # Note: NW may be a float number
        n_tapers = np.floor(2*NW-1).astype(int)
    window_length = time_step if window_length is None else window_length
    n_winlen = int(window_length*fs)
    tapers,weights = _get_dpss_windows(n_winlen,NW,n_tapers,weight_type)
    return _spectrogram(data=data,fs=fs,time_step=time_step,win=tapers,weights=weights,freq_range=freq_range,detrend=detrend,nfft=nfft,db_scale=db_scale,p_ref=p_ref,boundary_pad=boundary_pad)

plot_spectrogram(times, freqs, psd, ax=None, **kwargs)

Plot the spectrogram.

Note: Convert the spectrogram to dB scale (set db_scale to True in the spectrogram functions, or convert it manually) before plotting, otherwise the plot may not be very informative.

Parameters:

Name	Type	Description	Default
times	(n_frames,)	Time points of each frame	required
freqs	(n_freqs,)	Frequency points of the spectrogram	required
psd	(n_freqs, n_frames)	PSD spectrogram	required
ax	Optional[Axes]	The Axes object to plot the spectrogram. If None, a new figure will be created. Defaults to None.	None
**kwargs		Additional arguments to ax.pcolormesh	{}

Returns:

Name	Type	Description
fig	Figure	The figure object
ax	Axes	The Axes object

Examples:

>>> f,ax = plt.subplots(1,1)
>>> plot_spectrogram(times,freqs,psd,ax=ax,cmap="viridis")

Source code in src/pymultitaper/spectral.py

def plot_spectrogram(times:NDArray,freqs:NDArray,psd:NDArray,ax:Optional[plt.Axes]=None,**kwargs)-> tuple:
    """
    Plot the spectrogram.

    Note: Convert the spectrogram to dB scale (set `db_scale` to `True` in the spectrogram functions, or convert it manually) before plotting, otherwise the plot may not be very informative.

    Args:
        times (n_frames,): Time points of each frame
        freqs (n_freqs,): Frequency points of the spectrogram
        psd (n_freqs,n_frames): PSD spectrogram
        ax (Optional[plt.Axes], optional): The Axes object to plot the spectrogram. If `None`, a new figure will be created. Defaults to None.
        **kwargs: Additional arguments to `ax.pcolormesh`

    Returns:
        fig (plt.Figure): The figure object
        ax (plt.Axes): The Axes object

    Examples:
        >>> f,ax = plt.subplots(1,1)
        >>> plot_spectrogram(times,freqs,psd,ax=ax,cmap="viridis")
    """
    if ax is None:
        fig,ax = plt.subplots()
    else:
        fig = ax.figure
    mesh = ax.pcolormesh(times,freqs,psd,**kwargs)
    ax.set_xlabel("Time (s)")
    ax.set_ylabel("Frequency (Hz)")
    fig.colorbar(mesh,ax=ax)
    return fig,ax

plot_spectrum(times, freqs, psd, time, ax=None, **kwargs)

Plot the spectrum at a specific time point.

Parameters:

Name	Type	Description	Default
times	(n_frames,)	Time points of each frame	required
freqs	(n_freqs,)	Frequency points of the spectrogram	required
psd	(n_freqs, n_frames)	PSD spectrogram	required
time	float	The time point to plot the spectrum	required
ax	Optional[Axes]	The Axes object to plot the spectrum. If None, a new figure will be created. Defaults to None.	None
**kwargs		Additional arguments to ax.plot	{}

Returns:

Name	Type	Description
fig	Figure	The figure object
ax	Axes	The Axes object

Examples:

>>> f,ax = plt.subplots(1,1)
>>> plot_spectrum(times,freqs,psd,time=0.7,ax=ax)

Source code in src/pymultitaper/spectral.py

def plot_spectrum(times:NDArray,freqs:NDArray,psd:NDArray,time:float,ax:Optional[plt.Axes]=None,**kwargs)-> tuple:
    """

    Plot the spectrum at a specific time point.

    Args:
        times (n_frames,): Time points of each frame
        freqs (n_freqs,): Frequency points of the spectrogram
        psd (n_freqs,n_frames): PSD spectrogram
        time (float): The time point to plot the spectrum
        ax (Optional[plt.Axes], optional): The Axes object to plot the spectrum. If `None`, a new figure will be created. Defaults to None.
        **kwargs: Additional arguments to `ax.plot`

    Returns:
        fig (plt.Figure): The figure object
        ax (plt.Axes): The Axes object

    Examples:
        >>> f,ax = plt.subplots(1,1)
        >>> plot_spectrum(times,freqs,psd,time=0.7,ax=ax)
    """
    if ax is None:
        fig,ax = plt.subplots()
    else:
        fig = ax.figure
    idx = np.argmin(np.abs(times-time))
    ax.plot(freqs,psd[:,idx],**kwargs)
    ax.set_xlabel("Frequency (Hz)")
    ax.set_ylabel("PSD (dB)")
    ax.set_title(f"Spectrum at time {time}s")
    return fig,ax

spectrogram(data, fs, time_step, window_length=None, window_shape='hamming', freq_range=None, detrend='constant', nfft=None, db_scale=True, p_ref=2e-05, boundary_pad=False)

Compute the ordinary (single-taper) PSD of the input data. This is similar to scipy.signal.spectrogram.

Parameters:

Name	Type	Description	Default
data	NDArray	(n_samples,) Input data	required
fs	float	Sampling frequency	required
time_step	float	Time step between frames in seconds	required
window_length	float	Window length in seconds. If None, will be set to the same as time_step. Defaults to None.	None
window_shape	Union[str, tuple]	The shape of the window function. Defaults to "hamming".	'hamming'
freq_range	Optional[list]	The desired frequency range. If None, will be set to [0, fs/2]. Defaults to None.	None
detrend	Literal['constant', 'linear', 'off']	Whether and how to detrend the signal. Defaults to "constant".	'constant'
nfft	Optional[int]	The number of FFT points. If None, will be set to the smallest power of 2 that is larger than the window length. Defaults to None.	None
db_scale	bool	Whether convert the result to db scale, i.e. 10log10(psd/p_ref**2). Defaults to True.	True
p_ref	float	If db_scale is True, the p_ref value used in the dB conversion. Defaults to 2e-5.	2e-05
boundary_pad	bool	Whether to pad the data with zeros at the beginning and end. This is useful when the data is not evenly divisible by the window length and time step. By default False. If True, the data will be padded with zeros at the beginning and end, so that the first frame is centered on the first sample of data, and all samples are included in (at least) one frame. If False, the first frame is centered at window_length/2 seconds after the first sample, and samples after n_frames*time_step+window_length seconds are ignored.	False

Returns:

Name	Type	Description
freqs	NDArray	(n_freqs,) Frequency points of the spectrogram
times	NDArray	(n_frames,) Time points of each frame
psd	NDArray	(n_freqs,n_frames) PSD spectrogram

Examples:

>>> times,freqs,psd = spectrogram(data,fs,time_step=0.001,window_length=0.005)

Source code in src/pymultitaper/spectral.py

def spectrogram(data:NDArray,fs:float,time_step:float,window_length:Optional[float]=None,window_shape:Union[str,tuple]="hamming",freq_range:Optional[list]=None,detrend:Literal["constant","linear","off"]="constant",nfft:Optional[int]=None,db_scale:bool=True,p_ref:float=2e-5,boundary_pad:bool=False)-> Tuple[NDArray,NDArray,NDArray]:
    """
    Compute the ordinary (single-taper) PSD of the input data. This is similar to `scipy.signal.spectrogram`.

    Args:
        data (NDArray): (n_samples,) Input data
        fs (float): Sampling frequency
        time_step (float): Time step between frames in seconds
        window_length (float, optional): Window length in seconds. If `None`, will be set to the same as `time_step`. Defaults to None.
        window_shape (Union[str,tuple], optional): The shape of the window function. Defaults to "hamming".
        freq_range (Optional[list], optional): The desired frequency range. If `None`, will be set to [0, fs/2]. Defaults to None.
        detrend (Literal["constant","linear","off"], optional): Whether and how to detrend the signal. Defaults to "constant".
        nfft (Optional[int], optional): The number of FFT points. If `None`, will be set to the smallest power of 2 that is larger than the window length. Defaults to None.
        db_scale (bool, optional): Whether convert the result to db scale, i.e. 10log10(psd/p_ref**2). Defaults to True.
        p_ref (float, optional): If `db_scale` is `True`, the `p_ref` value used in the dB conversion. Defaults to 2e-5.
        boundary_pad (bool, optional): Whether to pad the data with zeros at the beginning and end. This is useful when the data is not evenly divisible by the window length and time step. By default `False`.

            - If `True`, the data will be padded with zeros at the beginning and end, so that the first frame is centered on the first sample of data, and all samples are included in (at least) one frame.
            - If `False`, the first frame is centered at `window_length/2` seconds after the first sample, and samples after `n_frames*time_step+window_length` seconds are ignored.

    Returns:
        freqs (NDArray): (n_freqs,) Frequency points of the spectrogram
        times (NDArray): (n_frames,) Time points of each frame
        psd (NDArray): (n_freqs,n_frames) PSD spectrogram

    Examples:
        >>> times,freqs,psd = spectrogram(data,fs,time_step=0.001,window_length=0.005)
    """
    window_length = time_step if window_length is None else window_length
    n_winlen = int(window_length*fs)
    win,weights = _get_1d_window(window_shape,n_winlen)
    return _spectrogram(data=data,fs=fs,time_step=time_step,win=win,weights=weights,freq_range=freq_range,detrend=detrend,nfft=nfft,db_scale=db_scale,p_ref=p_ref,boundary_pad=boundary_pad)