4

「a」と「b」という 2 つの条件があるとします。ニューロンは、条件「a」では平均 40 スパイク/秒 (Hz)、条件「b」では 80 スパイク/秒で発火します。条件「a」への応答は 20 回表示され、条件「b」は 10 回表示され、各表示は 1000 ミリ秒です。

AB <- rbind(
    ldply( 1:20, 
        function(trial) { 
          data.frame( 
              trial=trial, 
              cond=factor('a',c('a','b')), 
              spiketime = runif(40,0,1000))
        }
    ), ldply(21:30, 
        function(trial) {
          data.frame(
              trial=trial, 
              cond=factor('b',c('a','b')), 
              spiketime = runif(80,0,1000))
        }
  )
)

単純なヒストグラムは、次のようにプロットできます。

qplot(spiketime, data=AB, geom='line',stat='bin',y=..count.., 
      xlim=c(0,1000), colour=cond, binwidth=100,xlab='Milliseconds')

ただし、これはプレゼンテーション全体の平均ではないため、y 軸の値はほぼ同じです。y 軸に沿ってスパイク率 (スパイク/秒) をプロットしたいと思います。これは、条件「b」が 1 秒あたり約 2 倍のスパイクを誘発することを示しています。プレゼンテーションの数が増えてもスパイク率は増加せず、ノイズが少なくなるだけです。データフレーム AB を前処理せずにこれを行う方法はありますか?

言い換えれば、次の行に沿って何かを行うことができますか?

qplot(spiketime, data=AB, geom='line',stat='bin',
      y=..count../num_presentations*1000/binwidth, ylab='Spikes per second',
      xlim=c(0,1000), colour=cond, binwidth=100,xlab='Milliseconds')

ここで、num_presentations は、条件「a」の場合は 20、条件「b」の場合は 10 であり、1000/binwidth は、単位を正しく取得するための定数にすぎませんか?

4

2 に答える 2

2

条件の平均ではありません。それはそれらを合計します。条件 a には 20x40 = 800 ポイントがあり、条件 b には 10*80 = 800 ポイントがあるため、これらの「ヒストグラム」の下の「面積」は同じになります。条件内の各試行の重みを等しくする必要がありますが、各ポイントの重みを等しくする必要はありません。これは、前処理ステップとして実行する必要があります。

trial.group <- unique(AB[,c("trial","cond")])
hists <- dlply(AB, .(trial), function(x) {hist(x$spiketime, breaks=10, plot=FALSE)})
hists.avg <- ddply(trial.group, .(cond), function(x) {
  hist.group <- ldply(hists[x$trial], function(y) {
    data.frame(mids=y$mids, counts=y$counts)
  })
  ddply(hist.group, .(mids), summarise, counts=mean(counts))
})

ggplot(data=hists.avg, aes(x=mids, y=counts, colour=cond)) + geom_line()

これはhist、各トライアルのデータを個別にビン化するために使用し、トライアル グループのカウントを平均化します。これにより、各条件に等しい重みが与えられ、各条件内の各試行に等しい重みが与えられます。

編集1：

@kohske ソリューションを使用しますが、試行回数を明示的に入力する代わりに計算します。

tmp <- as.data.frame(table(unique(AB[,c("trial","cond")])["cond"]))
names(tmp) <- c("cond","ntrial")
AB <- merge(AB, tmp)

ggplot(AB, aes(spiketime, ntrial=ntrial, colour=cond)) + 
  stat_bin(aes(y=..count../ntrial*1000/100), binwidth=100, geom="line", position="identity") +
  xlim(0,1000) + 
  labs(x='Milliseconds', y="Firing rate [times/s]")
于 2011-08-15T21:55:56.860 に答える
2

ここに解決策があります：

AB$ntrial <- ifelse(AB$cond=="a", 20, 10)
ggplot(AB, aes(spiketime, ntrial=ntrial, colour=cond)) + 
  stat_bin(aes(y=..count../ntrial*1000/100), binwidth=100, geom="line", position="identity") +
  xlim(0,1000) + 
  labs(x='Milliseconds', y="Firing rate [times/s]")

于 2011-08-16T00:00:38.463 に答える