微信公众号搜"智元新知"关注
微信扫一扫可直接关注哦!

如何计算R中给定范围内高于某个值的实例数?

如何解决如何计算R中给定范围内高于某个值的实例数?

我有一个相当大的数据集,查看整个基因组中的 SNP。我正在尝试生成一个热图,该热图基于在整个基因组中 x 碱基对的滑动窗口内有多少 SNP 具有超过 50 的 BF(贝叶斯因子)值。例如,在前 1,000,000 个碱基对中可能有 5 个感兴趣的 SNP,然后在接下来的 1,000 个碱基对中有 3 个,依此类推,直到我到达基因组的末端,这将用于生成单行热图。目前,我的数据如下所示:

SNP BF  BP
0001_107388 11.62814713 107388
0001_193069 2.333472447 193069
0001_278038 51.34452334 278038
0001_328786 5.321968927 328786
0001_523879 50.03245434 523879
0001_804477 -0.51777189 804477
0001_990357 6.235452787 990357
0001_1033297    3.08206707  1033297
0001_1167609    -2.427835577    1167609
0001_1222410    52.96447989 1222410
0001_1490205    10.98099565 1490205
0001_1689133    3.75363951  1689133
0001_1746080    3.519987207 1746080
0001_1746450    -2.86666016 1746450
0001_1777011    0.166999413 1777011
0001_2114817    3.266942137 2114817
0001_2232084    50.43561123 2232084
0001_2332903    -0.15022324 2332903
0001_2347062    -1.209000033    2347062
0001_2426273    1.230915683 2426273

其中 SNP = SNP ID,BF = 贝叶斯因子,而 BP = 基因组上的位置(我在其中捏造了几个 > 50 个值,以便数据适合此示例)。>

问题是我没有每个基因组位置的 SNP,否则我可以简单地根据行数拆分感兴趣的窗口,然后计算 BF 列中的许多行超过 50。有什么办法我可以计算基因组位置的不同窗口内感兴趣的 SNP 的数量吗?最好使用 R,但如果能完成工作,则使用其他语言(如 Python 或 Bash)没有问题。

谢谢!

解决方法

library(slider); library(dplyr)
my_data %>%
  mutate(count = slide_index(BF,BP,~sum(.x > 50),.before = 999999))

这里统计BP最后1M的窗口中有多少BF>50。

            SNP         BF      BP count
1   0001_107388 11.6281471  107388     0
2   0001_193069  2.3334724  193069     0
3   0001_278038 51.3445233  278038     1
4   0001_328786  5.3219689  328786     1
5   0001_523879 50.0324543  523879     2
6   0001_804477 -0.5177719  804477     2
7   0001_990357  6.2354528  990357     2
8  0001_1033297  3.0820671 1033297     2
9  0001_1167609 -2.4278356 1167609     2
10 0001_1222410 52.9644799 1222410     3
11 0001_1490205 10.9809957 1490205     2
12 0001_1689133  3.7536395 1689133     1
13 0001_1746080  3.5199872 1746080     1
14 0001_1746450 -2.8666602 1746450     1
15 0001_1777011  0.1669994 1777011     1
16 0001_2114817  3.2669421 2114817     1
17 0001_2232084 50.4356112 2232084     1
18 0001_2332903 -0.1502232 2332903     1
19 0001_2347062 -1.2090000 2347062     1
20 0001_2426273  1.2309157 2426273     1

版权声明:本文内容由互联网用户自发贡献,该文观点与技术仅代表作者本人。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容, 请发送邮件至 dio@foxmail.com 举报,一经查实,本站将立刻删除。

相关推荐


Selenium Web驱动程序和Java。元素在(x,y)点处不可单击。其他元素将获得点击?
Python-如何使用点“。” 访问字典成员?
Java 字符串是不可变的。到底是什么意思?
Java中的“ final”关键字如何工作?(我仍然可以修改对象。)
“loop:”在Java代码中。这是什么,为什么要编译?
java.lang.ClassNotFoundException:sun.jdbc.odbc.JdbcOdbcDriver发生异常。为什么?
这是用Java进行XML解析的最佳库。
Java的PriorityQueue的内置迭代器不会以任何特定顺序遍历数据结构。为什么?
如何在Java中聆听按键时移动图像。
Java“Program to an interface”。这是什么意思?