-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathREADME.Rmd
155 lines (117 loc) · 8.17 KB
/
README.Rmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
---
output: github_document
---
<!-- README.md is generated from README.Rmd. Please edit that file -->
```{r, include = FALSE}
knitr::opts_chunk$set(
collapse = TRUE,
comment = "#>",
fig.path = "man/figures/README-",
out.width = "100%"
)
```
# preciseTAD: A transfer learning framework for 3D domain boundary prediction at base-pair resolution
<!-- [](https://travis-ci.com/stilianoudakis/preciseTAD) -->
<!-- badges: start -->
[](https://www.tidyverse.org/lifecycle/#maturing)
<!-- badges: end -->
**preciseTAD: A transfer learning framework for 3D domain boundary prediction at base-pair resolution** Spiro C. Stilianoudakis, Maggie A. Marshall, Mikhail G. Dozmorov. bioRxiv 2020.09.03.282186; doi: https://doi.org/10.1101/2020.09.03.282186
Predicted preciseTAD boundary points (PTBPs) and regions (PTBRs) for 60 cell lines are available [here](https://drive.google.com/drive/folders/15Rc6PhrrBjThwE-5dSyNX-ILELaUu6uG?usp=sharing).
## Overview
preciseTAD provides functions to predict the location of boundaries of topologically associated domains (TADs) and chromatin loops at base-level resolution. As an input, it takes BED-formatted genomic coordinates of domain boundaries detected from low-resolution Hi-C data, and coordinates of high-resolution genomic annotations from ENCODE or other consortia. preciseTAD employs several feature engineering strategies and resampling techniques to address class imbalance, and trains an optimized random forest model for predicting low-resolution domain boundaries. Translated on a base-level, preciseTAD predicts the probability for each base to be a boundary. Density-based clustering and scalable partitioning techniques are used to detect precise boundary regions and summit points. Compared with low-resolution boundaries, preciseTAD boundaries are highly enriched for CTCF, RAD21, SMC3, and ZNF143 signal and more conserved across cell lines. The pre-trained model can accurately predict boundaries in another cell line using CTCF, RAD21, SMC3, and ZNF143 annotation data for this cell line.
The main functions (in order of implementation) are:
- `extractBoundaries()` accepts a 3-column data.frame or matrix with the chromosomal coordinates of user-defined domains and outputs the unique boundaries. The second and third columns are the domain anchor centers.
- `bedToGRangesList()` accepts a filepath containing BED files representing the coordinates of ChIP-seq defined functional genomic annotations
- `createTADdata()` accepts a set of unique boundaries and genomic annotations derived from `extractBoundaries()` and `bedToGRangesList()`, respectively, to create the data matrix used to build a model to predict domain boundary regions
- `TADrandomForest()` a wrapper of the `randomForest` package which implements a random forest binary classification algorithm on domain boundary data
- `preciseTAD()` which leverages a domain boundary prediction model (i.e., random forest) and density-based clustering to predict TAD boundary coordinates at a base-level resolution
## Installation
`preciseTAD` can be installed from Bioconductor:
```{r}
# if (!requireNamespace("BiocManager", quietly=TRUE))
# install.packages("BiocManager")
# BiocManager::install("preciseTAD")
library(preciseTAD)
```
The latest version of `preciseTAD` can be directly installed from Github:
```{r eval=FALSE}
devtools::install_github("dozmorovlab/preciseTAD", build_vignettes = TRUE)
library(preciseTAD)
```
## Usage
Below is a brief workflow of how to implement `preciseTAD` on binned data from CHR1 to get precise base pair coordinates of TAD boundaries for a 10mb section of CHR 22. For more details, including the example how to use the pre-trained model, see `vignette("preciseTAD")`
First, you need to obtain called TAD boundaries using an established TAD-caller. As an example, consider the [Arrowhead](https://github.com/aidenlab/juicer/wiki/Arrowhead) TAD-caller, a part of the juicer suite of tools developed by the Aiden Lab. Arrowhead outputs a .txt file with the chromosomal start and end coordinates of their called TADs. As an example, we have provided Arrowhead TADs for GM12878 at 5kb resolution.
```{r}
data("arrowhead_gm12878_5kb")
head(arrowhead_gm12878_5kb)
```
The unique boundaries for CHR1 and CHR22 can be extracted as:
```{r}
bounds <- extractBoundaries(domains.mat = arrowhead_gm12878_5kb, filter = FALSE, CHR = c("CHR1", "CHR22"), resolution = 5000)
bounds
```
Next, you will need to download cell line-specific ChIP-seq data in the form of BED files from [ENCODE](https://www.encodeproject.org/chip-seq-matrix/?type=Experiment&replicates.library.biosample.donor.organism.scientific_name=Homo%20sapiens&assay_title=TF%20ChIP-seq&status=released). Once, you have downloaded your preferred list of functional genomic annotations, store them in a specific file location. These files can then be converted into a GRangesList object and used for downstream modeling using the following command:
```{r eval=FALSE}
path <- "pathToBEDfiles"
tfbsList <- bedToGRangesList(filepath = path, bedList = NULL, bedNames = NULL, pattern = "*.bed", signal = 4)
```
As an example, we have already provided a GRangesList object with a variety of transcription factor binding sites specific to the GM12878 cell line. Once you load it in, you can see the list of transcription factors using the following:
```{R}
data("tfbsList")
names(tfbsList)
```
For the purposes of this example, let's focus only on CTCF, RAD21, SMC3, and ZNF143 transcription factors.
```{r}
tfbsList_filt <- tfbsList[names(tfbsList) %in% c("Gm12878-Ctcf-Broad", "Gm12878-Rad21-Haib", "Gm12878-Smc3-Sydh", "Gm12878-Znf143-Sydh")]
```
Now, using the “ground-truth” boundaries and the following TFBS, we can build the data matrix that will be used for predictive modeling. The following command creates the training data from CHR1 and reserves the testing data from CHR22. We specify 5kb sized genomic bins (to match the resolution used to call the original TADs), a distance-type feature space, and apply random under-sampling (RUS) on the training data only.
```{r}
set.seed(123)
tadData <- createTADdata(bounds.GR = bounds,
resolution = 5000,
genomicElements.GR = tfbsList_filt,
featureType = "distance",
resampling = "rus",
trainCHR = "CHR1",
predictCHR = "CHR22"
)
```
We can now implement our machine learning algorithm of choice to predict TAD-boundary regions. Here, we opt for the random forest algorithm.
```{r}
set.seed(123)
tadModel <- TADrandomForest(trainData = tadData[[1]],
testData = tadData[[2]],
tuneParams = list(mtry = 2, ntree = 500, nodesize = 1),
cvFolds = 3,
cvMetric = "Accuracy",
verbose = TRUE,
model = TRUE,
importances = TRUE,
impMeasure = "MDA",
performances = TRUE)
# The model itself
tadModel[[1]]
# Variable importances (mean decrease in accuracy)
tadModel[[2]]
# Model performance metrics
tadModel[[3]]
```
Lastly, we take our TAD-boundary region predictive model and use it to make predictions on a 10mb section of CHR22:35,000,000-45,000,000.
```{r}
# Run preciseTAD
set.seed(123)
pt <- preciseTAD( genomicElements.GR = tfbsList_filt,
featureType = "distance",
CHR = "CHR22",
chromCoords = list(35000000, 45000000),
tadModel = tadModel[[1]],
threshold = 1.0,
verbose = FALSE,
parallel = 2,
DBSCAN_params = list(30000, 100),
slope = 5000,
genome = "hg19")
# View preciseTAD predicted boundary coordinates between CHR22:35mb-45mb
pt[[2]]
```