R:优化尖峰修剪功能
R: Optimise spike pruning function
由于我还没有找到用于分析电生理数据的 R 包,所以我使用了我组的尖峰 p运行ing 函数:
prune.spikes <- function(spikes, min.isi) {
# copy spike matrix
prunedspikes <- spikes
# initialise index of last spike: infinitely before the first one.
for (i in 1:ncol(spikes)) {
last <- -Inf
for (j in 1:nrow(spikes)) {
if (spikes[j, i] == 1) {
if (j - last < min.isi) {
prunedspikes[j, i] <- 0; # remove the spike
}
else {
last <- j
}
}
}
}
return(prunedspikes)
}
该函数采用由 0 和 1 值组成的峰值向量或矩阵,如果出现在最小间隔内,则删除任何 1。
由于两个嵌套循环,它需要很长时间才能 运行。为了优化它,我想出了这个解决方案(删除一个循环):
prune.cols <- function(spikes, min.isi) {
prunedspikes <- apply(spikes, 2, FUN = prune.rows, min.isi = min.isi)
return(prunedspikes)
}
prune.rows <- function(spikes, min.isi) {
prunedspikes <- spikes
last <- -Inf
for (i in 1:length(spikes)) {
if (spikes[i] == 1) {
if (i - last < min.isi) {
prunedspikes[i] <- 0; # remove the spike
}
else {
last <- i
}
}
}
return(prunedspikes)
}
与原始版本相比,在大型数据集上调用 prune.cols
的速度明显加快(约 60 次)。不过,一个循环仍然存在。到目前为止,我还没有想出一个好的简单的解决方案。如何进一步完善功能?
如果速度差异还不是问题,可能保留循环而不是使用Rcpp会更好
根据 Hadley Wickham 的文章 Loops that should be left as is,使用此循环并不是一个坏主意,因为它可以归类为 递归关系 情况。
一旦速度成为瓶颈,那么求助于 Rcpp 或 this page(文章也建议)可能是解决方案。
像@Khashaa提出的那样,我借助Rcpp实现了功能:
NumericMatrix prunespikes(NumericMatrix spikes, double minisi) {
NumericMatrix prunedspikes = spikes;
int ncol = spikes.ncol();
int nrow = spikes.nrow();
for (int i = 0; i < ncol; i++) {
int last = 0;
while (spikes(last, i) == 0) {
last++;
}
for (int j = last + 1; j < nrow; j++) {
if (spikes(j, i) == 1) {
if (j - last < minisi) {
prunedspikes(j, i) = 0;
} else {
last = j;
}
}
}
}
return prunedspikes;
}
由于我还没有找到用于分析电生理数据的 R 包,所以我使用了我组的尖峰 p运行ing 函数:
prune.spikes <- function(spikes, min.isi) {
# copy spike matrix
prunedspikes <- spikes
# initialise index of last spike: infinitely before the first one.
for (i in 1:ncol(spikes)) {
last <- -Inf
for (j in 1:nrow(spikes)) {
if (spikes[j, i] == 1) {
if (j - last < min.isi) {
prunedspikes[j, i] <- 0; # remove the spike
}
else {
last <- j
}
}
}
}
return(prunedspikes)
}
该函数采用由 0 和 1 值组成的峰值向量或矩阵,如果出现在最小间隔内,则删除任何 1。 由于两个嵌套循环,它需要很长时间才能 运行。为了优化它,我想出了这个解决方案(删除一个循环):
prune.cols <- function(spikes, min.isi) {
prunedspikes <- apply(spikes, 2, FUN = prune.rows, min.isi = min.isi)
return(prunedspikes)
}
prune.rows <- function(spikes, min.isi) {
prunedspikes <- spikes
last <- -Inf
for (i in 1:length(spikes)) {
if (spikes[i] == 1) {
if (i - last < min.isi) {
prunedspikes[i] <- 0; # remove the spike
}
else {
last <- i
}
}
}
return(prunedspikes)
}
与原始版本相比,在大型数据集上调用 prune.cols
的速度明显加快(约 60 次)。不过,一个循环仍然存在。到目前为止,我还没有想出一个好的简单的解决方案。如何进一步完善功能?
如果速度差异还不是问题,可能保留循环而不是使用Rcpp会更好
根据 Hadley Wickham 的文章 Loops that should be left as is,使用此循环并不是一个坏主意,因为它可以归类为 递归关系 情况。
一旦速度成为瓶颈,那么求助于 Rcpp 或 this page(文章也建议)可能是解决方案。
像@Khashaa提出的那样,我借助Rcpp实现了功能:
NumericMatrix prunespikes(NumericMatrix spikes, double minisi) {
NumericMatrix prunedspikes = spikes;
int ncol = spikes.ncol();
int nrow = spikes.nrow();
for (int i = 0; i < ncol; i++) {
int last = 0;
while (spikes(last, i) == 0) {
last++;
}
for (int j = last + 1; j < nrow; j++) {
if (spikes(j, i) == 1) {
if (j - last < minisi) {
prunedspikes(j, i) = 0;
} else {
last = j;
}
}
}
}
return prunedspikes;
}