内容

1简介

preciseTAD是一个R包，旨在将域边界调用转换为有监督的机器学习框架。preciseTAD提供完整的功能，建立预测-响应数据和选择最佳模型，精确预测边界区域。这个功能可以分为两个主要用途:

1．1输入数据

用于建模的训练/测试数据表示为一个矩阵，其中行是基因组区域，列是基因组注释，单元格包含它们之间的关联度量。用户可以选择连接来自多条染色体的基因组区域。

要创建数据的逐行维度，preciseTAD使用改变了装箱，这是一种通过将线性基因组分割为不重叠区域来确定用于建模的数据矩阵维度的策略。这个步骤对用户来说是透明的。要创建移位的箱子，染色体特定的箱子从分辨率的一半开始r，并以相等的长度间隔继续r直到染色体的末端(国防部r+r / 2)，使用hg19基因组坐标。移位的基因组箱被定义为边界区域(Y = 1)，如果它们包含一个称为边界的区域，而非边界区域(Y = 0)，从而建立二进制响应向量(Y)用于分类。直观地说，移位的箱子以原始箱子之间的边界为中心，从而捕获潜在的边界。我们发现移位分箱策略提高了模型的性能。

列维由基因组注释构成，如转录因子结合位点(TFBSs)、组蛋白修饰标记、染色质状态。的（\ (log_ {2} \))的距离，列举了从每个基因组仓中心到最近的ChIP-seq峰值区域中心的基因组距离，形成特征空间。其他特性类型选项包括二进制重叠- ChIP-seq区域是否与基因组bin重叠的指示器，数重叠-每个bin中重叠的数量，和重叠的百分比-任何bin之间的重叠百分比和所有ChIP-seq区域重叠的总宽度。提供的定制培训/测试数据格式preciseTAD允许用户实现任意二进制分类机器学习算法。

1．2preciseTAD功能和输出

preciseTAD实现了一个随机森林(RF)模型，允许调整超参数和应用特征约简。主要输入是训练和测试数据、超参数值列表、要使用的交叉验证折叠次数(如果提供了超参数值网格)以及用于优化的度量。输出包括模型对象(下游边界预测所必需的)，变量重要性值，以及在测试数据上验证模型时的性能度量列表(见表1)。然后，该模型用于预测基本级别的精确边界位置。

为了预测域边界的基本位置，可以使用preciseTAD模型应用于用上述基因组注释注释的每个碱基。使用基于密度的聚类和可伸缩的数据分区技术对基本水平的概率进行聚类，以缩小边界区域和点。首先，概率向量，\ (p_ {n_{我}}\)，则提取(\ (n_{我}\)表示染色体的长度我\ \ ())．接下来,DBSCAN(基于密度的带噪声应用的空间聚类)(Ester et al. 1996;Hahsler, Piekenbrock, and Doran 2019)应用于碱基之间成对基因组距离的矩阵\(p_{n_{i}} \ge t\),在那里\ \ (t)由用户确定的阈值。由DBSCAN识别的高度预测碱基的结果集群被称为精确的边界区域(PTBR)。为了精确地识别每个PTBR中的一个碱基，preciseTAD在每个集群中实现围绕中间体(PAM)的分区。对应的集群中介id定义为a精确的边界点(PTBP)，使其成为每个聚类PTBR中最具代表性和生物学意义的碱基。输出包括PTBPs和PTBSs的基因组坐标列表。

preciseTAD允许我们使用预训练的模型来预测细胞类型中结构域边界的精确位置，没有Hi-C数据，但有基因组注释数据。具体来说，只需要CTCF、RAD21、SMC3和ZNF143转录因子结合位点的细胞类型特异性ChIP-seq数据(BED格式)。我们提供了使用GM12878和K562基因组注释数据预训练的模型，以及用于其他细胞类型精确域边界的染色体特异性预测的Arrowhead-和peakachu -预测边界。

2开始

2.1安装

# if (!requireNamespace("BiocManager", quietly=TRUE)) # install.packages("BiocManager") BiocManager::install("preciseTAD") # For the latest version install from GitHub # BiocManager::install("dozmorovlab/preciseTAD")
library(knitr) library(e1071) library(preciseTAD)

3.实现

3．1模型建立

3.1.1数据矩阵的构造

preciseTAD要求用户提供“ground-truth”域边界的基因组坐标，以建立响应向量(Y)．作为一个例子，我们考虑从流行衍生的TAD边界箭头TAD调用器，由Aiden实验室开发的Juicebox套件工具的一部分(Durand et al. 2016)．为了得到边界，箭头使用5 kb GM12878 Hi-C数据应用于常染色体(Rao et al. 2014)（GSE63525)．一个命令行脚本的例子，用来调用TADs箭头，以及结果TAD坐标的前三列。一个更详细的实现教程箭头可以找到在这里．

箭头-c 1 \ #染色体上调用TADs对-r 5000 \ #HiC数据分辨率~/GSE63525_GM12878_insitu_primary。hic \ # .HIC文件的位置~/arrowhead_output #存储输出的位置
数据(“arrowhead_gm12878_5kb”)头(arrowhead_gm12878_5kb)
## v1 v2 v3 ## 1 1 49375000 50805000 ## 2 1 16830000 17230000 ## 31 163355000 164860000 ## 4 1 231935000 233400000 ## 5 1 149035000 149430000 ## 6 1 3995000 5505000

方法将此TAD坐标数据转换为具有唯一边界坐标的GRanges对象preciseTAD: extractBoundaries ()函数。请注意，该函数的输入是一个类似bed的数据帧，其中第二和第三列是域锚点中点的基因组坐标。这里我们只提取CHR1和CHR22的边界(CHR = c("CHR1"， "CHR22")参数)。我们指定预处理= FALSE因为我们只对所有边界感兴趣，而不是根据长度过滤tad。最后，我们指定分辨率= 5000匹配所使用的分辨率箭头(尽管这个参数会被忽略预处理= FALSE)．如下所示，总共有1901个独特的TAD边界报告箭头对于染色体1和22。

bounds <- extractBoundaries(域。mat = arrowhead_gm12878_5kb, filter = FALSE, CHR = c("CHR1"， "CHR22")， resolution = 5000) #查看唯一边界边界
## seqnames ranges strand ##    ## 787 chr1 815000-815001 * ## 9196 chr1 8900000 -890001 * ## 37 chr1 915000-915001 * ## 8923 chr1 1005000-1005001 * ## 275 chr1 1015000-1015001 * ## ... ... ... ...## 8379 chr22 50815000-50815001 * ## 16788 chr22 50935000-50935001 * ## 8382 chr22 50945000-50945001 * ## 16791 chr22 51060000-51060001 * ## 16670 chr22 51235000-51235001 * ## ------- ## seqinfo:来自未指定基因组的22个序列;没有seqlengths

为了确定最能预测TAD边界的基因组特征，preciseTAD需要功能基因组注释数据。它们用于建立特征空间(\(\textbf{X}=\{X_{1}， X_{2}， \cdots, X_{p} \}\))．细胞类型特异性基因组注释数据可从编码数据门户作为BED文件。下载了首选的功能基因组注释列表后，将它们存储在特定的文件位置(即“。/path/to/BEDfiles”)。对象将这些文件转换为GRangesList对象，并用于创建特性空间preciseTAD: bedToGRangesList ()函数。的信号参数是BED文件中包含峰值信号强度值的列，用于将元数据分配给相应的grange(仅用于下游绘图)。我们已经提供了一个包含26个GM12878细胞类型的转录因子结合位点(TFBS)的GRangesList对象。加载后，可以使用以下命令查看转录因子列表。

# path <- "。*. /path/to/BEDfiles" # tfbsList <- bedToGRangesList(filepath=path, bedList=NULL, bedNames=NULL, pattern = "*. /path/to/BEDfiles"bed"， signal=4) data("tfbsList") names(tfbsList)
##[1]“Gm12878-Atf3-Haib”“gm12878 - cfs - sydh”“Gm12878-Cmyc-Uta”“gm12878 - ctfc - broad”“Gm12878-Egr1-Haib”“Gm12878-Gabp-Haib”“gm12878 - jdd - sydh”“Gm12878-Mazab85725-Sydh”“Gm12878-Mef2a-Haib”“Gm12878-Mxi1-Sydh”“Gm12878-P300-Sydh”“Gm12878-Pol2-Haib”“Gm12878-Pu1-Haib”##[16]“Gm12878-Rad21-Haib”“Gm12878-Rfx5-Sydh”“Gm12878-Sin3a-Sydh”“gm12878 - six75 - haib”“gm12878 - sp - haib”“Gm12878-Srf-Haib”“Gm12878-Taf1-Haib”“Gm12878-Taf1-Haib”“Gm12878-Yy1-Sydh”“Gm12878-Znf143-Sydh”
tfbsList
# # GRangesList对象长度26:# # $“Gm12878-Atf3-Haib”与1677年# #农庄组织对象范围和1元数据列:# # seqnames链|覆盖范围# # < Rle > < IRanges > < Rle > | <数字> # # [1]chr1 1167429 - 1167429 * | 28.5019 # # [2] chr1 1510215 - 1510215 * | 64.7440 # # [3] chr1 1655730 - 1655730 * | 25.1686 # # [4] chr1 1709843 - 1709843 * | 15.5525 # # [5] chr1 2343938 - 2343938 * 19.6961 |  ## ... ... ... ... . ...## [1673] chrX 118699292-118699462 * | 17.6970 ## [1674] chrX 119947928-119948098 * | 13.8884 ## [1675] chrX 128481457-128481627 * | 26.1075 ## [1676] chrX 132399481-132399651 * | 20.7070 ## [1677] chrX 136520393-136520563 * | 24.4197 ## ------- # seqinfo:来自未指定基因组的24个序列;No seqlength## ##…## <25个更多的元素>

使用“ground-truth”边界和以下TFBS，我们可以构建用于预测建模的数据矩阵。的preciseTAD: createTADdata ()函数可以创建训练和测试数据。这里，我们指定在染色体1上进行训练，在染色体22上进行测试。此外，我们指定分辨率= 5000构建5kb移位基因组箱(以匹配Hi-C数据分辨率)，featureType = "distance"对于一个\(log_2(距离+ 1)\)-type特征空间，和重新采样= "rus"在训练数据上应用随机欠采样(RUS)来平衡边界区域与非边界区域的类别。我们还指定了一个种子，以确保执行重新采样时的可重复性。结果是一个包含两个数据帧的列表:(1)重新采样(如果指定)的训练数据，(2)测试数据。

set.seed(123) tadData <- createaddata (bounds. seed)GR = bounds，分辨率= 5000，基因组元素。GR = tfbsList, featureType = "distance"， resampling = "rus"， trainCHR = "CHR1"， predictCHR = "CHR22") #训练数据查看子集tadData[[1]][1:5,1:4]
## y Gm12878-Atf3-Haib gm12878 - cfs - sydh Gm12878-Cmyc-Uta ## 5799 No 18.73205 17.73674 14.60901 ## 32930 No 21.23327 15.32696 15.03635 ## 5767 Yes 18.07127 17.99360 13.22415 ## 8504 Yes 19.24515 19.26026 18.31016 ## 39486 No 21.23396 17.28550 19.09578
#查看表是否平衡(tadData[[1]]$y)
## ##否是## 1595 1595
#查看测试数据的子集tadData[[2]][1:5,1:4]
## y Gm12878-Atf3-Haib Gm12878-Cfos-Sydh Gm12878-Cmyc-Uta ## 1 No 24.08978 24.1100 24.12327 ## 2 No 24.08937 24.1096 24.12287 ## 3 No 24.08897 24.1092 24.12248 ## 4 No 24.08857 24.1088 24.12208 ## 5 No 24.08816 24.1084 24.12169

3.1.2使用递归特征消除的特征选择

我们现在可以实现我们选择的机器学习算法来预测边界区域。在这里，我们选择随机森林算法进行二进制分类。preciseTAD提供了执行递归特征消除(RFE)的功能，作为通过使用的特征减少的一种形式preciseTAD: TADrfe ()函数的包装器rfe在脱字符号包库恩(2012)．preciseTAD: TADrfe ()使用5倍交叉验证，在数据中从2到2的最大特征数量的最佳特征子集上实现随机森林模型。我们指定精度作为性能度量。下面是一个示例。

set.seed(123) rfe_res <- TADrfe(trainData = tadData[[1]]， tuneParams = list(ntree = 500, nodesize = 1)， cvFolds = 5, cvMetric = "Accuracy"， verbose = TRUE)
rfe_res[[1]]
##变量MCC ROC Sens Spec Pos Pred Value Neg Pred Value Accuracy Kappa MCCSD ROCSD SensSD SpecSD Pos Pred ValueSD Neg Pred ValueSD AccuracySD KappaSD ## 12 0.3491082 0.7420849 0.6526646 0.6952978 0.6822854 0.6667258 0.6739812 0.3479624 0.02627631 0.010098192 0.009508308 0.02993689 0.019873426 0.008273672 0.013160545 0.02632109 ## 24 0.4613042 0.8026680 0.6984326 0.7455349 0.7163087 0.7297806 0.4595611 0.01904174 0.009974699 0.015099178 0.0227584018 0.0079491990.008769029 0.01753806 ## 38 0.507629440 0.7122257 0.7931034 0.5053292 0.01922275 0.010486133 0.022648767 0.01063061 0.008681254 0.014664805 0.010711191 0.02142238 ## 416 0.5150537 0.8263765 0.7197492 0.7943574 0.7779659 0.7394767 0.7570533 0.5141066 0.01881406 0.01489678 0.022886160 0.01696849 0.012896280 0.010802546 0.02160509 ## 526 0.5307963 0.8294258 0.7210031 0.8075235 0.7894162 0.7439457 0.5285266 0.02605025 0.0152340240.036287890 0.01480341 0.009062094 0.022029556 0.013957766 0.02791553
#查看每个CV折头中前n个特征的重要性变量(rfe_res[[2]])
##整体var变量Resample ## 1 26.76923 gm12878 - smc1 - sydh 26 Fold1 ## 2 23.70939 Gm12878-Rad21-Haib 26 Fold1 ## 3 20.08605 Gm12878-Ctcf-Broad 26 Fold1 ## 4 17.65042 Gm12878-Znf143-Sydh 26 Fold1 ## 5 14.17874 Gm12878-Pol2-Haib 26 Fold1 ## 6 13.25481 Gm12878-Gabp-Haib 26 Fold1

递归特征消除结果表明，仅考虑前四个转录因子构建随机森林预测模型时，模型精度开始趋于稳定(图1A)。在每个交叉褶皱中聚合前四个TFBS的可变重要性值后，我们看到五个褶皱中每一个最重要的前四个TFBS是SMC3、RAD21、CTCF和ZNF143(图1B)。这些是环状挤压模型的已知组成部分，该模型已被提出作为人类基因组三维结构的一种机制(Sanborn et al. 2015;Fudenberg et al. 2016;Hansen et al. 2018)．

图1。(A)递归特征消除(RFE)表明，仅考虑前4个最具预测性的转录因子结合位点(TFBS)时，性能稳定。(B) 5个交叉折叠中每个前4个TFBS的总平均变量重要性(使用精度的平均下降)。

图1．(A)递归特征消除(RFE)表明，仅考虑前4个最具预测性的转录因子结合位点(TFBS)时，性能稳定。(B) 5个交叉折叠中每个前4个TFBS的总平均变量重要性(使用精度的平均下降)。

3.1.3实现边界预测的随机森林

现在我们已经适当地减少了我们的特征空间，我们可以实现一个简单地建立在上面提到的TFBS (SMC3, RAD21, CTCF和ZNF143)上的随机森林算法。我们可以利用preciseTAD: TADrandomForest ()函数的包装器randomForest包(Breiman 2001;lilaw, Wiener和其他人，2002)．我们指定了训练和测试数据、超参数值、交叉验证折叠的数量、要考虑的性能度量(这里是精度)、初始化的可再现性种子、保留模型对象的指标、保留变量重要性的指标、要考虑的变量重要性度量(这里是平均精度下降(MDA))，以及根据测试数据保留模型性能的指标。函数返回一个包含以下内容的列表:tadModel [[1]]2) adata.frame模型中包含的每个特征的重要性各不相同(tadModel [[2]])， 3) adata.frame各种模型的性能指标。

#限制数据矩阵只包含SMC3、RAD21、CTCF和ZNF143 tfbsList_filt <- tfbsList[names(tfbsList) %in% c("Gm12878-Ctcf-Broad"， "Gm12878-Rad21-Haib"， "Gm12878-Smc3-Sydh"， "Gm12878-Znf143-Sydh")] set.seed(123) tadData <- createTADdata(bounds. dbdata)GR = bounds，分辨率= 5000，基因组元素。GR = tfbsList_filt, featureType = "distance"， resampling = "rus"， trainCHR = "CHR1"， predictCHR = "CHR22") #运行RF set.seed(123) tadModel <- TADrandomForest(trainData = tadData[[1]]， testData = tadData[[2]]， tuneParams = list(mtry = 2, ntree = 500, nodesize = 1)， cvFolds = 3, cvMetric = "Accuracy"， verbose = FALSE, model = TRUE, importances = TRUE, impMeasure = "MDA"，性能= TRUE)
#查看模型性能性能<- tadModel[[3]]性能$Performance <- round(性能$Performance, digits = 2) rownames(性能)<-性能$Metric kable(t(性能)，标题= "在CHR22测试数据上验证基于CHR1的RF时的模型性能列表。")

表1: 在CHR22测试数据上验证基于CHR1的RF时的模型性能列表。
TN FN 《外交政策》 TP 总计 灵敏度 特异性 卡巴 精度 BalancedAccuracy 精度 玻璃钢 FNR 净现值 世纪挑战集团 F1 AUC Youden AUPRC
度规 TN FN 《外交政策》 TP 总计 灵敏度 特异性 卡巴 精度 BalancedAccuracy 精度 玻璃钢 FNR 净现值 世纪挑战集团 F1 AUC Youden AUPRC
性能 6400.00 67.00 2954.00 239.00 9660.00 0.78 0.68 0.08 0.69 0.73 0.07 0.32 0.01 0.99 0.17 0.14 0.80 0.47 0.09

正如您可能知道的那样，在这种设置中可以使用其他机器学习二进制分类器。例如，假设我们选择实现一个支持向量机(SVM)。这很容易做到preciseTAD: createTADdata ()方便地设置训练和测试数据集。我们使用e1071包，使用下面的示例命令运行具有径向内核、cost = 1且gamma = 0.5的SVM。我们看到SVM模型的准确率只有0.67，而随机森林的准确率为0.69(表1)。

svmModel <- svm(y ~ .， data = tadData[[1]]， kernel = "radial"， cost = 1, gamma = 0.5) svmPreds <- predict(svmModel, tadData[[2]][， -1]， positive = "Yes") #视图混淆矩阵表(svmPreds, tadData[[2]][， 1])
## ## svmPreds否是##否6278 49 ##是3076 257

3.2精确的边界预测

回想一下，我们的模型对边界进行了分类地区，其中每个预测指的是一个宽度为5000个碱基的基因组仓。为了更精确地预测基本分辨率的边界坐标，我们可以通过使用preciseTAD: preciseTAD ()函数。从概念上讲，我们用所选的基因组注释和featureType注释每个碱基，而不是基因组箱。然后，我们将我们的模型应用于该注释矩阵，以预测每个碱基是给定相关基因组注释的边界的概率。为了使计算成本最小化，基本水平预测应该在选定的区域执行。

3.2.1之上运行preciseTAD

假设我们想要精确地预测CHR22:25,500,000-27,500,000的2Mb部分的域边界坐标。为此，我们指定chromCoords = list(25500000, 27500000)．此外，我们设置了一个概率阈值范围，包括0.975、0.99和1.0，用于构建ptbr。对于DBSCAN，我们分配了一个范围\ε(\ \)-邻域值包括5000、10000、15000、20000和30000，并保持3为MinPts价值。我们指定verbose = TRUE所以这个函数会打印出一些小的结果对于每个组合\ ((t) \ε)\)以及整体最佳组合

#运行preciseTAD set.seed(123) pt <- preciseTAD(基因组元素。GR = tfbsList_filt, featureType = "distance"， CHR = "CHR22"， chromCoords = list(17000000,18000000)， tadModel = tadModel[[1]]， threshold = 1, verbose = TRUE, parallel = NULL, DBSCAN_params = list(eps = c(30000)， # c(10000, 30000, 50000) MinPts = c(100))， # c(3,10,100) slope = 5000, genome = "hg19"， BaseProbs = TRUE)
##[1]“使用距离类型特征空间建立bp分辨率测试数据”##[1]“建立概率向量”##[1]“确定epsilon邻域(eps)和最小点数(MinPts)的最佳组合”##[1]“初始化DBSCAN为MinPts=100和eps=30000”##[1]“preciseTAD确定2个PTBRs”##[1]“建立PTBPs”##[1]“集群1出2”##[1]“集群2出2”##[1]“MinPts和eps的最佳组合是=(100,30000)”##[1]“preciseTAD识别出总共1678个碱基对，其预测概率等于或超过阈值1”“初始化DBSCAN for MinPts = 100和eps = 30000”“preciseTAD识别出2个ptbr”##[1]“建立PTBPs”##[1]“集群2中的1”##[1]“集群2中的2”
#什么被返回?ptbr -精确的TAD边界区域和ptbp -边界点名称(pt)
## [1] "preciseTADparams" "PTBR" "PTBP" "摘要"
#查看结果

当指定verbose = TRUE该函数为所有函数提供简要输出\(t \乘以\)组合，包括每个组合产生多少ptbr / ptbp。注:建议提供verbose = FALSE在实现preciseTAD使用多个\ \ (t)而且\ε(\ \)在整个染色体上。此外，我们从详细语句中看到的最佳组合\((t， \epsilon) = (1,30000)\)．

的完整输出preciseTAD: preciseTAD ()函数是一个包含4个元素的列表:1)对于t和eps的每个组合具有平均(和标准偏差)归一化富集(NE)值的数据帧(仅当为至少参数提供多个值时);所有后续总结应用于(t, eps))， 2)跨越每个精确etad预测区域的基因组坐标(PTBR)， 3)精确etad预测边界点的基因组坐标(PTBP)的最佳组合。4)包括以下统计数据的名单:PTBRWidth - PTBR宽度，PTBRCoverage -概率超过阈值的基础水平坐标与PTBRWidth的比值，DistanceBetweenPTBR -前一个PTBR结束与后续PTBR开始之间的基因组距离，NumSubRegions -每个PTBR集群中的元素数量，SubRegionWidth -跨越与每个PTBR相关的子区域的基因组坐标，DistBetweenSubRegions -前一个ptbr特定区域的结束与后续ptbr特定区域的开始之间的基因组距离，normmilizedenrichment -在模型中使用的围绕侧边ptbp的基因组注释的归一化富集，以及BaseProbs -每个对应基坐标的概率数字向量。归一化富集计算方法为与侧翼预测边界重叠的峰值区域总数除以预测边界的数量。图2给出了如何计算每个摘要的示意图。

图2。说明如何在preciseTAD()输出中计算每个诊断摘要的示意图。

图2．中如何计算每个诊断摘要的示意图preciseTAD ()输出。

作为边界的基本特定概率向量(BaseProbs)可能很大，视乎chromCoords设置。在我们的例子中，它有200万个条目。设置时要小心BaseProbs = TRUE避免内存问题。为了看看前500个碱基的概率是怎样的，我们这样做:

情节(pt总结BaseProbs美元[1∶])

我们可以看到，成为边界的概率似乎在200和400个碱基左右达到峰值。它远低于1，但可能值得进一步研究。

3.2.2使用preciseTAD与Juicebox

Juicebox是Aiden Lab提供的一个接口，允许将边界坐标叠加到Hi-C接触图上。要可视化被预测边界包围的域，必须首先选择Hi-C映射。以Rao et al. 2014导出的GM12878的联系矩阵为例，可以通过选择导入文件->打开…-> Rao and Huntley等-> GM12878 -> in situ Mbol -> primary．格式preciseTAD结果使用在Juicebox，用户可以利用preciseTAD: juicer_func ()方法中提供的函数preciseTADR包，它将GRanges对象转换为如下所示的数据帧。

#转换pt_juice <- juicer_func(pt$PTBP)

然后需要将ptbp保存为BED文件使用write.table如下图所示。一旦保存，导入到Juicebox使用显示注释面板-> 2D注释->添加Local -> pt_juice.bed．

Filepath = "~/path/to/store/ptbps"写入。表(pt_juice、文件。path(filepath， "pt_juice.bed")， quote = FALSE, col.names = FALSE, row.names = FALSE, sep = "\t")

3.2.3Cross-cell-type预测

preciseTAD允许我们使用预训练的模型来预测细胞类型中结构域边界的精确位置，没有Hi-C数据，但有基因组注释数据。具体来说，只需要CTCF、RAD21、SMC3和ZNF143转录因子结合位点的细胞类型特异性ChIP-seq数据(BED格式)。这些数据需要使用上述方法组装到GRanges对象中bedToGRangesList ()函数(例如，创建另一个单元格类型特定的tfbsList_filt对象)。

预训练模型tadModel可以通过使用bedToGRangesList ()，extractBoundaries ()，createTADdata (),TADrandomForest ()函数，如所述。或者，染色体特异性预训练模型可以从OneDrive存储 ．我们推荐使用GM12878基因组注释数据和Peakachu峰预训练的模型。注意，每个模型都带有“CHRi_holdout”注释，表示训练中持有的染色体，因此用于预测边界。需要注意的是，预训练preciseTAD在癌细胞系(GM12878和K562)上训练模型。尽管我们期望在癌细胞中保留域边界形成的基本规则(CTCF等结合)，但在应用预训练的模型来预测正常细胞系的边界时要谨慎。

作为一个例子，我们考虑Hela-S3细胞系。我们将下载并使用CHR22_GM12878_5kb_Peakachu模型来预测Hela-S3细胞系中chr22的边界。

#读取建立在chr1-chr21上的RF模型tadModel <- readRDS("~/Downloads/CHR22_GM12878_5kb_Peakachu.rds")

hela - s3特定的bed格式的TFBS数据可以在编码．所选的TFBS数据可从辅助存储库．它们可以转换为GRangesList。

#建立hela特异性基因组注释#从存储库中读取CTCF, RAD21, SMC3和ZNF143，并存储为列表CTCF <- as.data.frame(read.delim("http://raw.githubusercontent.com/stilianoudakis/preciseTAD_supplementary/master/data/bed/hela/ctcf.bed"， sep = "\t"， header = FALSE)) RAD21 <- as.data.frame(read.delim("http://raw.githubusercontent.com/stilianoudakis/preciseTAD_supplementary/master/data/bed/hela/rad21.bed"， sep = "\t"，header = FALSE)) smc3 <- as.data.frame(read.delim("http://raw.githubusercontent.com/stilianoudakis/preciseTAD_supplementary/master/data/bed/hela/smc3.bed"， sep = "\t"， header = FALSE)) znf143 <- as.data.frame(read.delim("http://raw.githubusercontent.com/stilianoudakis/preciseTAD_supplementary/master/data/bed/hela/znf143.bed"， sep = "\t"， header = FALSE)) helaList <- list(ctcf, rad21, smc3, znf143) hela。GR <- bedToGRangesList(bedList=helaList, bedNames=c("Ctcf"， "Rad21"， "Smc3"， "Znf143")， signal=5)

现在我们可以预测Hela-S3的chr22边界了。

#运行preciseTAD set.seed(123) pt <- preciseTAD(基因组元素。GR = hela。GR, featureType = "distance"， CHR = "CHR22"， chromCoords = list(25500000, 27500000)， tadModel = tadModel[[1]]， threshold = 1.0, verbose = FALSE, parallel = NULL, DBSCAN_params = list(30000, 100)， slope = 5000, BaseProbs = FALSE, savetobed = FALSE)

参考文献

布雷曼，利奥，2001年。“随机森林。”机器学习45(1): 5-32。

杜兰德，Neva C, Muhammad S Shamim, Ido Machol, Suhas SP Rao, Miriam H Huntley, Eric S Lander, Erez Lieberman Aiden. 2016。榨汁机提供一键式系统分析循环分辨率Hi-c实验电池系统3(1): 95-98。

Ester, Martin, Hans-Peter Kriegel, Jörg Sander，徐晓伟等。1996.“一种基于密度的在有噪声的大型空间数据库中发现聚类的算法。”在知识发现(Kdd), 96:226-31。34.

Fudenberg, Geoffrey, Maxim Imakaev, Carolyn Lu, Anton Goloborodko, Nezar Abdennur, Leonid A Mirny. 2016。环状挤压形成染色体结构域细胞的报道15(9): 2038-49。

哈斯勒，迈克尔，马修·皮肯布洛克，德里克·多兰，2019年。Dbscan:基于密度的快速聚类统计软件杂志25: 409 - 16。

汉森，安德斯·S，克劳迪娅·卡特里奥，泽维尔·达扎克，罗伯特·吉安。2018。“最近的证据表明，Tads和染色质环是动态结构。”核9(1): 20-32。

马克斯·库恩，2012年。“插入式包装。”R统计计算基金会，维也纳，奥地利。URL Https://Cran。r项目。Org/Package =插入符号．

lilaw, Andy, Matthew Wiener等人。2002.随机森林的分类与回归R新闻2(3): 18-22。

Rao, Suhas SP, Miriam H Huntley, Neva C Durand, Elena K Stamenova, Ivan D Bochkov, James T Robinson, Adrian L Sanborn，等。2014。“千碱基分辨率的人类基因组3D地图揭示了染色质循环的原理。”细胞159(7): 1665-80。

Sanborn, Adrian L, Suhas SP Rao, Huang Su-Chen, Neva C Durand, Miriam H Huntley, Andrew I Jewett, Ivan D Bochkov，等。2015。“染色质挤压解释了野生型和工程基因组中环和结构域形成的关键特征。”美国国家科学院院刊112 (47): e6456-e6465。