r - 计算 R 中的单例数

Question

我有一些 RNA-seq 数据，我需要计算单例的数量。我们将单例定义为在附近(距离任一侧 100 个碱基的距离内)没有任何其他读取的读取。

我有一个数据框，其中包含每次读取的开始坐标和结束坐标。我正在使用 R 来执行此操作。

我暂时写了这段代码，但应用不正确，因此出现错误。

begin_end <- data.frame(begin_coordinate, final_coordinate)
apply(begin_end, 1, function(x) x[,1]-(x-1)[,2])

数据框的第一行是:

> head(begin_end)

    begin   final
1   60507   60551
2   60790   60840
3   62004   62051
4   62819   62868
5   65141   65187

第一个似乎是单例，因为下一个读取在结束后开始超过 100 个碱基，数据集第一行中的其余部分也是如此。但是数据框很长，我希望不是所有的读取都是单例的。

Answer 1

这与@jeremycg 对 dplyr 的 lag 和 lead 所做的相同，但在 data.table 中:

library(data.table)
setDT(begin_end)

begin_end[{
  d = begin - shift(final, type="lag")
  pmin(d, shift(d, type="lead"), na.rm=TRUE) > 100
}]

---

评论。基本的data.table语法是DT[i,j]。 i 用于过滤输入，而 j 用于修改输出。

我们在上面使用了 i，但是为了检查它是如何工作的，我们可以将相关向量折腾到 j 中:

begin_end[,{
  d       = begin - shift(final, type="lag")
  d_lead  = shift(d, type="lead")
  my_pmin = pmin(d, d_lead, na.rm=TRUE)
  c(.SD, list(d = d, d_lead = d_lead, my_pmin = my_pmin))
}]

#    begin final    d d_lead my_pmin
# 1: 60507 60551   NA    239     239
# 2: 60790 60840  239   1164     239
# 3: 62004 62051 1164    768     768
# 4: 62819 62868  768   2273     768
# 5: 65141 65187 2273     NA    2273

.SD 是表中已有的列向量列表，是 Subset of Data 的缩写。