#' Compute bands powers from signal.
#' Divides signal into frequency bands and performs Fast Fourier Transforms (FFT) to compute bands powers.
#' @param x Vector containing signal.
#' @param Fs Signal sample rate.
#' @param bands List of frequency bands with lower and upper limits. Defaults to \code{list(delta = c(0.5,3.5), theta = c(3.5,8), alpha = c(8,12), beta = c(12,30), gamma1 = c(30,40))}
#' @param normalize If set
#' @param butter butter order or false.
#' @return A named vector containing bands powers.
#' @export
bands_power <- function(x, Fs,
bands = list(delta = c(0.5,3.5),
theta = c(3.5,8),
alpha = c(8,12),
beta = c(12,30),
gamma1 = c(30,40)),
butter = FALSE,
normalize = NA){
#x[is.na(x)] <- 0
if(butter != FALSE){
filt <- signal::butter(butter, 0.1)
x <- signal::filtfilt(filt, x)
}
powers <- c()
for(band_label in names(bands)){
power <- sum(eegkit::eegfft(x,Fs,bands[[band_label]][1],bands[[band_label]][2])$strength)
if(length(normalize) == 2){
power <- power/sum(eegkit::eegfft(x,Fs,normalize[1],normalize[2])$strength)
}
powers <- c(powers,power)
}
names(powers) <- names(bands)
powers
}
#' Compute bands powers from a full hypnogram.
#'
#' @param record record.
#' @param channel channel to split and hyp.
#' @param bands bands. Default "delta","theta","alpha","beta", "gamma1", "denominator", "broadband".
#' @param normalize normalize.
#' @param labels labels.
#' @param butter butter order
#' @param resampling sRate to resample by.
#' @return a df.
#' @export
hypnogram_band_powers <- function(record,
channel,
bands=list(delta = c(0.5,3.5),
theta = c(3.5,8),
alpha = c(8,12),
beta = c(12,30),
gamma1 = c(30,40)),
normalize = c(4,40),
labels = c("N3","N2","N1","REM","AWA"),
butter = FALSE,
resampling = FALSE){
sRate <- record[["channels"]][[channel]][["metadata"]][["sRate"]]
hypnogram <- hypnogram(record[["events"]],labels)
signal <- record[["channels"]][[channel]][["signal"]]
if(resampling != FALSE){
signal <- signal::resample(signal, 1,(sRate/resampling) )
sRate <- resampling
}
signal <- sleepr::split_signal(signal = signal,
hypnogram = hypnogram,
sRate = sRate)
pw <- dplyr::bind_rows(lapply(signal,function(x){
as.list(sleepr::bands_power(x = x, Fs = sRate, bands = bands, normalize = normalize, butter = butter))
}))
pw$stage <- hypnogram$event
pw$epoch <- c(1:nrow(pw))
return(pw)
}
#' Compute EEMD mean,var,skew,kur
#'
#' @param record record.
#' @param channel channel to split and hyp.
#' @param num_imfs numinfs
#' @param labels lablssl
#' @param butter butter
#' @return a df.
#' @export
hypnogram_eemd <- function(record,
channel,
num_imfs = 8,
labels = c("N3","N2","N1","REM","AWA"),
butter = FALSE){
# Rlibeemd::eemd(c3a2.period$signal,
# num_imfs = 8)
sRate <- record[["channels"]][[channel]][["metadata"]][["sRate"]]
hypnogram <- hypnogram(record[["events"]],labels)
signal <- sleepr::split_signal(signal = record[["channels"]][[channel]][["signal"]],
hypnogram = hypnogram,
sRate = sRate)
eemddf <- data.frame()
for(i in c(1:length(signal))){
x <- signal[[i]]
if(butter != FALSE){
filt <- signal::butter(butter, 0.1)
x <- signal::filtfilt(filt, x)
}
eemd <- Rlibeemd::eemd(x,num_imfs = num_imfs)
epdf <- data.frame(epoch = i)
for(j in c(1:num_imfs)){
epdft <- data.frame(x1 = mean(eemd[,j]),
x2 = stats::var(eemd[,j]),
x3 = e1071::skewness(eemd[,j]),
x4 = e1071::kurtosis(eemd[,j]))
colnames(epdft) <- c(paste0("eemd_",j,"_mean"),
paste0("eemd_",j,"_var"),
paste0("eemd_",j,"_skew"),
paste0("eemd_",j,"_kur"))
epdf <- cbind(epdf,epdft)
}
eemddf <- dplyr::bind_rows(eemddf,epdf)
}
return(eemddf)
}
#' aggregate_band_powers
#'
#' @param hypnogram_band_powers hypnogram_band_powers.
#' @return Wide dataframe.
#' @export
aggregate_band_powers <- function(hypnogram_band_powers){
colnames <- colnames(hypnogram_band_powers)[!(colnames(hypnogram_band_powers) %in% c("denominator","broadband","epoch"))]
bandnames <- colnames[!(colnames %in% c("stage"))]
means <- stats::aggregate(hypnogram_band_powers[,bandnames],by=list(hypnogram_band_powers$stage),FUN=mean,na.rm=TRUE)
hypnogram_band_powers_nrem <- hypnogram_band_powers[hypnogram_band_powers$stage %in% c("N1","N2","N3"),]
if (nrow(hypnogram_band_powers_nrem) > 0){
hypnogram_band_powers_nrem$stage <- "NREM"
meansnrem <- stats::aggregate(hypnogram_band_powers_nrem[,bandnames],by=list(hypnogram_band_powers_nrem$stage),FUN=mean,na.rm=TRUE)
means <- rbind(means,meansnrem)
}
means$id <- 1
means <- stats::reshape(means, idvar = "id", timevar = "Group.1", direction = "wide")
means$id <- NULL
colnames(means) <- tolower(gsub("\\.", "_", colnames(means)))
return(means)
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.