BNBC 1.18.0
这BNBC软件包提供了执行归一化和批化校正的功能跨样品关于从HI-C获得的数据(Lieberman-Aiden等,2009)实验。
在此软件包中,我们实施了用于从HI-C接触矩阵集的一般子集和数据提取的工具,以及接触矩阵的平滑,跨样本归一化和跨样本批处理效果校正方法。
BNBC期望输入
格兰格
代表基因组测定的对象,单个范围的宽度等于用于划分基因组的bin尺寸。列表
(平方)接触矩阵。数据框架
包含样本级协变量的对象(即性别,处理日期等)。如果您使用此包,请引用(Fletez-Brant等人,2017年)。
非常理解的是,HI-C接触矩阵在观察到的接触数量中表现出指数衰减,这是一对相互作用基因座之间的距离的函数。在这项工作中,我们像最近一样运作(即(Yang等人,2017年)),在特定距离的所有基因座相互作用的集合中,一次一个染色体。对于给定的距离\(k \),相关的基因座集合在每个联系矩阵中列出,包括条目\(k \)- 矩阵对角线(主对角线称为第一个对角线)。我们将这些对角线称为矩阵“频段”。
该文档具有以下依赖关系
图书馆(BNBC)
BNBC使用联络组
代表给定基因组LOCIS相互作用集的一组接触矩阵的类。该课程有3个插槽:
Rowdata
: 一个格兰格
对象为分区基因组的每个bin具有1个元素。酷塔
: 一个数据框架
每个样本中包含信息的对象(即性别)。联系人
: 一个列表
接触矩阵。我们希望每个人联络组
表示来自1个染色体的数据。我们正在考虑一个全基因组容器。包装提供了一个示例数据集。
数据(CGEX)CGEX
##类的对象`contactGroup`代表与## 1282垃圾箱的联系矩阵## 40 kb平均宽度平均宽度## 14样本
创建一个联络组
对象需要指定上面的3个插槽:
CGEX < - ContactGroup(Rowdata = Rowdata(CGEX),Contacts = Contacts(CGEX),Coldata = Coldata(CGEX))
请注意,在此示例中,我们为每个插槽使用了访问者方法。也有相应的“设置器”方法,例如rowdata(cgex)< -
。
打印联络组
对象给出由Rowdata
插槽,基因组距离(即100KB)和样品数量的宽度宽度:
CGEX
##类的对象`contactGroup`代表与## 1282垃圾箱的联系矩阵## 40 kb平均宽度平均宽度## 14样本
将数据输入BNBC您需要一个接触矩阵的列表,每个样本一个。我们假设接触矩阵是正方形的,没有缺少值。我们不要求通过各种偏置校正软件对数据进行了转换或预处理,并为某些简单的预处理技术提供了方法。
当前没有标准的HI-C数据格式。相反,不同的组通常以文本文件的形式产生自定义格式。因为接触矩阵是正方形的,所以通常只分配上部或下三角矩阵。在这种情况下,您可以使用以下技巧制作矩阵正方形:
垫子<-matrix(1:9,nrow = 3,ncol = 3)垫[lower.tri(mat)] <-0 mat
## [,1] [,2] [,3] ## [1,] 1 4 7 ## [2,] 0 5 8 ## [3,] 0 0 9
##现在,我们用上三角垫填充下部三角矩阵[power.tri(mat)] <-mat [upper.tri(mat)]垫子
## [,1] [,2] [,3] ## [1,] 1 4 7 ## [2,] 4 5 8 ## [3,] 7 8 9
下面,我们证明将一组假设的接触矩阵转换为联络组
目的。物体upper.mats.list
应该是接触矩阵的列表,每个矩阵表示为上三角矩阵。我们也假设定位
成为一个基因组机
包含接触矩阵的基因座的对象,样本数据
成为一个数据框架
每样本信息(即单元类型,示例名称等)的信息。我们首先将所有接触矩阵转换为对称矩阵,然后使用构造方法方法ContactGroup()
创建对象。
##示例不运行##将上三角转换为对称矩阵垫<-lapply(upper.mats.list,function(m){m [lower.tri(m)] <-m [upper.tri(m)})##使用ContactGroup构造方法CG < - ContactGroup(rowdata = locidata,contacts = matslist,coldata = sampledata)
为此,联系人
列表必须具有与Rownames相同的名称酷塔
。
*.Cooler
文件这.Cooler
文件格式被广泛采用和支持BNBC
。我们假设一个简单的冷却器文件格式(请参阅?getchrcgfromcools
完整描述;重要的是,我们假设在所有样本中都观察到相同的交互冷却器
程序。我们的进入点是目录哪些交互存储在冷却器文件中。我们通过使用函数生成文件位置的索引来做到这一点getGenomeIdx()
。
coolerdir <-system.file(“ cooler”,package =“ bnbc”)cools cools <-list.files(coolerdir,pattern =“ cool $”,full.names = true)step <-4e4 bin.ixns.lists.list <--bnbc ::: getGenomeIdx(冷却[1],步骤)bin.ixns <-bin.ixns.list $ bin.ixns
作为输出,我们有一个列表bin.ixns.list
,两个要素。第一个是元素bin.ixns
,这是一个Data.Table
对象列出了在冷却器文件中观察到的一组交互。另一个元素是data.frame
目的垃圾箱
,其中列出了一组基因组垃圾箱及其坐标。这些是某些功能中使用的便利性,通常对最终用户不感兴趣。
有了我们的索引,我们可以一次将数据加载到内存中,这是一个染色体的数据(此时我们的方法无法处理 - 互相)。我们强调的是,通过在内存中一个染色体上的基因座之间的相互作用中进行所有观察,我们的算法非常有效,具有用于矩阵更新的自定义例程,并且只需要传递数据。
data(cgex)cool.cg <-bnbc ::: getchrcgfromcools(bin.ixns,files,files = cools = cool =“ chr22”,coldata = coldata = coldata(cgex)[1:2,])[[1]],联系人(cool.cg)[[1]))
## [1]
在此示例中,我们加载联络组
目的CGEX
在记忆中与IT的表示相比凉爽的
由冷却器
程序。然后我们使用该方法getchrcgfromcools()
将整个染色体的相互作用矩阵(在所有受试者上观察)加载到内存中。此时,用户有效联络组
对象,并且可以按照后续部分中所述进行分析。
我们提供了操纵单个矩阵频段以进行接触矩阵的设置器和GETER方法。首先,我们具有与单个矩阵谱带(无BNBC相关)的功能:
MAT.1 < - 联系人(CGEX)[[1]] MAT.1 [1000:1005,1000:1005]
## [,1] [,2] [,3] [,4] [,5] [,6] ## [1,] 5 30 36 32 53 20 ## [2,] 30 7 36 38 50 16## [3,] 36 36 3 44 51 15 ## [4,] 32 38 44 4 89 39 ## [5,] 53 50 51 89 36 55 ## [6,] 20 16 16 15 39 55 5
b1 < - band(mat = mat.1,band.no = 2)band(mat = mat.1,band.no = 2)<-b1 + 1 Mat.1 [1000:1005,1000:1005]
## [,1] [,2] [,3] [,4] [,5] [,6] ## [1,] 5 31 31 36 32 53 20 ## [2,] 31 7 37 38 38 50 16## [3,] 36 37 3 45 51 15 ## [4,] 32 38 45 4 90 39 ## [5,] 53 50 51 90 36 56 ## [6,] 20 16 16 15 39 56 5
在此示例中,接触矩阵的主对角也是上面印刷示例的主要对角线。同样,第二个乐队也是第一个偏离对角线,也是印刷示例的第一个非对角线。从打印的示例可以看出,更新矩阵频段是对称操作,并更新了矩阵上部和下三角形中的第一个偏外式操作。
为了在联系矩阵列表中利用它,我们有cgapply()
将相同函数应用于每个矩阵的函数。它支持使用平行。
为了调整深度测序的差异,我们首先应用logcpm
转换(Law等,2014)每个触点矩阵。这种转换将每个触点矩阵除以该矩阵的上三角的总和(在每个矩阵中添加0.5,将1个添加到上三角形的总和中),将所得矩阵缩放\(10^6 \)最后,获取缩放矩阵的日志(在接受对数之前,添加了一个福吉因子和分母)。
cgex.cpm <-logcpm(cgex)
此外,我们使用平滑的内核平滑每个接触矩阵,以减少垃圾箱宽度选择的伪像。我们支持Box和Gaussian Smoother。
cgex.smooth < - boxsmoother(cgex.cpm,h = 5)##或## cgex.smooth <-gausssmoother(cgex.cpm,radius = 3,sigma = 4)
BNBC分别在每个矩阵频段上操作。对于每个矩阵频段\(k \),我们提取每个样品在该频段上的观察结果并形成一个矩阵\(M \)从那些乐队如果乐队\(k \)有\(d \)条目,然后之后logcpm
转型,\(m \ in \ mathbb {r}^{n \ times d} \)。对于每个此类矩阵,我们首先应用分位数归一化(Bolstad等,2003)纠正分配差异,然后战斗(Johnson,Li和Rabinovic 2007)纠正批处理效果。
在这里我们将使用BNBC要在前10个矩阵频段上进行批处理校正,从第二个矩阵频段开始,以第十一个结束。
cgex.bnbc <-bnbc(cgex.smooth,batch = coldata(cgex.smooth)$ batch,threshold = 1e7,step = 4e4,nbands = 11,verbose = false)
##在单个批处理中发现了70个具有均匀表达的基因(所有零);这些不会调整以批量调整。##在单个批处理中发现了95个具有均匀表达的基因(所有零);这些不会调整以批量调整。##在单个批次(所有零)中找到了103个具有均匀表达的基因;这些不会调整以批量调整。##在单个批处理中发现了77个具有均匀表达的基因(所有零);这些不会调整以批量调整。##在单个批处理中找到了63个具有均匀表达的基因(所有零);这些不会调整以批量调整。 ## Found 76 genes with uniform expression within a single batch (all zeros); these will not be adjusted for batch. ## Found 66 genes with uniform expression within a single batch (all zeros); these will not be adjusted for batch. ## Found 82 genes with uniform expression within a single batch (all zeros); these will not be adjusted for batch. ## Found 79 genes with uniform expression within a single batch (all zeros); these will not be adjusted for batch. ## Found 77 genes with uniform expression within a single batch (all zeros); these will not be adjusted for batch.
## R版本4.2.0 RC(2022-04-19 R82224)##平台:x86_64-pc-linux-gnu(64位)### blas:/home/biocbuild/bbs-3.15-bioc/r/lib/libblas.so ## lapack:/home/biocbuild/bbs-3.15-bioc/rib/lib/libb/librlapack.so ### ## ## locale:## [1] lc_ctype = en_us.utf-8 lc_numeric = c ## [3] lc_time = en_gb lc_collate = c ## [5] lc_us.us.utf-8 lc_messages = en_us.utf-8 ## [7]lc_paper = en_us.utf-8 lc_name = c ## [9] lc_address = c lc_telephone = c ## [11] lc_measurement = en_us.utf-8 lc_istientification = c ## ## ## ## ##附件:统计图形GRDEVICES UTILS方法## [8]基础## ##其他附件:## [1] BNBC_1.18.0 SummarizedExperiment_1.26.0 ## [3] BioBase_2.56.0 Genomicranges_1.1.1.48.48.48.0 ## [5] GenomeInomeInfodb_32。0iranges_2.30.0 ## [7] s4Vectors_0.34.0 matrixgenerics_1.8.0 ## [9] matrixstats_0.62.0 biocgenerics_0.42.0.42.0 ## [11] Biocstyle_2.24.01] HTTR_1.4.4 SASS_0.4。1edgeR_3.38.0 ## [4] bit64_4.0.5 jsonlite_1.8.0 splines_4.2.0 ## [7] bslib_0.3.1 BiocManager_1.30.17 blob_1.2.3 ## [10] GenomeInfoDbData_1.2.8 tiff_0.1-11 yaml_2.3.5 ## [13] RSQLite_2.2.12 lattice_0.20-45 limma_3.52.0 ## [16] digest_0.6.29 XVector_0.36.0 htmltools_0.5.2 ## [19] preprocessCore_1.58.0 Matrix_1.4-1 XML_3.99-0.9 ## [22] genefilter_1.78.0 bookdown_0.26 zlibbioc_1.42.0 ## [25] fftwtools_0.9-11 xtable_1.8-4 jpeg_0.1-9 ## [28] BiocParallel_1.30.0 annotate_1.74.0 KEGGREST_1.36.0 ## [31] mgcv_1.8-40 EBImage_4.38.0 cachem_1.0.6 ## [34] cli_3.3.0 survival_3.3-1 magrittr_2.0.3 ## [37] crayon_1.5.1 memoise_2.0.1 evaluate_0.15 ## [40] nlme_3.1-157 data.table_1.14.2 tools_4.2.0 ## [43] stringr_1.4.0 Rhdf5lib_1.18.0 locfit_1.5-9.5 ## [46] DelayedArray_0.22.0 AnnotationDbi_1.58.0 Biostrings_2.64.0 ## [49] compiler_4.2.0 jquerylib_0.1.4 rlang_1.0.2 ## [52] rhdf5_2.40.0 grid_4.2.0 RCurl_1.98-1.6 ## [55] rhdf5filters_1.8.0 htmlwidgets_1.5.4 bitops_1.0-7 ## [58] rmarkdown_2.14 abind_1.4-5 DBI_1.1.2 ## [61] R6_2.5.1 knitr_1.38 fastmap_1.1.0 ## [64] bit_4.0.4 stringi_1.7.6 parallel_4.2.0 ## [67] sva_3.44.0 Rcpp_1.0.8.3 vctrs_0.4.1 ## [70] png_0.1-7 xfun_0.30
Bolstad,B M,R A Irizarry,M Astrand和T P Speed。2003年。“基于方差和偏置的高密度寡核苷酸阵列数据的归一化方法的比较。”生物信息学19:185–93。https://doi.org/10.1093/bioinformatics/19.2.185。
Fletez-Brant,Kipper,Yunjiang Qiu,David U Gorkin,Ming Hu和Kasper D Hansen。2017年。“在HI-C实验中消除样品之间不需要的变化。”生物xiv,214361。https://doi.org/10.1101/214361。
约翰逊,W埃文,成李和阿里尔·拉比诺维奇。2007年。“使用经验贝叶斯方法调整微阵列表达数据中的批处理效应” 8:118-27。https://doi.org/10.1093/biostatistics/kxj037。
Law,慈善机构W,Yunshun Chen,Wei Shi和Gordon K Smyth。2014年。“ VOOM:精确权重解锁RNA-Seq读取计数的线性模型分析工具。”基因组生物学15:R29。https://doi.org/10.1186/GB-2014-15-2-R29。
Lieberman-Aiden,Erez,Nynke L Van Berkum,Louise Williams,Maxim Imakaev,Tobias Ragoczy,Agnes Tell,Ido Amit等。2009年。“远程相互作用的全面映射揭示了人类基因组的折叠原理。”科学326:289–93。https://doi.org/10.1126/science.1181369。
Yang,Tao,Feipeng Zhang,Galip Gurkan Yardimci,歌迷歌曲,Ross C Hardison,William Stafford Noble,Feng Yue和Qunhua Li。2017年。“打hicrep:使用层调整后的相关系数评估HI-C数据的可重复性。”基因组研究。https://doi.org/10.1101/gr.220640.117。