Power Spectra

[Kevan Hashemi]

A customer writes:

I’m looking to extract the ECoG power in my GBM-epilepsy mice (single lead cortical electrode/transmitter) for different lengths of time e.g. 24 hours, 1 hour etc. Specifically power in the typical frequency bands (delta, theta, alpha, beta, gamma (<80Hz)) for the epochs. I know you had a .tcl script on your website at one point for power, do you have an updated version of this or any other way of exporting this kind of info within your software? I’m also looking to do the same kind of thing for coastline if that’s something that’s possible! Let me know what you think.

[Kevan Hashemi]

We generate long-term power spectra in two stages. We do address these stages in our Event Detection page, so each stage has a link to a section of our Event Detection page.

Stage One: Calculate Interval Spectra

We use interval processing in the Neuroplayer to calculate the power spectrum of each interval in a selection of frequency bands that we specify. The processing produces a characteristics file. Each line in the characteristics file provides the spectrum of all the channels we selected for processing. The interval can be anything from 1/16 s to 32 s. For twenty-four hour spectrogram, we might as well choose 8-s intervals. They let us see everything down to 1 Hz clearly, and they provide ample resolution.

To perform "interval processing" you need an "interval processor", which is a text file containing TclTk code that the Neuroplayer will run. Here is an example processor script, including comments. It produces power in four bands, and appends the coastline as well.

# Obtain the power spectrum in bands 1-4, 4-12, 12-30, 30-80 Hz, where each band
# includes the lower boundary but not the upper boundary. Units are ksqcnt.
# Follow these band powers with the signal coastline in units of kcnt. Set all
# values to zero if signal loss is greater than 20%.
set max_loss 20.0
append result "$info(channel_num) "
set f_lo 1.0
foreach f_hi {4.0 12.0 30.0 80.0} {
if {$info(loss) <= $max_loss} {
set bp [Neuroplayer_band_power [expr $f_lo + 0.01] $f_hi 0]
} else {
set bp 0.0
}
append result "[format %.2f [expr 0.001 * $bp]] "
set f_lo $f_hi
}
if {$info(loss) <= $max_loss} {
set cl [lwdaq coastline_x $info(values)]
} else {
set cl "0.0"
}
append result "[format %.1f [expr 0.001 * $cl]] "

I save this file to PB.tcl, a text file on my hard driver. I select it as my processor in the Neuroplayer. I check Enable and Save. I set the interval to 8 s, I enter "88 108" into the "select" string. I open the sample LWDAQ/Images/M1631648948.ndf archive and press Play. I get a file M1631648948_PB.txt that contains the following lines.

M1631648948.ndf 0.0 88 72.9 10.52 6.43 2.00 1.52 108 68.2 1.84 4.55 1.84 1.20
M1631648948.ndf 8.0 88 70.2 9.76 8.52 1.96 1.23 108 68.0 2.80 5.90 1.84 1.23
M1631648948.ndf 16.0 88 76.2 2.93 5.98 2.59 1.55 108 63.1 3.32 9.04 1.75 0.94
M1631648948.ndf 24.0 88 75.8 11.77 8.83 2.47 1.38 108 63.4 2.80 6.92 1.90 1.00
M1631648948.ndf 32.0 88 74.7 5.21 5.71 2.04 1.58 108 69.8 4.36 5.72 1.49 1.54
M1631648948.ndf 40.0 88 73.5 4.12 5.85 2.42 1.34 108 67.5 1.28 4.61 1.95 1.31
M1631648948.ndf 48.0 88 78.4 5.70 8.36 2.52 1.75 108 70.8 5.09 10.75 1.83 1.27

Here we see the archive name, archive time, channel number 88, band powers, channel number 108, band powers. The band power are in "kilo square counts", or "ksqcnt". If the dynamic range of your transmitter is 27 mV, one count is 0.4 µV, so 1 ksqcnt is 160 µV^2. Coastline is in kcnt, where 1 kcnt is 400 µV.

Stage Two: Calculate Average Power

Once we have power from 8-s intervals, we can obtain the average power for each band over longer intervals with a characteristics averaging script such as PBA. We edit the PBA script to select one or more channels that are present in each characteristics line, specify the averaging interval, and tell it the number of "bands" that are present for each channel. In the above case we have four frequency bands, but also a coastline value, so we would set num_bands to 5. When we are satisfied, we save our version of PBA, load it into the LWDAQ Toolmaker, and Execute. We will be asked to select all the characteristics files. The script will then produce our output file with average powers.

Now that the average band powers are in a text file, it's up to you to plot them in the way you like.

Best, Kevan