# Event Classification Processor, Version 1.4
#
# This classifier detects dominant artifacts in EEG in an EEG signal, such as
# epileptic seizures, hiss, intermittent hiss, grooming, eating, and transient
# jumps caused by imperfect electrodes and animal movement. Its detection metric
# is event_pwr, which is the power in a band event_lo to event_hi. We intend for
# this band to exclude transients caused by electrodes and include power from the
# events we are interested in. The processor includes a baseline power calibrator,
# which we should run on the EEG recording for a while before we go back and start
# storing characteristics to disk.
#
# The processor produces six metrics. The event power metric is a measure of the
# ratio of the event-band power to the baseline power. The transient power is the
# power below the event band, and this we use to detect problems with electrodes
# or that result from signal loss. The transient power metric is a measure of the
# ratio of transient power to baseline power. When the transient power exceeds the
# baseline power by more than a factor of ten or so, this is usually an indication
# of a problem with the electrode.
#
# The high-frequency metric is a measure of the ratio of the high-frequency power
# to the event-band power. To be clear: the hf_pwr metric is NOT the ratio of hf
# power to baseline power, but the ratio of hf power to event power. This
# processor is designed to detect dominant events. If there is high frequency in
# the signal, this must be the dominant event, so that its power is a significant
# portion of the total event-band power. We are not trying to study high-frequency
# power in isolation, as we might if we were looking for HFOs in EEG.
#
# The spikiness and asymmetry metrics are measurs of how much the signal power is
# concentrated in short spikes and in upward or downward spikes. The asymmetry
# metric is below one half for negative-going spikes.
#
# The intermittency metric is a measure of the power of the demodulated
# high-frequency signal. If we look at the amplitude envelope of the
# high-frequency signal, and take its fourier transform so as to filter the
# envelope, we obtain the intermittency signal. The power of this signal compared
# to the power of the original high-frequency signal is an indication of how
# intermittent the high-frequency power is in the interval.
#
# Set up the Neuroarchiver's fourier transform and playback interval. Remove or alter
# this code if you want to change the interval the processor looks at.
set config(window_fraction) "0.1"
set config(play_interval) "1.0"
set config(glitch_threshold) "200"
set config(match_limit) "0.1"
# We set up the Neuroclassifier.
set config(classifier_types) "\
Transient brown \
Eating brown \
Grooming yellow \
Hiss_Short blue \
Hiss_Spike red \
Hiss_Large darkgreen \
Hiss_Small lightgreen"
set config(classifier_metrics) "event_pwr \
transient_pwr \
hf_pwr \
spikiness \
asymmetry \
intermittency"
# We begin this channel's section of the characteristics line with the
# channel number.
append result "$info(channel_num) "
# Define the event band in Hertz.
set event_lo 4.0
set event_hi 160.0
# After rectifying the event-band signal, which amounts to a demodulation
# of its amplitude envelope, we take the fourier transform of the de-modulated
# signal and calculate the power of this transform in a range of frequencies.
# We define this range below, in Hertz.
set int_lo 4.0
set int_hi 16.0
# Define the limits for what we consider to be high-frequency power, or hiss.
set hf_lo 60.0
set hf_hi 160.0
# If we have too much signal loss, we ignore this interval so as to avoid
# false event detection due to bad messages and jumps in signal value after
# a period of loss.
set max_loss 20.0
# The transient band extends from just above zero frequency to just below
# the event band. When we see excessive transient power, we ignore the
# playback interval during event detection.
set transient_pwr [expr 0.001*[Neuroarchiver_band_power 0.1 [expr $event_lo - 0.1] 0 0]]
# We choose our event band to include all activity that might contribute
# to the events we are looking for, and to exclude both transient noise and
# very high frequency quantization noise. In order to calculate the spikiness
# of a signal, which we do below, we must look at the sample values of the
# signal itself rather than its spectrum. If we look at the unfiltered signal,
# we will find our calculation of spikiness is disturbed by transients. Here
# we instruct the band power routine to overwrite the info(values) string with
# the inverse transform of the event band spectrum.
set event_pwr [expr 0.001*[Neuroarchiver_band_power $event_lo $event_hi 0 1]]
# The spikiness of a signal is the ratio of its range to its standard deviation
# after filtering to the event frequency band. We divide by the the fourth root
# of the number of samples in order to overcome the effect of a normal distribution
# of samples.
scan [lwdaq ave_stdev $info(values)] %f%f%f%f ave stdev max min
if {$stdev > 0} {
set spikiness [expr ($max-$min)/($stdev)]
} {
set spikiness 1.0
}
# The asymmetry of a signal is a measure of how much it protrudes downwards
# or upwards from its average value. We construct a measure that is 1.0 when
# the signal is symmetric, less than 1.0 when it protrudes downwards, and greater
# than 1.0 for upwards.
set count 0
set upcount 0
set downcount 0
foreach v $info(values) {
if {($v - $ave) >= 2.0*$stdev} {incr upcount}
if {($ave - $v) >= 2.0*$stdev} {incr downcount}
if {abs($v - $ave) >= $stdev} {incr count}
}
if {$count>0} {
set asymmetry [expr 1.001+1.0*($upcount-$downcount)/$count]
} {
set asymmetry 1.0
}
# Obtain the high-frequency power. We choose this band to include the
# power we see in high-frequency bursts, but to exclued the harmonics of
# eating and other rhythms. In order to calculate the intermittency
# of a signal, we overwrite the info(values) string with the inverse
# transform of the high-frequency band spectrum.
set hf_pwr [expr 0.001*[Neuroarchiver_band_power $hf_lo $hf_hi 0 1]]
# The intermittency of a signal is the low-frequency power of the rectified
# high-frequency signal. Thus we take the high-frequency signal, rectify it,
# take the transform, and calculate the power in the low-frequency band.
set rectified_values ""
foreach v $info(values) {
append rectified_values "[expr abs($v)] "
}
set spectrum [lwdaq_fft $rectified_values]
set f 0
set intermittency_pwr 0.0
foreach {a p} $spectrum {
if {$f >= $int_lo} {
set intermittency_pwr [expr $intermittency_pwr + 0.001*$a*$a]
}
set f [expr $f + $info(f_step)]
if {$f > $int_hi} {break}
}
# Here we perform baseline power calibration, which works by finding the
# minimum event band power. We increase the baseline power by a small
# fraction every interval, in order to accommodate an increase in electrode
# sensitivity. We take care to avoid using corrupted intervals for our
# baseline. We won't use an interval with poor reception. We go to some
# effort to avoid being fooled by the glitch filter, which can happen
# under rare circumstances with bad messages. .
upvar #0 Neuroarchiver_info(bp_$info(channel_num)) baseline_pwr
set baseline_pwr [expr 1.0001*$baseline_pwr]
if {($event_pwr < $baseline_pwr) && ($event_pwr > 0.0) && ($info(loss) < $max_loss)} {
set info(values) [Neuroarchiver_values $info(signal)]
set info(spectrum) [lwdaq_fft $info(values) -glitch 0 \
-window [expr round([llength $info(values)] * $config(window_fraction))]]
set pwr [expr 0.001*[Neuroarchiver_band_power $event_lo $event_hi 0 0]]
if {$pwr < $baseline_pwr} {
set baseline_pwr $event_pwr
}
}
# We calculate the characteristics we use for determining similarity
# of events. If we have signal loss, we set these metrics to zero.
if {($info(loss)<$max_loss) && ($event_pwr>0.0)} {
set transient_m [expr 1.0 / (1.0 + pow(5.0/($transient_pwr/$baseline_pwr),0.75))]
set event_m [expr 1.0 / (1.0 + pow(5.0/($event_pwr/$baseline_pwr),1.0))]
set hf_m [expr 1.0 / (1.0 + pow(0.1/($hf_pwr/$event_pwr),1.0))]
set spikiness_m [expr 1.0 / (1.0 + pow(8.0/$spikiness,4))]
set asymmetry_m [expr 1.0 / (1.0 + pow(1.0/$asymmetry,2))]
set intermittency_m [expr 1.0 / (1.0 + pow(0.1/($intermittency_pwr/$hf_pwr),4.0))]
} {
foreach a "event_m transient_m hf_m spikiness_m asymmetry_m intermittency_m" {
set [set a] 0.0
}
}
# We add a single-letter code that indicates the type of the interval, so
# far as this processor can deduce it. We have L for loss, U for unclassified
# and N for no event.
if {($info(loss)>=$max_loss)} {
append result "L "
} elseif {$event_m>=0.5} {
append result "U "
} else {
append result "N "
}
# We add the baseline power to the characteristics. If we calculate
# the baseline power later, when looking through an archive, we could miss
# the defining low-power interval that sets the baseline.
append result "[format %.1f $baseline_pwr] "
# We record in the result line all the characteristics of the playback
# interval. The first one we record is the one we will use as a threshold
# for the occurance of an event, this threshold being 0.5 by default.
foreach a "event_m transient_m hf_m spikiness_m asymmetry_m intermittency_m" {
append result "[format %.3f [set $a]] "
}