+
+
+
+
+
+In [1]:
+
+
+
+
+
+
+suppressMessages(library(ComplexHeatmap))
+suppressMessages(library(RColorBrewer))
+suppressMessages(library(tidyverse))
+suppressMessages(library(dplyr))
+suppressMessages(library(ggplot2))
+suppressMessages(library(cowplot))
+
+
+
+
+
+
+In [2]:
+
+
+
+
+
+
+df <- read.table("./DiffFootprinting/All_statistics.txt", header = TRUE) %>%
+ as.data.frame()
+
+## sum up the TF and protection score for each cell type
+df$CLP <- df$Protection_Score_CLP + df$TC_CLP
+df$CMP <- df$Protection_Score_CMP + df$TC_CMP
+df$GMP <- df$Protection_Score_GMP + df$TC_GMP
+df$HSC <- df$Protection_Score_HSC + df$TC_HSC
+df$LMPP <- df$Protection_Score_LMPP + df$TC_LMPP
+df$MEP <- df$Protection_Score_MEP + df$TC_MEP
+df$MPP <- df$Protection_Score_MPP + df$TC_MPP
+df$pDC <- df$Protection_Score_pDC + df$TC_pDC
+
+
+
+
+
+
+In [3]:
+
+
+
+
+
+
+# For some TFs, JASPAR database includes different variants of motifs, here we remove them
+df <- df %>%
+ dplyr::filter(!grepl("var", Motif)) %>%
+ textshape::column_to_rownames("Motif")
+
+
+# We only keep TFs with more than 1000 binding sites
+df <- subset(df, Num > 1000)
+
+df <- subset(df, select = c("CLP", "CMP", "GMP",
+ "HSC", "LMPP", "MEP",
+ "MPP", "pDC"))
+
+
+
+
+
+
+In [4]:
+
+
+
+
+
+
+# we can select the top 100 highly variable TFs
+df$Var <- apply(df, 1, sd)
+
+df <- df %>%
+ top_n(50, wt = Var) %>%
+ select(-Var)
+
+df <- t(scale(t(df)))
+
+options(repr.plot.width = 6, repr.plot.height = 8)
+
+p <- Heatmap(as.matrix(df),
+ name = "TF Activity",
+ cluster_columns = TRUE,
+ cluster_rows = TRUE,
+ show_row_names = TRUE,
+ rect_gp = gpar(col = "black", lwd = 0.5)
+ )
+
+p
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
diff --git a/docs/_static/hint/zeb1.png b/docs/_static/hint/zeb1.png
new file mode 100644
index 000000000..cef9280f4
Binary files /dev/null and b/docs/_static/hint/zeb1.png differ
diff --git a/docs/_static/motif_analysis/MA0098.2.Ets1.png b/docs/_static/motif_analysis/MA0098.2.Ets1.png
new file mode 100644
index 000000000..0e657acf8
Binary files /dev/null and b/docs/_static/motif_analysis/MA0098.2.Ets1.png differ
diff --git a/docs/_static/motif_analysis/logo.png b/docs/_static/motif_analysis/logo.png
new file mode 100644
index 000000000..cf2cdb897
Binary files /dev/null and b/docs/_static/motif_analysis/logo.png differ
diff --git a/docs/_static/rgt-viz/boxplot_ex.png b/docs/_static/rgt-viz/boxplot_ex.png
new file mode 100644
index 000000000..76be20e5b
Binary files /dev/null and b/docs/_static/rgt-viz/boxplot_ex.png differ
diff --git a/docs/_static/rgt-viz/heatmap_ex.png b/docs/_static/rgt-viz/heatmap_ex.png
new file mode 100644
index 000000000..e9dc6b12b
Binary files /dev/null and b/docs/_static/rgt-viz/heatmap_ex.png differ
diff --git a/docs/_static/rgt-viz/intersection_ex.png b/docs/_static/rgt-viz/intersection_ex.png
new file mode 100644
index 000000000..9fba3e12b
Binary files /dev/null and b/docs/_static/rgt-viz/intersection_ex.png differ
diff --git a/docs/_static/rgt-viz/jaccard_ex.png b/docs/_static/rgt-viz/jaccard_ex.png
new file mode 100644
index 000000000..c7215491f
Binary files /dev/null and b/docs/_static/rgt-viz/jaccard_ex.png differ
diff --git a/docs/_static/rgt-viz/lineplot_ex.png b/docs/_static/rgt-viz/lineplot_ex.png
new file mode 100644
index 000000000..73fadbb08
Binary files /dev/null and b/docs/_static/rgt-viz/lineplot_ex.png differ
diff --git a/docs/_static/rgt-viz/logo.png b/docs/_static/rgt-viz/logo.png
new file mode 100644
index 000000000..0f03bdfc1
Binary files /dev/null and b/docs/_static/rgt-viz/logo.png differ
diff --git a/docs/_static/rgt-viz/projection_ex.png b/docs/_static/rgt-viz/projection_ex.png
new file mode 100644
index 000000000..7c5db3b52
Binary files /dev/null and b/docs/_static/rgt-viz/projection_ex.png differ
diff --git a/docs/_static/rgt-viz/rgt-viz.png b/docs/_static/rgt-viz/rgt-viz.png
new file mode 100644
index 000000000..3f0363596
Binary files /dev/null and b/docs/_static/rgt-viz/rgt-viz.png differ
diff --git a/docs/_static/rgt_logo.png b/docs/_static/rgt_logo.png
new file mode 100644
index 000000000..eee973e61
Binary files /dev/null and b/docs/_static/rgt_logo.png differ
diff --git a/docs/_static/tdf/Genomic_region_test_1.png b/docs/_static/tdf/Genomic_region_test_1.png
new file mode 100644
index 000000000..2aa982a4b
Binary files /dev/null and b/docs/_static/tdf/Genomic_region_test_1.png differ
diff --git a/docs/_static/tdf/Genomic_region_test_2.png b/docs/_static/tdf/Genomic_region_test_2.png
new file mode 100644
index 000000000..989abacfe
Binary files /dev/null and b/docs/_static/tdf/Genomic_region_test_2.png differ
diff --git a/docs/_static/tdf/Promoter_test_tutorial_1.png b/docs/_static/tdf/Promoter_test_tutorial_1.png
new file mode 100644
index 000000000..1b0fd4139
Binary files /dev/null and b/docs/_static/tdf/Promoter_test_tutorial_1.png differ
diff --git a/docs/_static/tdf/Promoter_test_tutorial_2.png b/docs/_static/tdf/Promoter_test_tutorial_2.png
new file mode 100644
index 000000000..f5fe57ded
Binary files /dev/null and b/docs/_static/tdf/Promoter_test_tutorial_2.png differ
diff --git a/docs/_static/tdf/kuo_triplex_overview.png b/docs/_static/tdf/kuo_triplex_overview.png
new file mode 100644
index 000000000..c6260966f
Binary files /dev/null and b/docs/_static/tdf/kuo_triplex_overview.png differ
diff --git a/docs/_static/tdf/logo.png b/docs/_static/tdf/logo.png
new file mode 100644
index 000000000..358d5fbf5
Binary files /dev/null and b/docs/_static/tdf/logo.png differ
diff --git a/docs/_static/thor/exampledata_THOR-600x309.png b/docs/_static/thor/exampledata_THOR-600x309.png
new file mode 100644
index 000000000..a6f383ff9
Binary files /dev/null and b/docs/_static/thor/exampledata_THOR-600x309.png differ
diff --git a/docs/_static/thor/overview_THOR.png b/docs/_static/thor/overview_THOR.png
new file mode 100644
index 000000000..423960c13
Binary files /dev/null and b/docs/_static/thor/overview_THOR.png differ
diff --git a/docs/_static/thor/rgt_logo.png b/docs/_static/thor/rgt_logo.png
new file mode 100644
index 000000000..eee973e61
Binary files /dev/null and b/docs/_static/thor/rgt_logo.png differ
diff --git a/docs/conf.py b/docs/conf.py
index a1c4acb50..0be2fd035 100644
--- a/docs/conf.py
+++ b/docs/conf.py
@@ -1,279 +1,180 @@
-# -*- coding: utf-8 -*-
+# Configuration file for the Sphinx documentation builder.
#
-# RGT documentation build configuration file, created by
-# sphinx-quickstart on Wed Jan 6 15:14:11 2016.
-#
-# This file is execfile()d with the current directory set to its
-# containing dir.
-#
-# Note that not all possible configuration values are present in this
-# autogenerated file.
-#
-# All configuration values have a default; values that are commented out
-# serve to show the default.
+# This file only contains a selection of the most common options. For a full
+# list see the documentation:
+# https://www.sphinx-doc.org/en/master/usage/configuration.html
-import sys
-import os
+# -- Path setup --------------------------------------------------------------
# If extensions (or modules to document with autodoc) are in another directory,
# add these directories to sys.path here. If the directory is relative to the
# documentation root, use os.path.abspath to make it absolute, like shown here.
-sys.path.insert(0, os.path.abspath('../rgt/'))
+#
+import inspect
+import logging
+import os
+import sys
+from datetime import datetime
+from pathlib import Path
+
+from sphinx.application import Sphinx
+from sphinx.ext import autosummary
+# sys.path.insert(0, os.path.abspath('.'))
# -- General configuration ------------------------------------------------
-# If your documentation needs a minimal Sphinx version, state it here.
-#needs_sphinx = '1.0'
+needs_sphinx = "1.7"
-# Add any Sphinx extension module names here, as strings. They can be
-# extensions coming with Sphinx (named 'sphinx.ext.*') or your custom
-# ones.
extensions = [
- 'sphinx.ext.autodoc',
+ "sphinx.ext.autodoc",
+ "sphinx.ext.doctest",
+ "sphinx.ext.coverage",
+ "sphinx.ext.mathjax",
+ "sphinx.ext.autosummary",
+ "sphinx.ext.napoleon",
+ "sphinx.ext.intersphinx",
+ "sphinx.ext.githubpages",
+ "sphinx_autodoc_typehints",
+ "nbsphinx",
+ "edit_on_github",
]
-# Add any paths that contain templates here, relative to this directory.
-templates_path = ['_templates']
-
-# The suffix of source filenames.
-source_suffix = '.rst'
-
-# The encoding of source files.
-#source_encoding = 'utf-8-sig'
-# The master toctree document.
-master_doc = 'index'
+# Generate the API documentation when building
+autosummary_generate = True
+napoleon_google_docstring = False
+napoleon_numpy_docstring = True
+napoleon_include_init_with_doc = False
+napoleon_use_rtype = False
+napoleon_custom_sections = [("Params", "Parameters")]
+# -- Project information -----------------------------------------------------
# General information about the project.
-project = 'RGT'
-copyright = '2016, Manuel Allhoff, Joseph C.C. Kuo, Eduardo G. Gusmao, Ivan G. Costa'
+project = "RGT"
+author = "RGT develop team"
+title = "RGT - Regulatory Genomics Toolbox"
+copyright = f"{datetime.now():%Y}, {author}"
-# The version info for the project you're documenting, acts as replacement for
-# |version| and |release|, also used in various other places throughout the
-# built documents.
-#
-# The short X.Y version.
-version = '0.0.1'
-# The full version, including alpha/beta/rc tags.
-release = '0.0.1'
-# The language for content autogenerated by Sphinx. Refer to documentation
-# for a list of supported languages.
-#
-# This is also used if you do content translation via gettext catalogs.
-# Usually you set "language" from the command line for these cases.
-language = None
+# -- General configuration ---------------------------------------------------
+
+# Add any Sphinx extension module names here, as strings. They can be
+# extensions coming with Sphinx (named 'sphinx.ext.*') or your custom
+# ones.
+extensions = ['myst_parser']
-# There are two options for replacing |today|: either, you set today to some
-# non-false value, then it is used:
-#today = ''
-# Else, today_fmt is used as the format for a strftime call.
-#today_fmt = '%B %d, %Y'
+# Add any paths that contain templates here, relative to this directory.
+templates_path = ['_templates']
+master_doc = "index"
# List of patterns, relative to source directory, that match files and
# directories to ignore when looking for source files.
-exclude_patterns = ['_build']
+# This pattern also affects html_static_path and html_extra_path.
+exclude_patterns = ['_build', 'Thumbs.db', '.DS_Store']
-# The reST default role (used for this markup: `text`) to use for all
-# documents.
-#default_role = None
-# If true, '()' will be appended to :func: etc. cross-reference text.
-#add_function_parentheses = True
-
-# If true, the current module name will be prepended to all description
-# unit titles (such as .. function::).
-#add_module_names = True
-
-# If true, sectionauthor and moduleauthor directives will be shown in the
-# output. They are ignored by default.
-#show_authors = False
-
-# The name of the Pygments (syntax highlighting) style to use.
-pygments_style = 'sphinx'
-
-# A list of ignored prefixes for module index sorting.
-#modindex_common_prefix = []
-
-# If true, keep warnings as "system message" paragraphs in the built documents.
-#keep_warnings = False
-
-
-# -- Options for HTML output ----------------------------------------------
+# -- Options for HTML output -------------------------------------------------
# The theme to use for HTML and HTML Help pages. See the documentation for
# a list of builtin themes.
-html_theme = 'default'
-
-# Theme options are theme-specific and customize the look and feel of a theme
-# further. For a list of options available for each theme, see the
-# documentation.
-#html_theme_options = {}
-
-# Add any paths that contain custom themes here, relative to this directory.
-#html_theme_path = []
-
-# The name for this set of Sphinx documents. If None, it defaults to
-# "
+
+
+
+
+In [ ]:
+
+
+
+
+
+
+
+
+
+## Method
+HINT (**H**mm-based **I**de**N**tification of **T**ranscription factor footprints) is a framework that uses open chromatin data to identify the active transcription factor binding sites. This method is originally proposed to model the active binding sites by simultaneous analysis of DNase-seq and the ChIP-seq profiles of histone modifications on a genome-wide level (paper).
+
+The HMM has as input a normalized and a slope signal of DNase-seq and one of the histone marks. It can, therefore, detect the increase, top, and decrease regions of either histone modification and DNase signals. And we next modified HINT to allow only DNase-seq data by removing the three histone-level states and the use of bias-corrected DNase-seq signal before normalization steps.
+
+Recently, we extended HINT to ATAC-seq, a new assay to identify accessible DNA regions, taking the protocol-specificity into consideration.
+
+## Basic Usage
+We describe here how to detect footprints using HINT for ATAC-seq, DNase-seq, and histone modifications data. To perform footprinting, you need at least two files, one with the aligned reads of your chromatin data and another describing the regions to detect footprints. You can use a peak caller, such as [MACS2](https://github.com/macs3-project/MACS), to define these regions of interest.
+
+
+### Footprinting for ATAC-seq data
+Download [here](https://costalab.ukaachen.de/open_data/hint/tutorial/HINT_ATACTest.tar.gz) the example data for ATAC-seq based on chromosome 1 of the GM12878 cell. Execute the following commands to extract the data from the download file:
+
+```shell
+tar xvfz HINT_ATACTest.tar.gz
+cd HINT_ATACTest
+```
+and the below command to perform footprinting:
+
+```shell
+rgt-hint footprinting --atac-seq ATAC.bam ATACPeaks.bed
+```
+For simplicity, we use only the first 1000 peaks from chromosome 1. The above commands will output a BED file containing the footprints in your **current folder** with **footprints** as the prefix. Moreover, You can set the below arguments
+
+```shell
+--output-location=your_directory --output-prefix=your_prefix
+```
+
+to tell HINT your preferred output directory and name. Each footprint, i.e. each line of the BED file, will contain information regarding the tag-count score (number of reads) of each footprint. This score can be used as a footprint quality assessment (the higher values indicates better candidates). In addition, a file including the details of reads and footprints will also be written in the same folder of BED file.
+
+If your data is paired-end, you may want to try another model which is optimized for paired-end sequencing data:
+
+```shell
+rgt-hint footprinting --atac-seq --paired-end --output-prefix=fp_paired ATAC.bam ATACPeaks.bed
+```
+
+**Note**: HINT performs bias correction for ATAC-seq by default, so you must download the genomes following these instructions and correctly specify the genome references with the following command before footprinting:
+
+```shell
+--organism=genome_version
+```
+Currently, the default setting is **hg19**. Find here for more information.
+
+### Footprinting for DNase-seq
+You can find [here](http://134.130.18.8/open_data/hint/tutorial/HINT_DNaseTest.tar.gz) example DNase-seq data. Execute the following commands to extract the data from a compressed file:
+```shell
+tar xvfz HINT_DNaseTest.tar.gz
+cd HINT_DNaseTest
+```
+
+and the following command to call the footprints:
+```shell
+rgt-hint footprinting --dnase-seq DNase.bam DNasePeaks.bed
+```
+
+We recommend you to use cleavage bias correction. This can be done by using the following command:
+
+```shell
+rgt-hint footprinting --dnase-seq --bias-correction DNase.bam DNasePeaks.bed
+```
+
+Don’t forget to define the proper genome references using :
+```shell
+--organism=genome_version
+```
+Currently, the default setting is **hg19**.
+
+### Footprinting for histone modification data
+Download [here](http://134.130.18.8/open_data/hint/tutorial/HINT_HistoneTest.tar.gz) the example data for histone modification. Execute the following commands to extract data:
+
+```shell
+tar xvfz HINT_HistoneTest.tar.gz
+cd HINT_HistoneTest
+```
+
+and call footprints
+```shell
+rgt-hint footprinting --histone histone.bam histonePeaks.bed
+```
+
+
+## Citation
+
+If you use HINT with ATAC-seq should cite the following publication:
+
+* Li, Z., Schulz, M. H., Look, T., Begemann, M., Zenke, M., & Costa, I. G. (2019). [Identification of transcription factor binding sites using ATAC-seq](https://genomebiology.biomedcentral.com/articles/10.1186/s13059-019-1642-2). Genome Biology, 20(1), 45.
+
+HINT with DNase with bias correction should cite:
+
+* Gusmao EG, Allhoff M, Zenke M and Costa IG. [Analysis of computational footprinting methods for DNase sequencing experiments](https://www.nature.com/articles/nmeth.3772). Nature Methods, 13(4):303-309, 2016.
+
+HINT with DNase-seq or histones cite the following publication:
+
+* Gusmao EG, Dieterich C, Zenke M and Costa IG. [Detection of active transcription factor binding sites with the combination of DNase hypersensitivity and histone modifications](https://academic.oup.com/bioinformatics/article/30/22/3143/2390674). Bioinformatics, 30(22):3143-3151, 2014.
+
+
+
diff --git a/docs/hint/tutorial-dendritic-cell.md b/docs/hint/tutorial-dendritic-cell.md
new file mode 100644
index 000000000..5e5622eb8
--- /dev/null
+++ b/docs/hint/tutorial-dendritic-cell.md
@@ -0,0 +1,65 @@
+# Tutorial – Differential Footprints on Dendritic Cell specification
+
+In this tutorial, we will show how to use HINT-ATAC to compare changes in the activity of transcription factors with differential footprinting. We will use this data presented at HINT-ATAC paper, where we contrast the footprints of two dendritic cells. In short, we cultured classical DC type1 (cDC1) and plasmacytoid DC (pDC) from mouse bone marrow and performed Omni ATAC-seq experiment (original data is deposited in [GEO](https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE118221)).
+
+We have performed low-level analysis steps including reads alignment with Bowtie2 and peaks calling with MACS2.
+
+1. Download the data [here](https://costalab.ukaachen.de/open_data/hint/tutorial/DendriticCells.tar.gz) and extract the pre-processed data with the following commands:
+
+```shell
+tar xvfz DendriticCells.tar.gz
+cd DendriticCells
+```
+
+2. Execute the following commands to call footprints for cDC1 and pDC cells:
+
+```shell
+rgt-hint footprinting --atac-seq --paired-end --organism=mm9 --output-location=./ --output-prefix=cDC1 cDC1.bam cDC1_peaks.narrowPeak
+rgt-hint footprinting --atac-seq --paired-end --organism=mm9 --output-location=./ --output-prefix=pDC pDC.bam pDC_peaks.narrowPeak
+```
+
+HINT requires a bam file (indexed and sorted) and peaks as input. HINT only considers peak regions for footprinting to speed up the analysis. Each command will generate an output file in the current folder, containing the genomic locations of the footprints. HINT-ATAC (pDC.bed and cDC1.bed for the example above). It also produces a file with ending “.info”, which has statistics from the libraries, i.e. the total number of reads and so on.
+
+When –paired-end is specified, HINT-ATAC filter signals by only considering cleavage events from paired-ends with a particular size range. This option is recommended given significant improvement in footprint prediction (see our paper). You also need to specify the genome/organism with the flag –organism. See [this](https://reg-gen.readthedocs.io/en/latest/rgt/setup_data.html) for supported genomes or perform the command bellow for a complete description of HINT parameters:
+
+3. HINT also outputs signals for visualization in a genome browser. For this, you can use a special HINT command to generate genomic profiles ([BigWig](https://genome.ucsc.edu/goldenpath/help/bigWig.html) files).
+
+```shell
+rgt-hint tracks --bc --bigWig --organism=mm9 cDC1.bam cDC1_peaks.narrowPeak --output-prefix=cDC1_BC
+rgt-hint tracks --bc --bigWig --organism=mm9 pDC.bam pDC_peaks.narrowPeak --output-prefix=pDC_BC
+```
+
+This bigwig file contains the number of ATAC-seq reads at each genomic position as estimated by HINT-ATAC after signal normalization and cleavage bias correction. Note that this command currently only supports bias correction and shifts of ATAC-seq. This is more accurate than simply looking a coverage profiles of a bam file.
+
+
+Open all bigwig and footprint files (bed) generated above in a genome browser, e.g. IGV. Remember these experiments are based on mouse genome **mm9**. Check for example the genomic profiles around the gene Zbtb46 in the below, which is expressed by cDC but not by pDC. We observe that this gene has several open chromatin regions for these two cell types, but one particular region has cDC1 specific footprints.
+
+
+
+4. One of the main applications of footprinting is to find TFs associated with a particular cellular condition. We can do this by first finding motifs overlapping with predicted footprints. RGT suite also offers a tool for finding motif predicted binding sites (MPBS) with JASPAR, UNIPROBE or HOCOMOCO motifs.
+
+Execute the following commands to do motif matching for footprints using motifs from [JASPAR](https://jaspar.genereg.net/) as default:
+
+```shell
+rgt-motifanalysis matching --organism=mm9 --input-files pDC.bed cDC1.bed
+```
+
+where pDC.bed and cDC1.bed are the footprint predictions obtained in step 2. Depending on the number of footprints and the size of the motif database, the previous command requires a long computation time. The above commands will make a folder named match and generate a BED file for each input file, containing the matched motif instances for each footprint region. The 4th column contains the motif name and the 5th column the bit-score of the motif match. Check [here](https://reg-gen.readthedocs.io/en/latest/motif_analysis/introduction.html) for more details about motif matching.
+
+5. Finally, we use HINT to generate average ATAC-seq profiles around binding sites of particular TF. This analysis allows us to inspect the chromatin accessibility for each particular TF. Moreover, by comparing the cleavage profiles from two ATAC-seq libraries (i.s. cDC1 vs pDC cells ), we can get insights on changes in binding in two cells. For this, execute the following commands:
+
+```shell
+rgt-hint differential --organism=mm9 --bc --nc 30 --mpbs-files=./match/cDC1_mpbs.bed,./match/pDC_mpbs.bed --reads-files=cDC1.bam,pDC.bam --conditions=cDC1,pDC --output-location=cDC1_pDC
+```
+
+The above command will read the motif matching files generated by step 4 and BAM files which contain the sequencing reads to perform the comparison. Note that here we specify –bc to use the bias-corrected signal (currently only ATAC-seq is supported). The command –nc allow parallel execution of the job.
+
+After the command is done, a new folder ‘cDC1_pDC’ will be created. Inside it, a scatter plot can be found, which shows the transcription factor (TF) activity dynamics between cDC1 and pDC as shown below. The y-axis represents the differences in TF activity and names of TFs with significant differential activity values are colored as red (x-axis is a random number for jittering purposes).
+
+
+
+Each point represents a factor, but we only label TFs with a significant change in activity score (p-value < 0.05). These statistics can be obtained from file cDC1_pDC_statistics.txt. From the above plot, we can see that motif BATF::JUN and Zeb1 show significant differences, suggesting that they are related to [cDC1](https://www.nature.com/articles/ni.3197) and [pDC](https://www.sciencedirect.com/science/article/pii/S0092867408011380), respectively.
+
+In addition, a folder called lineplot can be found, which contains the ATAC-seq profile for each of the motifs found in the mpbs bed files. Let’s check the profiles of these two factors, you will see that BATF has a higher number of ATAC-seq counts in cDC1 cells, while Zeb1 has more ATAC-seq in pDC cells. Higher ATAC-seq signal indicates the higher activity of the factor in that particular cell. Other known DC factors identified by HINT-ATAC include [Spi1](https://academic.oup.com/nar/article/43/20/9680/1396036?login=false) (also referred as PU.1 or SFPI1), which is a master regulator of DC differentiation and mostly active in cDC.
+
+
diff --git a/docs/hint/tutorial-single-cell.md b/docs/hint/tutorial-single-cell.md
new file mode 100644
index 000000000..e958cca61
--- /dev/null
+++ b/docs/hint/tutorial-single-cell.md
@@ -0,0 +1,106 @@
+# Tutorial – Differential Footprints on single-cell ATAC-seq
+
+Analysis of single-cell chromatin accessibility data allows us to detect novel cell clusters, while how to interpret these clusters remains a challenge. In this tutorial, we will show how to use HINT-ATAC to compare transcription factor binding activity based on differential footprinting analysis. We will use a public scATAC-seq dataset from [human hematopoietic cell](https://www.cell.com/cell/fulltext/S0092-8674(18)30446-X).
+
+## Biological backgroud
+Human hematopoietic stem cells give rise to all blood cells through a continuous differentiation which involves remarkable changes of gene regulation. The transcription factors (TFs) play an important role during this process. Below is a schematic showing the classical definition of human hematopoietic differentiation (Source: Buenrostro, Jason D., et al (2018) ).
+
+
+
+## Identifty cell-type-specfic TFs with footprinting analysis
+ATAC-seq (Assay for Transposase-Accessible Chromatin using sequencing) is a simple and sensitive protocol for profiling genome-wide chromatin accessibility based on Tn5 transposition. In this process, the DNA sequences directly bound by TFs are protected from transposition thus leaving a footprint. We will utilize this information to identify TF binding sites and compare TF activity across distinct biological conditions. Prior to this, please ensure that you have installed RGT and setup the data successfully. See [here](https://reg-gen.readthedocs.io/en/latest/rgt/installation.html) for instructions.
+
+1. Download the data [here](https://costalab.ukaachen.de/open_data/scOpen/HematopoieticCells/HematopoieticCells.tar.gz) and extract the pre-processed data with the following commands:
+
+```shell
+tar xvfz HematopoieticCells.tar.gz
+cd HematopoieticCells
+```
+
+Inside the directory, you will find a sub-folder called BAM which includes BAM filles for all the cell types that we would like to analyze.
+
+2. To predict TF binding sites, we need a peak file for each of the BAM files. Here we will do the peak calling using [MACS2](https://genomebiology.biomedcentral.com/articles/10.1186/gb-2008-9-9-r137):
+
+```shell
+## create a folder for peaks
+mkdir Peaks
+
+## loop through all BAM files and call the peaks
+for celltype in CLP CMP GMP HSC LMPP MEP MPP pDC;
+do
+macs2 callpeak -t ./BAM/${celltype}.bam -n ${celltype} --outdir ./Peaks -g hs --nomodel -f BAMPE -q 0.01 --keep-dup all
+sed -i '/random\|chrUn\|chrM/d'./Peaks/${celltype}_peaks.narrowPeakdone
+done
+```
+
+3. Execute the following commands to call footprints for all the cell types:
+```shell
+mkdir Footprints
+
+for celltype in CLP CMP GMP HSC LMPP MEP MPP pDC;
+do
+rgt-hint footprinting --atac-seq --paired-end --organism=hg19 \
+--output-location=./Footprints \
+--output-prefix=${celltype} \
+./BAM/${celltype}.bam ./Peaks/${celltype}_peaks.narrowPeak
+done
+```
+
+HINT only considers peak regions for footprinting to speed up the analysis. For each cell type, this will generate a BED file containing the genomic locations of the footprints under the Footprints folder. It also produces a file with ending “.info”, which has statistics from the libraries, i.e. the total number of reads and so on.
+
+When –paired-end is specified, HINT-ATAC filter signals by only considering cleavage events from paired-ends with a particular size range. This option is recommended given a significant improvement in footprint prediction (see our paper). You also need to specify the genome/organism with the flag –organism. See [this](https://reg-gen.readthedocs.io/en/latest/rgt/setup_data.html) for supported genomes or perform the command bellow for a complete description of HINT parameters:
+
+```shell
+rgt-hint footprinting --help
+```
+
+**NOTE**: this is the most time-consuming step, and to speed it up, you can run the above command for each cell type at the same time.
+
+4. One of the main applications of footprinting is to find TFs associated with a particular cellular condition. We can do this by first finding motifs overlapping with predicted footprints. RGT suite also offers a tool for finding motif predicted binding sites (MPBS) with JASPAR, UNIPROBE or HOCOMOCO motifs. Execute the following commands to do motif matching for footprints using motifs from [JASPAR](http://jaspar.genereg.net/) as default:
+
+```shell
+mkdir MotifMatching
+
+for celltype in CLP CMP GMP HSC LMPP MEP MPP pDC;
+do
+rgt-motifanalysis matching \
+--organism=hg19 \
+--output-location=./MotifMatching \
+--input-files ./Footprints/${celltype}.bed
+done
+```
+
+The above commands will generate a BED file for each input file, containing the matched motif instances for each footprint region. The 4th column contains the motif name and the 5th column the bit-score of the motif match. Check [here](https://reg-gen.readthedocs.io/en/latest/motif_analysis/introduction.html) for more details about motif matching.
+
+5. Finally, we use HINT to generate average ATAC-seq profiles around binding sites of particular TF. This analysis allows us to inspect the chromatin accessibility for each particular TF. Moreover, by comparing the cleavage profiles from different ATAC-seq libraries, we can get insights on changes in binding in two cells. For this, execute the following commands:
+
+```shell
+mkdir -p DiffFootprinting
+
+rgt-hint differential --organism=hg19 --bc --nc 64 \
+--mpbs-files=./MotifMatching/CLP_mpbs.bed,\
+./MotifMatching/CMP_mpbs.bed,\
+./MotifMatching/GMP_mpbs.bed,\
+./MotifMatching/HSC_mpbs.bed,\
+./MotifMatching/LMPP_mpbs.bed,\
+./MotifMatching/MEP_mpbs.bed,\
+./MotifMatching/MPP_mpbs.bed,\
+./MotifMatching/pDC_mpbs.bed \
+--reads-files=./BAM/CLP.bam,./BAM/CMP.bam,./BAM/GMP.bam,\
+./BAM/HSC.bam,./BAM/LMPP.bam,./BAM/MEP.bam,./BAM/MPP.bam,\
+./BAM/pDC.bam \
+--conditions=CLP,CMP,GMP,HSC,LMPP,MEP,MPP,pDC \
+--output-location=./DiffFootprinting --output-prefix=All
+```
+
+The above command will read the motif matching files generated by last step and BAM files which contain the ATAC-seq to estimate a cell-type-specific binding activity score for each TF. We can use this to identify important TFs as shown below. Note that here we specify –bc to use the bias-corrected signal (currently only ATAC-seq is supported). The command –nc allows parallel execution of the job.
+
+After the command is done, a new folder **DiffFootprinting** will be created. Under this directory, there is a sub-folder called Lineplots which includes all the line plots for visualization of the TF footprint profiles. For example, we can check the activity of GATA1 which is know to be relevant for megakaryocyte-erythroid progenitro (MEP) cells differentiation:
+
+
+
+6. To summarize all TF activity between distinct conditions, we can generate a heatmap using the data from **All_statistics**.txt. This file contains a tag count and protection score for each conditon and cell type.
+
+
+
+Based this result, we can identify more cell-type-specific TFs. The code for generating the heatmap is found [here](../_static/hint/plot_heatmap.html).
\ No newline at end of file
diff --git a/docs/index.rst b/docs/index.rst
index 5056a2790..d6e6b7617 100644
--- a/docs/index.rst
+++ b/docs/index.rst
@@ -1,56 +1,80 @@
.. RGT documentation master file, created by
- sphinx-quickstart on Wed Jan 6 15:14:11 2016.
+ sphinx-quickstart on Mon Jul 4 13:14:44 2022.
You can adapt this file completely to your liking, but it should at least
contain the root `toctree` directive.
-Welcome to RGT's documentation!
-===============================
+RGT - Regulatory Genomics Toolbox
+===================================
-Contents:
+.. image:: _static/rgt_logo.png
+ :width: 175px
+ :align: left
-.. toctree::
- :maxdepth: 2
-
-.. automodule:: GenomicRegion
- :members:
-
-.. automodule:: GenomicRegionSet
- :members:
-
-.. automodule:: CoverageSet
- :members:
+RGT is an open source Python 3.6+ library for analysis of regulatory genomics.
+It is programmed in an oriented object fashion and its core classes provide functionality for handling regulatory genomics data.
-.. automodule:: MultiCoverageSet
- :members:
+This project is developed and maintained by `Institute for Computational Genomics 
+
+In this section we will describe how to create the motif logos and the various ways they can be added to the final HTML report. We will also explain how to add additional motif repositories, both from known sources (JASPAR, Hocomoco) and your own hand-crafted ones.
+
+## Create Motif Logos
+```shell
+# move inside your RGTDATA folder
+cd ~/rgtdata
+
+# run this to generate logos for all available repositories
+python setupLogoData.py --all
+```
+
+After this, your directory tree should look somewhat like this:
+```shell
+├── fig
+├── fp_hmms
+├── hg19
+├── hg38
+├── lib
+├── logos
+│ ├── hocomoco
+│ ├── jaspar_vertebrates
+│ ├── uniprobe_primary
+│ └── uniprobe_secondary
+├── mm10
+├── mm9
+├── motifs
+│ ├── hocomoco
+│ ├── jaspar_vertebrates
+│ ├── uniprobe_primary
+│ └── uniprobe_secondary
+├── zv10
+└── zv9
+```
+
+Note how the directories within **motifs** match the newly created directories within **logos**.
+
+But what if you added a new repository, **jaspar_plants**?
+```shell
+├── motifs
+│ ├── hocomoco
+│ ├── jaspar_plants
+│ ├── jaspar_vertebrates
+│ ├── uniprobe_primary
+│ └── uniprobe_secondary
+```
+
+You wouldn’t want to re-create all the logos. So you run, instead, the following command:
+
+```shell
+# generate logos only for the jaspar_plants directory,
+# as long as it exists inside the motifs directory
+python setupLogoData.py jaspar_plants
+```
+
+Afterwards, your **logos** directory will look like this:
+```shell
+├── logos
+│ ├── hocomoco
+│ ├── jaspar_plants
+│ ├── jaspar_vertebrates
+│ ├── uniprobe_primary
+│ └── uniprobe_secondary
+```
+
+Inside each of these directories you’ll see the logo images, ready for use by the Motif Analysis tool when creating the HTML report.
+
+## Add a new Motif Repository
+There are two ways to add a motif repository: an easy one, and a more complex one. We will cover both, with their pros and cons, in this section.
+
+In short, the easy way requires only adding the PWM files into a directory. The harder way requires obtaining and formatting the motif annotation.
+
+### Easy way: no annotation
+
+Go inside your **rgtdata/motifs** folder and create a directory where you’ll put the PWM files:
+```shell
+├── motifs
+│ ├── hocomoco
+│ ├── jaspar_vertebrates
+│ ├── my_custom_repository
+│ ├── uniprobe_primary
+│ └── uniprobe_secondary
+```
+
+It is important that the PWM files be named as they are in the default repositories, with the “.pwm” extension, eg:
+
+```shell
+TF1.pwm
+```
+
+The format we use is the old JASPAR format (before 2016):
+```
+13 13 3 1 54 1 1 1 0 3 2 5
+13 39 5 53 0 1 50 1 0 37 0 17
+17 2 37 0 0 52 3 0 53 8 37 12
+11 0 9 0 0 0 0 52 1 6 15 20
+```
+
+It must have exactly four rows (A, C, G, T from top to bottom – although we are adding support for two more letters, m and 1, for methylation).
+The columns are the counts for each base of the motif, left to right.
+
+**NB: if you want to add one of the JASPAR repositories, you can follow these steps:**
+
+1. Go to http://jaspar.genereg.net/downloads and download the PFMs you want (one of the files marked as **Single batch file (txt)**)
+
+2. Run the createPwm.py script we provide inside your RGTDATA folder (**motifs** subfolder), like this:
+
+```shell
+python createPwm.py -f jaspar-2016 -i jaspar_file.txt -o my_custom_repository
+```
+
+You can now create the logos as explained in the previous section:
+```shell
+python setupLogoData.py my_custom_repository
+```
+
+To use it in your experiments, you can pass it as an argument to the **matching** and **enrichment** commands:
+```
+rgt-motifanalysis matching --motif-dbs ~/rgtdata/motifs/my_custom_repository [other arguments..]
+rgt-motifanalysis enrichment --motif-dbs ~/rgtdata/motifs/my_custom_repository [other arguments..]
+```
+
+There are two main issues with using this approach:
+
+* the repositories set in your **data.config** and **data.config.user** are completely ignored
+
+* the **–filter** argument utility is much reduced, because it relies on metadata (you can still filter on motif names, though)
+
+A more organic approach is to add the repository and also add its corresponding **annotation**.
+
+### Hard way: with annotation
+
+The initial steps are the same as before, but we need to write an annotation file: **my_custom_repository.mtf**, still inside the motifs directory.
+
+But what is an MTF file? It’s a format we use to unique all the scattered annotation formats each TF database has.
+
+It is strictly a **tab separated file**, so careful with whitespaces. Let’s look at the first line in the **hocomoco.mtf** file:
+
+```shell
+AHR AHR_HUMAN.H11MO.0.B 0.B AHR PAS domain factors P35869 Integrative vertebrates Homo sapiens 2.8475,6.6065,8.544,11.1115,11.6185,12.7235
+```
+
+* The first field simply the clean “name” of the motif.
+* The second field is the full, unique, original name of the motif.
+* The third field is the version of this motif (it changes slightly across different repositories)
+* The fourth field is the TF gene name, the so-called “gene symbol”
+* The fifth field is the so-called “motif family”, a description that varies a lot across repositories.
+* The sixth field is one or more Uniprot identifiers
+* The seventh field is the data source (Chip-Seq, Selex, Integrative)
+* The eight field is the taxonomic group
+* The ninth field is the species (usually the full name, eg “Home sapiens”, not “hs”)
+* The tenth field is a list of precomputed thresholds for several FPR values: 0.005, 0.001, 0.0005, 0.0001, 0.00005, 0.00001
+
+The good news is that you don’t need to set them all. Any field you do not want to set you can replace with a single dot.
+
+This is an acceptable mtf entry:
+
+```shell
+AHR AHR_HUMAN.H11MO.0.B 0.B AHR . . . . . .
+```
+
+After you have added the metadata for each of the PWMs inside **my_custom_repository**, you can set the repository inside your **data.config.user** file:
+```shell
+[MotifData]
+repositories: my_custom_repository
+```
+
+This will make the next runs of rgt-motifanalysis use only and exclusively your custom repository. You can also add more than one repository, eg:
+```shell
+[MotifData]
+repositories: my_custom_repository, jaspar_vertebrates, hocomoco
+```
+
+and then use the **--filter** parameter as appropriate to select one or more of those, as your experiment requires (this only needed during the **matching** stage):
+
+```shell
+rgt-motifanalysis matching --filter "database:my_custom_repository,hocomoco" [other arguments..]
+```
+
+This command would only use your custom repository and hocomoco, ignoring jaspar vertebrates without the need to modify the data.config.user file.
+
+
diff --git a/docs/motif_analysis/introduction.md b/docs/motif_analysis/introduction.md
new file mode 100644
index 000000000..48ce5a94d
--- /dev/null
+++ b/docs/motif_analysis/introduction.md
@@ -0,0 +1,44 @@
+# Introduction
+
+
+
+## Method
+We provide a flexible framework to perform operations related to motif analyses. Currently, the following analyses are available:
+
+### Motif Matching
+Matches a set of transcription factor motifs found on the main repositories (such as [Jaspar](http://jaspar.genereg.net/), [Uniprobe](http://the_brain.bwh.harvard.edu/uniprobe/) and [Hocomoco](http://autosome.ru/HOCOMOCO/)) against a set of genomic regions of interest.
+
+The search is performed using a PSSM (Position-Specific Scoring Matrix), which is a normalized/background-corrected version of a simple PFM (Position Frequency Matrix), which contains the frequency that each nucleotide happens in each position of the motif. Such structure (PSSM) is matched to every position within the regions of interest in the selected genome. This will generate a score for every position in which the PSSM was applied, which can be interpreted as the binding affinity of such cis-acting element to such positions.
+
+In order to define a cutoff in which to accept/reject motif instances based on the affinity score evaluated we use a False Positive Rate approach (FPR). In this criterion, the distribution of scores for a particular PSSM is estimated based on a dynamic programming algorithm. Then a cutoff can be established based on a user-selected p-value. Best results can be achieved with p-values between 10^(-5) to 10^(-3).
+
+### Motif Enrichment
+
+Evaluate the enrichment of transcription factors in certain genomic regions. After performing motif matching of transcription factors in certain regions, we can evaluate which transcription factors are more likely to occur in those regions than in “background regions”. There are three types of test available:
+
+1. Input regions vs. Background regions: In this test, all input regions are verified against background genomic regions. These background regions are either user-provided, or random genomic regions generated during “matching” with the same average length distribution of the original input regions. It is highly recommended that the number of background regions to be at least 10 times the number of original regions for statistical significance. A [Fisher’s exact test](https://en.wikipedia.org/wiki/Fisher's_exact_test) is performed for each transcription factor with the following criteria:
+* A = number of input regions with at least 1 occurrence of that particular transcription factor.
+* B = number of input regions with no occurrence of that particular transcription factor.
+* C = number of random background regions with at least 1 occurrence of that particular transcription factor.
+* D = number of random background regions with no occurrence of that particular transcription factor.
+
+After performing the Fisher’s exact test with the variables above, the results for all transcription factors are corrected for multiple testing with the [Benjamini-Hochberg procedure](https://en.wikipedia.org/wiki/False_discovery_rate#Benjamini.E2.80.93Hochberg_procedure).
+
+2. Gene-associated regions vs. Non-gene-associated regions: In some particular analyses, we want to check whether a group of regions that are associated with genes of interest (e.g. up-regulated genes) are enriched for some transcription factors vs. regions that are not associated to those genes. In this case we perform a gene-region association in order to divide our input regions into these two groups. This association considers promoter-proximal regions, gene body and distal regions. After the association, we perform a Fisher’s exact test followed by multiple testing correction as mentioned in the previous analysis type.
+
+3. Promoter regions of input genes vs. Background regions: If a gene list is provided outside of the context of the “gene association test” (ie, NOT in the experimental matrix) then a promoter test is performed. We take all provided genes, find their promoter regions in the target organism and create a “target regions” BED file from those. Motif matching is then performed on the target regions and the provided background. In this case, three ways of specifying background are available: as a normal BED file, randomly generated or with the **--make-background** options, a newly-generated background made of the promoter regions of all genes not included in the provided gene list. Note that is usually a big task so it might take a long time, and require a lot of memory.
+
+Which analysis to use is determined by both the input and the options used. Please refer to the [tutorial](https://reg-gen.readthedocs.io/en/latest/motif_analysis/tutorial.html) and the [tool usage](https://reg-gen.readthedocs.io/en/latest/motif_analysis/tool_usage.html) page to understand how to perform either of these three.
+
+## Installation
+
+You can follow the general installation instructions for the [RGT](https://reg-gen.readthedocs.io/en/latest/rgt/installation.html) suite, which will automatically install the MotifAnalysis tool with all its dependencies.
+
+After that, you must setup your [RGT Data folder](https://reg-gen.readthedocs.io/en/latest/rgt/setup_data.html), unless you have already done it.
+
+Specifically for Motif Analysis, you also need to create the Weblogos in the RGT Data folder, otherwise they will be missing from the enrichment files. More details [here](https://reg-gen.readthedocs.io/en/latest/motif_analysis/additional_motif_data.html), but the simple command is (within the rgtdata folder):
+```shell
+python setupLogoData --all
+```
+
+Further installation instructions, including installation without pip, are found [here](https://reg-gen.readthedocs.io/en/latest/rgt/installation.html).
\ No newline at end of file
diff --git a/docs/motif_analysis/tool_usage.md b/docs/motif_analysis/tool_usage.md
new file mode 100644
index 000000000..1ed639a44
--- /dev/null
+++ b/docs/motif_analysis/tool_usage.md
@@ -0,0 +1,64 @@
+# Tool Usage
+
+To run the motif analysis tool you have to specify the analysis type (matching or enrichment). Furthermore, in order to run the motif enrichment tool you should run the motif match in advance and then the motif enrichment.
+
+## Motif Matching
+
+Run this command to see all available options:
+```shell
+rgt-motifanalysis matching --help
+```
+
+### Input
+Inputs are genomic regions in BED or BigBed format that you wish to perform motif matching on. If an experimental matrix is provided, the inputs are ignored. If a Gene List is provided, then both the inputs and the experimental matrix are ignored, since the regions will be extracted from the promoters of the organism’s genes.
+
+### Output
+This analysis will populate the specified matching folder with the following files:
+
+* **X_mpbs.[bed or bb]**: BED or BigBed file containing the matched motif instances for each input region X (or for the selected genes’ promoter regions, if **--gene-list** is provided). MPBS stands for Motif-Predicted Binding Sites. If or **--make-background** are set, then either a **random_regions_mpbs.bed** or a **background_regions_mpbs.bed** file will also be present.
+
+* **random_regions.[bed or bb]**: BED or BigBed file containing the random background regions (this file is generated only if random regions are requested with **--rand-proportion**), plus the corresponding _mpbs file.
+
+* **background_regions.[bed or bb]**: BED or BigBed file containing the genomic regions corresponding to all genes of the target organism (this file is generated only if such a background is requested with **--make-background**), plus the corresponding _mpbs file.
+
+## Motif Enrichment
+Run this command to see all available options:
+
+```shell
+rgt-motifanalysis enrichment --help
+```
+
+### Input
+A required input is the file with the regions to be used as a “background” in the enrichment statistics. The file must be in either BED or BigBed format, and it must have a corresponding MPBS file as produced from motif matching. If the background file is named “background.bed”, then the MPBS file must be named “background_mpbs.bed” and be located either in the matching location, or in the same directory as the background.
+
+Input files can be specified just like for matching, in either BED or BigBed format. If an experimental matrix is provided, input files are ignored.
+
+Alternatively, you can specify an experimental matrix that defines the regions in which the motif matching will be applied to. Learn more about the Experimental Matrix format here.
+
+In case gene sets are being used, the experimental matrix may have one additional column. This column will group regions by gene sets. See the example below:
+```
+name type file genegroup
+# Regions
+H1-hESC regions ./Input/regions_H1hesc.bed geneset1
+HeLa-S3 regions ./Input/regions_HeLaS3.bed geneset1
+HepG2 regions ./Input/regions_HepG2.bed geneset2
+K562 regions ./Input/regions_K562.bed geneset2
+# Genes
+UP_REG genes ./Input/up_regulated_genes.txt geneset1
+DW_REG genes ./Input/down_regulated_genes.txt geneset2
+```
+
+In this example, H1-hESC and HeLaS3 cell types will be associated to the up-regulated (UP_REG) genes and HepG2 and K562 will be associated to down-regulated (DW_REG) genes. Commentary lines (starting with #) were added only for clarity but not required.
+
+### Output
+This analysis will populate the folder specified with **--output-location**. This folder will contain subfolders with names depending on the analysis type. If gene sets are provided, each folder will be named X__Y, where X is the region file name and Y is the gene set file name associated to that region. Else, if no gene sets are provided, each folder will be named after each region name only. Each of these folders will contain the following files:
+
+* **coord_association.[bed or bb]**: File containing the association between gene sets and regions. If no gene sets was provided, this file will associate all regions with all genes for the organism being used.
+
+* **mpbs_ev.[bed or bb]**: Contain all Motif-Predicted Binding Sites (MPBS) that occurred in regions that were matched to genes, in case gene sets were given; or in all regions in case no gene sets was provided.
+
+* **mpbs_nev.[bed or bb]**: Contain all Motif-Predicted Binding Sites (MPBS) that occured in regions that were NOT matched to genes, in case gene sets were given. In case gene sets were not given, such output is not available.
+
+* **fulltest_statistics.[html and txt]**: Contains the enrichment results in html and tab-separated text format. The results include, in order, the motif name, enrichment p-value, corrected enrichment p-value, Fisher’s exact test A value, Fisher’s exact test B value, Fisher’s exact test C value, Fisher’s exact test D value, foreground frequency (A/(A+B)), background frequency (C/(C+D)) and a list of genes that were matched to the regions in which that motif occurred separated by comma.
+
+* **genetest_statistics.[html and txt]**: Same as the above files, but regarding the gene association test. Available only if gene sets are specified in the experimental matrix.
\ No newline at end of file
diff --git a/docs/motif_analysis/tutorial.md b/docs/motif_analysis/tutorial.md
new file mode 100644
index 000000000..b41d71b6f
--- /dev/null
+++ b/docs/motif_analysis/tutorial.md
@@ -0,0 +1,55 @@
+# Tutorial
+
+## Full-site Test
+[Download here](http://www.regulatory-genomics.org/wp-content/uploads/2017/03/RGT_MotifAnalysis_FullSiteTest.tar.gz). Execute the following commands in order to perform a motif matching followed by motif enrichment:
+
+```shell
+cd RGT_MotifAnalysis_FullSiteTest
+rgt-motifanalysis matching --input-files input/regions_K562.bed input/background.bed
+rgt-motifanalysis enrichment input/background.bed input/regions_K562.bed
+```
+
+In the enrichment step, the order of the bed files matters. The background must always come first.
+
+You can also reduce the amount of motifs used via the **–filter** parameter:
+
+```shell
+rgt-motifanalysis matching --filter "species:sapiens;name:EG" [other arguments..]
+rgt-motifanalysis enrichment --filter "species:sapiens;name:EG" [other arguments..]
+```
+
+This will restrict the search to only the motifs with the string “EG” in their name (by default, the search is inexact; use **–filter-type** to explore different modes) and whose “species” metadata contains the “sapiens” string.
+
+See the help file for all the available keys:
+
+```
+rgt-motifanalysis matching --help
+```
+
+You should also look into the metadata files – [read here](https://reg-gen.readthedocs.io/en/latest/motif_analysis/additional_motif_data.html) for more information.
+
+## Promoter Test
+
+[Download here](http://www.regulatory-genomics.org/wp-content/uploads/2017/03/RGT_MotifAnalysis_PromoterTest.tar.gz). Execute the following commands in order to perform a motif matching followed by motif enrichment:
+```shell
+cd RGT_MotifAnalysis_PromoterTest
+rgt-motifanalysis matching --target-genes input/genes.txt --input-files input/background.bed
+rgt-motifanalysis enrichment input/background.bed match/target_regions.bed
+```
+
+In the enrichment step, the order of the bed files matters. The background must always come first.
+
+### Gene association Test
+[Download here](http://www.regulatory-genomics.org/wp-content/uploads/2017/03/RGT_MotifAnalysis_GeneAssocTest.tar.gz). Execute the following commands in order to perform a motif matching followed by motif enrichment:
+
+```
+cd RGT_MotifAnalysis_GeneAssocTest
+rgt-motifanalysis matching --input-matrix input_matrix.txt --rand-proportion 10
+rgt-motifanalysis enrichment --input-matrix input_matrix.txt match/random_regions.bed
+```
+
+The matching command will read the experimental matrix, which specifies the PATH to the genomic regions and the genes to make an association test on. It also creates a background made of random regions of size 10 times the biggest genomic region in input. It might take between 10 and 30 minutes to run. It if is taking too long, use **--rand-proportion 1**.
+
+The enrichment command will calculate the enrichment statistics for the all the input regions over the random background. This step should take about 5 minutes to complete.
+
+Further usage instructions are found [here](https://reg-gen.readthedocs.io/en/latest/motif_analysis/tool_usage.html).
\ No newline at end of file
diff --git a/docs/requirements.txt b/docs/requirements.txt
new file mode 100644
index 000000000..949a1b541
--- /dev/null
+++ b/docs/requirements.txt
@@ -0,0 +1,2 @@
+setuptools==57.5.0
+myst_parser
\ No newline at end of file
diff --git a/docs/rgt-viz/introduction.md b/docs/rgt-viz/introduction.md
new file mode 100644
index 000000000..0314a1315
--- /dev/null
+++ b/docs/rgt-viz/introduction.md
@@ -0,0 +1,22 @@
+# Introduction
+
+
+
+
+RGT-Viz is a collection of tests for association analysis and tools for visualizaiton of genomic data such as genomic regions and genomic signals. RGT-Viz allows the integrative analyses of the multiple regulatory genomic data, i.e. ChIP-Seq from several histone modifications, over distinct biological conditions. A crucial aspect is the definition of a experiment matrix describing the experimental design of the underlying data. RGT-Viz results are displayed by graphics and statistics, which allow you to find association between your genomic data.
+
+RGT-Viz’s functionality can be divided in two main categories as following:
+
+* Statistical tests comparing genomic regions, i.e. peaks of distinct factors. Among others RGT-Viz implements projection test, intersection test, Jaccard test, and combinatorial test.
+
+* Tools for association of genomic regions with genomic signals. Among others, RGT-Viz allows you to perform panels with boxplots, lineplots, and heatmaps.
+
+Moreover, RGT-Viz allows to easily group data, as defined in your experimental matrix, to contrast distinct aspects of your experimental design. RGT-Viz outputs a webpage with graphical visualization and statistics as well as vector based graphics.
+
+
+
+## Install
+
+If you have followed the generic instructions for the [RGT suite installation](https://reg-gen.readthedocs.io/en/latest/rgt/installation.html), then you can start using RGT-Viz.
+
+If you have any questions, comments, installation problems or bug reports, please access our [Github](https://github.com/CostaLab/reg-gen/issues).
diff --git a/docs/rgt-viz/method.md b/docs/rgt-viz/method.md
new file mode 100644
index 000000000..57f9b8e6f
--- /dev/null
+++ b/docs/rgt-viz/method.md
@@ -0,0 +1,65 @@
+# Method description
+
+## Category 1: Tests for regions versus regions
+
+Tests for regions versus regions require two types of inputs, one is the reference regions and the other is the query regions. Both of the reference and the query are in BED format. The common aim of the below tests is to evaluate the association between the reference and the query.
+
+### Projection Test
+
+The projection test uses binomial test to evaluate the association between the reference and the query by the following processes:
+
+* Convert the query regions into spatial points by the midpoint of each region.
+* Calculate the background proportion (P) by the coverage of the reference.
+* Count the number of overlaps between the reference regions and the query points.
+* Performing binomial test by background proportion as success probability, number of overlaps as success number, and number of regions in the query as number of trials.
+
+The projection test requires large reference regions which has large coverage on the genome, because the small background proportion, which is caused by small coverage, will result in bias in binomial test.
+
+
+
+### Intersection Test
+The intersection test provides a more intuitive way to measure the association between reference and query by counting the number of intersections directly.
+
+The reference and the query should have real overlap, otherwise nothing can be calculated. One more assumption is needed: no internal overlap within the reference itself (or query itself) exists. The internal overlap within reference or query may cause bias because of the overestimation of the intersection between reference and query.
+
+
+
+### Jaccard Test
+There is another way to measure the intersection between the reference and the query, which is the Jaccard index (also called Jaccard similarity coefficient). The basic concept of the Jaccard index, given two sets A and B, is the ratio of their intersection (A∩B) to their union (A∪B). Through this Jaccard index, the amount of intersection can be expressed by a value between zero to one.
+
+
+
+### Combinatorial Test
+The intersection test can draw the inference of the association between the reference and the query directly from their intersected regions, but it cannot find the association among various combinations of the different references. For example, if the query is associated to the presence of one reference and the absence of the other reference, this association cannot be approached by the intersection test because only one reference is considered.
+
+The criteria for the reference and the query selection for the combinatorial test is different than the previous tests. The previous tests examine the association between each reference and each query. By contrast, combinatorial test uses the different combinations of the references to compare their intersection numbers among the different queries. Therefore, the sizes among references as well as the sizes among queries should not vary too much.
+
+## Category 2: Tests for regions versus signals
+
+### Boxplot
+
+Boxplot is introduced to compare sets of reads from different ChIP-seq experiments on the given regions. Based on all the intersection vectors, boxplot is generated. Conceptually, the generation of boxplot is simply counting the number of reads within the each interested regions, then plot these counts in boxplot.
+
+
+
+### Lineplot & Heatmap
+The boxplot can detect the overlaps between reads and the given region, but not their spatial distribution around the region. In order to find the signal of reads for a particular position, boxplot does not provide enough information. Therefore, lineplot and heatmap are introduced here to display the spacial distribution of reads for the particular positions.
+
+The whole process of lineplot or heatmap is described as following:
+
+1. Extend each region with the given window size which defines the boundaries for plotting. Make all regions having the same length.
+
+2. Calculate the coverage of reads on the given regions according to the given bin size and step size. The bin is a smaller fragment which moves through the window according to the defined step size. In each step, the number of reads within the bin is counted. This approach is called sliding window read counting, which generate more smooth signal than direct counting for the regions in boxplot [3]. Through this strategy, the signal can be smoothed and still preserve its local distribution.
+
+3. Plotting the lineplot or heatmap. The lineplot takes the average of the counting numbers on all the extended regions (windows), but will lose the information of single data points. The heatmap displays all information even including the counting of single bin. The lineplot shows the genome-wide trend of the signal which gets obscured by the heatmap.
+
+
+
+
+### Reference
+
+1. Alexander Favorov, Loris Mularoni, Leslie M. Cope, Yulia Medvedeva, Andrey A. Mironov, Vsevolod J. Makeev, and Sarah J. Wheelan. Exploring massive, genome scale datasets with the genometricorr package. PLoS computational biology, 8(5):1002529, 2012.
+
+2. Raimundo Real and Juan M. Vargas. The probabilistic basis of jaccard’s index of similarity. Systematic biology, pages 380–385, 1996.
+
+3. Tarjei S Mikkelsen, Manching Ku, David B Jaffe, Biju Issac, Erez Lieberman, Georgia Giannoukos, Pablo Alvarez, William Brockman, Tae-Kyung Kim, Richard P Koche, et al. Genome-wide maps of chromatin state in pluripotent and lineage-committed cells. Nature, 448(7153):553–560, 2007.
diff --git a/docs/rgt-viz/tool_usage.md b/docs/rgt-viz/tool_usage.md
new file mode 100644
index 000000000..bff2f9b8d
--- /dev/null
+++ b/docs/rgt-viz/tool_usage.md
@@ -0,0 +1,15 @@
+# Tool usage
+
+## Category 1: Tests for regions versus regions
+
+### Projection Test
+
+Commands:
+```shell
+rgt-viz projection -r 
+
+## Method
+
+Triplex Domain Finder (TDF) characterizes the triplex-forming potential between RNA and DNA regions. Our method tests whether particular RNA regions are likely to form DNA binding domains (DBD) in the given lncRNA to the potential target DNA regions (e.g. promoters of differentially regulated genes after knockdown of the lncRNA).
+
+Moreover, DNA binding sites associated with the predicted DBDs are used to indicate potential target DNA regions, e.g. genes with high binding site coverage in their promoter. TDF provides two distinct statistical tests:
+
+* **promoter test** – used to evaluate the triple helix potential in promoters
+* **genomic region test** – used to test triple helix potential in a given set of genomic regions
+
+The command line tool provides results in a user friendly and graphical HTML interface.
+
+## Installation
+
+### Installation of TDF
+
+* If you have followed the generic instructions for the [RGT suite installation](https://reg-gen.readthedocs.io/en/latest/rgt/installation.html), then you can start using TDF.
+* If you have any questions, comments, installation problems or bug reports, please report in our [Github page](https://github.com/CostaLab/reg-gen/issues).
+
+### Configure genome sequence and annotation data
+
+Before performing the test, the genome sequence and annotation data should be configured. Previous installation of wget is required for Mac OS (see [complete installation](https://reg-gen.readthedocs.io/en/latest/rgt/setup_data.html)).
+```shell
+cd ~/rgtdata
+python setupGenomicData.py --mm9
+python setupGenomicData.py --hg38
+```
+
+Note that this operation will require a few minutes for completion and will download genome FASTA files.
+
+## Example Promoter Test
+
+We will describe the necessary steps to find potential DNA binding domains of Fendrr in the promoter region of genes UP/DOWN regulated after Fendrr siRNA as presented in the [TDF publication](https://academic.oup.com/nar/advance-article/doi/10.1093/nar/gkz037/5304312?searchresult=1).
+
+To run TDF with the **promoter test**, you need as input:
+
+* the FASTA ([What is FASTA?](http://de.wikipedia.org/wiki/FASTA-Format))
+* the gene list in form of gene symbol or ensembl id
+
+You can download data for executing TDF on Fendrr example [here](https://costalab.ukaachen.de/open_data/TDF/TDF_examples.zip) and find two files in **FENDRR_mm9** directory:
+
+* Fendrr.fasta: The sequence of Fendrr in FASTA format.
+* fendrr_gene_list.txt: A list of genes UP/DOWN regulated after Fendrr siRNA.
+
+### Run promoter test
+
+Go to the directory with these files (**FENDRR_mm9/**). The promoter test can be executed with the following command:
+```shell
+rgt-TDF promotertest -r FENDRR.fasta -de fendrr_gene_list.txt -organism mm9 -rn FENDRR -o promoter_test/FENDRR/ -l 20
+```
+
+After running the command, all the result files are in the output directory (**promoter_test/Fendrr/**) with html report: promoter_test/Fendrr/index.html (see [result webpage here](https://costalab.ukaachen.de/open_data/TDF/FENDRR_mm9/promoter_test/index.html)). Complete usage instructions and more descriptive examples are found here.
+
+### Example Genomic Region Test
+The genomic region test uses any set of genomic locations as target regions. This test should be used when target DNA regions are indicated by functional studies, such as ChIRP-Seq, CHART-Seq or ChIP-Seq.
+
+Genomic region test requires two files: a RNA sequence in FASTA format and a set of interested genomic regions in BED format.
+
+Before running the test, please make sure you configure the hg38 genome data ([here](https://reg-gen.readthedocs.io/en/latest/rgt/setup_data.html)). We demonstrate the genomic region test with the following files in MEG3_hg38 from example data obtained [here](http://costalab.org/files/tdf/TDF_examples.zip):
+
+* **MEG3_sequence.fa**: RNA sequence of MEG3 in FASTA format.
+* **MEG3_hg38_CHOP.bed**: Potential target regions of MEG3.
+
+### Run genomic region test
+Go to the directory with the files and execute the following command to run the test:
+```shell
+rgt-TDF regiontest -r MEG3_sequence.fa -bed MEG3_hg38_CHOP.bed -rn MEG3 -organism hg38 -o genomic_region_test/MEG3 -n 100 -l 14
+```
+
+Then all the results and graphics are stored in the output directory (**genomic_region_test/MEG3**) (check example results page [here](https://costalab.ukaachen.de/open_data/TDF/MEG3_hg38/genomic_region_test/index.html)).
+
+Complete usage instructions and more descriptive examples are found in tutorial for genomic region test.
+
+## Example – TRIPLEXES
+
+You can also use TDF command line to make use of functionally of TRIPLEXES. First, we show how to list all triple helices for a given pair of RNA and DNA regions similar as Triplexator. We use MEG3 as the example to show how to run TRIPLEXES to find all the triplexes between lncRNA MEG3 and MEG3 CHOP regions.
+
+```shell
+cd TDF_examples/MEG3_hg38/
+rgt-TDF get_TTS -i MEG3_hg38_CHOP.bed -tts MEG3_TTSs.bed -r MEG3_sequence.fa -organism hg38
+```
+
+Then TRIPLEXES reports:
+```shell
+Total binding events: 118
+```
+
+with two files as below:
+
+* **MEG3_TTSs.bed**: A BED file for all predicted TTSs.
+* **MEG3_TTSs.txp**: A tab-delimited text file for statistics of triplexes.
+
+Another unique feature of TRIPLEXES is the detection of auto binding sites in a given genome. In this case you only need to specify the DNA regions, as TRIPLEXES will look for all DNA position forming a triple helix to a transcript potentially expressed at the same loci. Note that we allow the TFO and TTS to be shifted by a small number of bases, 5bps as default. Here we take mm9 genome as the example:
+
+```shell
+rgt-TDF autobinding -i /data/rgtdata/mm9/genome_mm9.fa
+```
+
+## Citation
+Chao-Chung Kuo, Sonja Hänzelmann, Nevcin Sentürk Cetin, Stefan Frank, Barna Zajzon, Jens-Peter Derks, Vijay Suresh Akhade, Gaurav Ahuja, Chandrasekhar Kanduri, Ingrid Grummt, Leo Kurian, Ivan G Costa; Detection of RNA–DNA binding sites in long noncoding RNAs, Nucleic Acids Research, gkz037, https://doi.org/10.1093/nar/gkz037
\ No newline at end of file
diff --git a/docs/tdf/method_description.md b/docs/tdf/method_description.md
new file mode 100644
index 000000000..4dcf35062
--- /dev/null
+++ b/docs/tdf/method_description.md
@@ -0,0 +1,25 @@
+# Method Description
+
+## General Concept
+
+The brief overview of triplexes in molecular, sequence, and celluar levels with further explanation of the computational problems:
+
+
+
+* (A) Triplexes are formed via Hoogsteen hydrogen bonds between individual nucleobases of the RNA strand (blue) and the Watson-Crick base pairs of the DNA double strand (green and orange). A A:AT triplet is shown here as example.
+* (B) The presence of homopurine-homopyrimidine sequences and the orientation of RNA dictates triplex formation. Here we show base paring rules for the anti-parallel triplex motif as an extension of the previous example (A:AT triplet).
+* (C) LncRNA can bind to DNA region via triplex formation and recruit chromatin remodeling proteins to repress (or activate) target genes. Here lncRNA Fendrr is shown as an example for binding its target gene via triplex formation and recruiting PRC2 (polycomb repressive complex 2) to trimethylates histone H3 on lysine 27 (i.e. H3K27me3).
+* (D) The binding of one single triplex between Fendrr and DNA is illustrated here to introduce Triplex Forming Oligonucleotide (TFO) and Triplex Target Site (TTS). TFO binds to TTS according to the matching motif.
+* (E) Given a lncRNA, multiple triplexes can be formed throughout the whole genome. The TFOs often cluster together and these overlapping TFOs are defined as DNA Binding Domain (DBD) where lncRNA actively form triplexes with DNA.
+
+Triplex Domain Finder (TDF) uses [TRIPLEXES](https://github.com/CostaLab/Triplexes) to find TFOs (Triplex Forming Oligonucleotides) and TTSs (Triplex Target Sites) of a lncRNA in DNA target regions. One set of regions comprises potential lncRNA targets (target regions), e.g. promoters of genes differentially expressed after FENDRR knockdown. A second set defined as non-target regions comprises promoters of genes not differentially expressed. Next, we list all candidate DNA binding domains (DBDs) by finding regions in the RNA with overlapping TFOs (see figure below). We only consider candidate DBDs with more than k TTSs in target regions. Finally, we test whether the number of target regions with at least one TTSs is higher than the number of non-target regions with at least one TTS for a given DBD.
+
+## Promoter Test
+The promoter test evaluates the triple helix binding potential of a lncRNA in the promoters of candidate genes, e.g. genes differentially expressed in a particular functional study. For this, it compares binding statistics on the promoters of candidate genes (target regions) with the binding statistics in all other promoters (non-target regions).
+
+## Genomic Region Test
+The genomic region test evaluates the triple helix binding potential of a lncRNA in a set of regions in the DNA, as indicated by genome-wide essays measuring the interaction of lncRNAs with DNA (ChIRP-Seq and CHART-Seq). Due to lack of a non-target region set, we select random genomic regions as non-target regions and use an empirical statistical test. This test evaluates whether the number of target regions with at least one DBS is larger than the number of non-target regions with at least one DBS.
+
+Exact details are provided here:
+
+Chao-Chung Kuo, Sonja Hänzelmann, Nevcin Sentürk Cetin, Stefan Frank, Barna Zajzon, Jens-Peter Derks, Vijay Suresh Akhade, Gaurav Ahuja, Chandrasekhar Kanduri, Ingrid Grummt, Leo Kurian, Ivan G Costa; Detection of RNA–DNA binding sites in long noncoding RNAs, Nucleic Acids Research, gkz037, https://doi.org/10.1093/nar/gkz037
\ No newline at end of file
diff --git a/docs/tdf/tool_usage.md b/docs/tdf/tool_usage.md
new file mode 100644
index 000000000..98bafea75
--- /dev/null
+++ b/docs/tdf/tool_usage.md
@@ -0,0 +1,140 @@
+# Tool Usage
+
+TDF includes the following submodes:
+
+* **promotertest**: Promoter test evaluates the association between the given lncRNA to the target promoters.
+* **regiontest**: Genomic region test evaluates the association between the given lncRNA to the target regions by randomization.
+* **get_dbss**: Get TTSs in BED format from the single BED file
+* **integrate**: Integrate the project’s links and generate project-level statistics.
+
+Here we introduce the common parameters for two main tests (promoter test and genomic region test) first and then describe their test-specific parameters. The last three scripts as the tools are introduced afterward.
+
+## Common Inputs for both tests
+TDF can be executed with the following command:
+```shell
+rgt-TDF {promotertest,regiontest} [required inputs] [options]
+```
+
+Where:
+
+* **{promotertest,regiontest}**: Define the applying test, either promoter test, or genomic region test.
+* **required inputs**: Required inputs files and paths.
+* **options**: Additional input parameters or output options.
+There are some inputs common for both tests shown below:
+
+### Required Input for both tests
+| **Option Name** | **Type** | **Description** |
+| :-------------- |:-------- | :--------------- |
+| -h, –help | | Show the help message and exit |
+| -r | PATH | Input file name for RNA sequence (in fasta format)|
+| -rn|String|Define the RNA name|
+|-o|PATH|Output directory name for all the results and temporary files|
+|-organism|String|Define the organism (hg19, hg38, mm9, mm10… etc)|
+
+### Options
+| **Option Name** | **Type** | **Default** | **Description** |
+| :-------------- |:-------- | :--------------- | :--------------- |
+|-t | String | RNA name |Define the title name for the results under the Output name.
+| -a | Float | 0.05 |Define significance level for rejection null hypothesis
+| -ccf |Integer | 20| Define the cut off value for promoter counts
+|-rt | Boolean | False |Remove temporary files (fa, txp…etc)
+| -log |Boolean | False |Set the plots in log scale
+|-ac |PATH | None |Input file for RNA accecibility
+|-accf |Integer | 500 |Define the cut off value for RNA accecibility
+|-obed |Boolean | False |Output the BED files for DNA binding sites.
+| -showpa | Boolean| False| Show parallel and antiparallel bindings in the plot separately
+|-filter_havana | Boolean| False| Apply filtering to remove HAVANA entries.
+|-protein_coding| Boolean| False| Apply filtering to get only protein coding genes.
+|-known_only | Boolean | False| Apply filtering to get only known genes.
+|-nofile | Boolean | False| Don’t save any files in the output folder, except the statistics.
+
+### Options for TRIPLEXES
+The arguments of the TRIPLEXES can be adjusted by the options below.
+| **Option Name** | **Type** | **Default** | **Description** |
+| :-------------- |:-------- | :--------------- | :--------------- |
+| -l | Integer | 20| Define the minimum length of triplex
+| -e | Integer | 20| Set the maximal error-rate in % tolerated
+| -c |Integer | 2 |Sets the tolerated number of consecutive errors with respect to the canonical triplex rules as were found to greatly destabilize triplexes in vitro.
+|-fr | String | off |Activates the filtering of low complexity regions and repeats in the sequence data
+|-fm |Integer | 0| Method to quickly discard non-hits (Default 0).’0′ = greedy approach; ‘1’ = q-gram filtering.
+|-of | Integer | 1 |Define output formats of Triplexator
+| -mf | Boolean | False| Merge overlapping features into a cluster and report the spanning region.
+|-rm | Boolean | False| Set the multiprocessing
+|-par | String | False |Define other parameters for TRIPLEXES. Please ignore the first “-” and replace space with underline. For example, when you want to add “-G 80 -g 20”, please do “-par G_80_-g_20”.
+
+## Particular Inputs for promoter test
+
+### Required Input for promoter test
+The target promoters can be defined in two ways:
+
+* A gene list, which contains gene symbols or Ensembl IDs, one gene per line in plain text format. The argument, -de, should be used;
+* Two BED files containing the regions of target promoters and non-target promoters (background). Two arguments, -bed and -bg, should be used together.
+
+| **Option Name** | **Type** | **Description** |
+| :-------------- |:--------| :--------------- |
+|-de | PATH| Input file for gene list (gene symbols or Ensembl ID)
+| -bed | PATH |Input BED file of the promoter regions of genes
+|-bg | PATH |Input BED file of the promoter regions of background genes
+
+### Options for promoter test
+
+| **Option Name** | **Type** | **Default** | **Description** |
+| :-------------- |:-------- | :--------------- | :--------------- |
+| -pl | Integer | 1000 |Define the promotor length
+|-score | Boolean | False| Load score column from input gene list of BED file for analysis.
+|-scoreh | Boolean | False| Use the header of scores from the given gene list or BED file.
+
+## Particular Inputs for region set test
+
+### Required Input for region set test
+| **Option Name** | **Type** | **Description** |
+| :-------------- |:--------| :--------------- |
+| -bed |PATH |Input BED file for interesting regions on DNA
+
+### Options for region set test
+-mp Integer 0Define the number of threads for multiprocessing.
+
+| **Option Name** | **Type** | **Default** | **Description** |
+| :-------------- |:-------- | :--------------- | :--------------- |
+| -n |Integer | 10000| Number iterations (randomization)
+| -f |PATH | None| Input BED file as mask for randomization
+|-score |Boolean | False| Load score column from input BED file
+
+### get_ttss
+Get TTSs of the given RNA sequence with the single BED file.
+```shell
+rgt-TDF get_ttss [options]
+```
+| **Option Name** | **Type** | **Default** | **Description** |
+| :-------------- |:-------- | :--------------- | :--------------- |
+|-h,| –help || show this help message and exit
+|-i |PATH || Input BED file of the target regions
+| -tts| PATH || Output BED file of the TTSs
+|-tfo| PATH || Output BED file of the TFOs
+|-tfo| PATH || Output BED file of the TFOs
+|-r |PATH || Input FASTA file of the RNA
+| -organism| PATH || Define the organism
+| -l | Integer | 20| **Triplexes** Define the minimum length of triplex
+|-e | Integer | 20| **Triplexes** Set the maximal error-rate in % tolerated
+| -c | Integer | 2 |**Triplexes** Sets the tolerated number of consecutive errors with respect to the canonical triplex rules as such were found to greatly destabilize triplexes in vitro
+| -fr | on/off | off |**Triplexes** Activates the filtering of low complexity regions and repeats in the sequence data
+| -fm | Integer| 0| **Triplexes** Method to quickly discard non-hits (default: 0).’0′ = greedy approach; ‘1’ = q-gram filtering.
+| -of | Integer | 1 |**Triplexes** Define output formats of Triplexes
+| -mf | Boolean | False |**Triplexes** Merge overlapping features into a cluster and report the spanning region.
+| -rm |Integer | 1| **Triplexes** Set the multiprocessing
+
+### integrate
+
+Integrate the project’s links and generate project-level statistics.
+```shell
+rgt-TDF integrate [options]
+````
+| **Option Name** | **Type** | **Default** | **Description** |
+| :-------------- |:-------- | :--------------- | :--------------- |
+| -h, –help ||| show this help message and exit
+|-path | PATH || Define the path of the project.
+|-exp | PATH || Include expression score for ranking.
+
+
+
+
diff --git a/docs/tdf/tutorial_for_genomic_region_test.md b/docs/tdf/tutorial_for_genomic_region_test.md
new file mode 100644
index 000000000..77a6dbde9
--- /dev/null
+++ b/docs/tdf/tutorial_for_genomic_region_test.md
@@ -0,0 +1,88 @@
+# Tutorial for genomic region test
+
+## Basic
+
+The genome data of hg38 should be already configured ([How to configure genome data](https://reg-gen.readthedocs.io/en/latest/rgt/setup_data.html)).
+
+Genomic region test requires one RNA sequence and a set of interesting target regions:
+
+* **RNA sequence**: A file in FASTA format. For known lncRNAs, you can download their sequences from the UCSC Genome Table Browser.
+* **Target regions**: A BED file, which contains the regions of interest.
+
+### Example data
+To demonstrate the genomic region test, we use MEG3 as an example. In the example files ([download here](https://costalab.ukaachen.de/open_data/TDF/TDF_examples.zip)), you can find the folder “MEG3_hg38“. There are several files:
+
+* **MEG3_sequence.fa**: RNA sequence of MEG3 in FASTA format.
+* **MEG3_hg38_CHOP.bed**: Potential target regions of MEG3.
+
+You can use the following command to see all the arguments related to the genomic region test:
+
+```shell
+rgt-TDF regiontest
+```
+
+### Analysis by default setting
+Given the RNA sequence and interested regions, it is easy to perform the test with the default arguments.
+Command:
+```shell
+rgt-TDF regiontest -r MEG3_sequence.fa -bed MEG3_hg38_CHOP.bed -rn MEG3 -organism hg38 -o genomic_region_test/MEG3
+```
+
+where **-r** is the lncRNA sequence, **-bed** is the region of interest, **-rn** defines the lncRNA name, **–organism** defines the organism, and **-o** is the output directory.
+
+Then, the result webpages and graphics are stored in **genomic_region_test/MEG3/**
+
+You can simply open **genomic_region_test/MEG3/index.html** to see all the results and graphics. (see Demo [here](https://costalab.ukaachen.de/open_data/TDF/MEG3_hg38/genomic_region_test/MEG3/index.html))
+
+### Result Interface
+
+There are three main pages for each test.
+
+* **RNA page**: Shows statistics and graphics of candidate DBDs.
+* **Target region page**: Shows DBSs statistics for all target regions and rankings.
+* **Parameters page**: Shows the parameters used by TDF.
+
+The main page is the RNA, which will display the main statistics regarding the candidate DBDs. Here, you can find which DBD binds significantly to the target regions.
+
+
+
+The second page shows the statistics of target regions potentially forming triplexes with a given RNA. It is possible to define the region with the highest binding possibility and its associated gene.
+
+
+
+## Advanced options for the genomic region test
+
+We will describe a few relevant options of TDF. See the tool usage for description of all options.
+
+### How to change the number of randomization?
+The default randomization processes is performed for 10,000 times. You can change it by the argument -n.
+```shell
+rgt-TDF regiontest -r MEG3_sequence.fa -bed MEG3_hg38_CHOP.bed -rn MEG3 -organism hg38 -o genomic_region_test/MEG3 -n 100
+```
+
+### How to change minimum length and error tolerance of triple helix binding sites
+When the predicted binding sites are too less to analyze, you can decrease the minimum length (-l) or increase the percentage error tolerance (-e) as follows:
+```shell
+rgt-TDF regiontest -r MEG3_sequence.fa -bed MEG3_hg38_CHOP.bed -rn MEG3 -organism hg38 -o genomic_region_test/MEG3 -n 100 -l 14 -e 20
+```
+
+### How to output target regions and DBSs in BED format
+
+Sometimes we want to output the interested regions in BED format for further investigation. TDF can output target promoters and TTSs on target promoters by adding an argument, -obed:
+```shell
+rgt-TDF regiontest -r MEG3_sequence.fa -bed MEG3_hg38_CHOP.bed -rn MEG3 -organism hg38 -o genomic_region_test/MEG3 -n 100 -obed
+```
+
+Three BED files will appear in the output directory:
+
+* **MEG3_target_region_dbs.bed**: The target regions which have at least one TTS within it.
+* **MEG3_dbss.bed**: All TTSs between lncRNA and target regions.
+T
+he file name MEG3 is defined by output directory name.
+
+### Remove temporary files
+
+TDF generates some temporary files such as .fa and .txp. They can be removed in order to save memory usage by the argument, -rt:
+```shell
+rgt-TDF regiontest -r MEG3_sequence.fa -bed MEG3_hg38_CHOP.bed -rn MEG3 -organism hg38 -o genomic_region_test/MEG3 -n 100 -obed -rt
+```
diff --git a/docs/tdf/tutorial_promoter_test.md b/docs/tdf/tutorial_promoter_test.md
new file mode 100644
index 000000000..2ca2ceeda
--- /dev/null
+++ b/docs/tdf/tutorial_promoter_test.md
@@ -0,0 +1,90 @@
+# Tutorial for promoter test
+
+## Basic
+
+The genome data of mm9 should already be configured ([How to configure genome data](file:///Users/zhijianli/Github/reg-gen/docs/_build/html/rgt/setup_data.html)).
+
+The promoter test requires an RNA sequence and a list of differentially expressed genes. Their format should be as following:
+
+* RNA sequence: A file in FASTA format. For known lncRNAs, you can download their sequence from databases, such as [UCSC Genome Table Browser](http://genome.ucsc.edu/cgi-bin/hgTables?command=start).
+* Gene list: A plain text file with one gene symbol (or Ensembl ID) per line. You may get them from your experiment or interested targets.
+
+### Example data
+For demonstrating promoter test, we use FENDRR as an example. In the example files (download [here](http://costalab.org/files/tdf/TDF_examples.zip)), you can find the folder “FENDRR_mm9“. There are three files:
+
+* **FENDRR.fasta**: RNA sequence of FENDRR in FASTA format.
+* **fendrr_gene_list.txt**: A plain text file with differentially expressed gene names.
+* **fendrr_gene_list_ford_change.txt**: A tab separated file containing gene name (first column) and absolute fold-change (second column) . The fold change information will be used as user provided experimental evidence to rank candidate target promoters.
+
+You can use the following command to see all the arguments related to promoter test.
+```shell
+rgt-TDF promotertest
+```
+
+### Analysis with default setting
+Given the example sequence and gene list, it is easy to perform the test with the default arguments.
+
+Command:
+```shell
+rgt-TDF promotertest -r FENDRR.fasta -de fendrr_gene_list_fold_change.txt -rn FENDRR -organism mm9 -o promoter_test/FENDRR/
+```
+
+where **-r** is the lncRNA sequence, **-de** is the gene list, **-rn** is the name of the lncRNA, **-organism** defines the organism and **-o** is the path to output result files. Here, we use the gene list with fold change as example.
+
+The result webpages and graphics are stored in **promoter_test/FENDRR/**
+
+You can simply open **promoter_test/FENDRR/index.html** to see all the results and graphics (or see [here](https://costalab.ukaachen.de/open_data/TDF/FENDRR_mm9/promoter_test/FENDRR_FC/index.html) an HTML output).
+
+### Result Interface
+There are three main pages for each test.
+
+* **RNA page**: Show statistics and graphics of candidate DBDs within the given RNA.
+* **Sig promoters**: Show all the statistics of triplexes formed by the significant DBDs.
+* **All promoters**: Show all the statistics of triplexes formed by the all DBDs.
+* **Autobinding**: Show the statistics of autobinding triplexes.
+* **Parameters page**: Show the parameters used by TDF.
+
+The main page is “RNA”, which will display main statistics regarding the candidate DBDs.
+
+
+
+You can see the Sig Promoters Page by selecting the corresponding tab. The page shows all target promoters ranked by a criteria of choice. Simply click at the header of each column to resort columns.
+
+
+
+## Advanced options for TDF
+We will describe a few relevant options of TDF. See the tool usage for a description of all options.
+
+### How to change the length of promoters
+The default length of promoter is 1000 bp. If you want to change it into 2000, you can add an argument, -pl, as below:
+```shell
+rgt-TDF promotertest -r FENDRR.fasta -rn FENDRR -de fendrr_gene_list.txt -organism mm9 -o promoter_test/FENDRR_2000/ -pl 2000
+```
+
+### How to change minimum length and error tolerance of binding sites
+You can also change relevant Triplexator parameters, as the size of the triple helix binding size (-l) and the maximin error tolerance (-e) as following:
+```shell
+rgt-TDF promotertest -r FENDRR.fasta -rn FENDRR -de fendrr_gene_list.txt -organism mm9 -o promoter_test/FENDRR_l15e25/ -l 15 -e 15
+```
+
+### How to output target promoters and DBSs in BED format
+TDF also generated BED files with DNA binding sites in BED format for further investigation. TDF can output target promoters and DBSs on target promoters by adding the argument, -obed:
+
+```shell
+rgt-TDF promotertest -r FENDRR.fasta -rn FENDRR -de fendrr_gene_list.txt -organism mm9 -o promoter_test/FENDRR_obed/ -obed
+```
+
+Three BED files will be generated in the output directory:
+
+* FENDRR_target_promoters.bed – target promoter generated in TDF
+* FENDRR_target_promoters_dbs.bed – target promoters with at least on DBS
+* FENDRR_dbss.bed – DNA binding sites
+
+### How to remove temporary files
+TDF generates some temporary files such as .fa and .txp. They can be removed in order to save disk space using the argument, -rt:
+
+```shell
+rgt-TDF promotertest -r FENDRR.fasta -rn FENDRR -de fendrr_gene_list.txt -organism mm9 -o promoter_test/FENDRR_rt/ -rt
+```
+
+
diff --git a/docs/thor/introduction.md b/docs/thor/introduction.md
new file mode 100644
index 000000000..9c1a5b32a
--- /dev/null
+++ b/docs/thor/introduction.md
@@ -0,0 +1,66 @@
+# Introduction
+
+
+
+## Method
+THOR is an Hidden Markov Model based approach to detect and analyze differential peaks in two sets of ChIP-seq data from distinct biological conditions with replicates. THOR performs genomic signal processing, peak calling and p-value calculation in an integrated framework.-
+
+## Basic Usage
+You can find [here](http://www.regulatory-genomics.org/wp-content/uploads/2015/07/THOR_example_data.tar.gz) example files containing reads aligned to selected genomic regions from a B Cell lymphoma [study](http://www.sciencedirect.com/science/article/pii/S1074761314004889). These files are from two donors with follicular lymphoma (FL) and two control samples (CC). We are interested in the differential peaks between the biological conditions FL and CC. The files are based on our experiment LYMP-FL-CC described in our paper.
+
+### Run THOR
+The following command finds differential peaks in the two conditions (here FL and CC):
+
+```shell
+rgt-THOR THOR.config
+```
+
+THOR.config contains all necessary files for THOR to run:
+
+```shell
+#rep1
+FL5_H3K27ac.100k.bam
+FL8_H3K27ac.100k.bam
+#rep2
+CC4_H3K27ac.100k.bam
+CC5_H3K27ac.100k.bam
+#chrom_sizes
+hg19.chrom.sizes
+```
+
+THOR needs 3 types of files:
+
+* BAM files for the first biological condition
+* BAM files for the second biological condition
+* The chromosome sizes
+
+The file hg19.chrom.sizes is contained in the example files as well. It is a tab limited file with the chromosome name and the chromosome size. See here how to get the chromosome sizes for further organisms.
+
+In the config file, we can also incorporate the genome (with header #genome) and control input-DNA files (with headers #inputs1 and #inputs2). This improves the accuracy of the differential peak estimates. See here for a more advanced example.
+
+### Results
+THOR creates several files. If we do not specify anything, THOR uses THOR-exp-
+
+## Citation
+
+Allhoff, M., Sere K., Freitas, J., Zenke, M., Costa, I.G. [Differential Peak Calling of ChIP-seq Signals with Replicates with THOR](http://nar.oxfordjournals.org/content/early/2016/08/01/nar.gkw680.abstract). Nucleic Acids Research, 44(20):e153, 2016.
diff --git a/docs/thor/tool_usage.md b/docs/thor/tool_usage.md
new file mode 100644
index 000000000..a5d99497b
--- /dev/null
+++ b/docs/thor/tool_usage.md
@@ -0,0 +1,164 @@
+# Tool Usage
+
+Command:
+
+```shell
+rgt-THOR [option] CONFIG
+```
+
+## Required Input
+| **Option Name** | **Type** | **Description** |
+| :-------------- |:-------- | :--------------- |
+| CONFIG | File | Configuration file for the experiment |
+
+### Options
+
+| **Option Name** | **Type** | **Default** | **Description** |
+| :--------------- | :-------- | :----------- | :--------------- |
+| \-\-name | string | none | Experiment’s name and prefix for all files that are created. |
+| \-\-merge | boolean | False | Merge peaks which have a distance less than the estimated mean fragment size (recommended for histone data). |
+| \-\-no-merge-bin | boolean | True | Merge the overlapping bin before filtering by p-value. |
+| \-\-housekeeping-genes | file (bed format) | none | Define housekeeping genes (BED format) used for normalizing. |
+| \-\-output-dir | string | none | Store files in output directory. |
+| \-\-report | boolean | False | Generate HTML report about experiment. |
+| \-\-deadzones | file (bed format) | none | Define blacklisted genomic regions avoided for analysis (BED format). |
+| \-\-no-correction | boolean | False | Do not use multipe test correction for p-values (Benjamini/Hochberg). |
+| \-\-pvalue | float | 0.1 | P-value cutoff for peak detection. Call only peaks with p-value lower than cutoff. |
+| \-\-exts | interger1,interger2,... | estimate | Read’s extension size for BAM files (comma separated list for each BAM file in config file). If option is not chosen, estimate extension sizes from reads. |
+| \-\-factors-inputs | float1,float2, ... | estimate | Normalization factors for input-DNA (comma separated list for each BAM file in config file). If option is not chosen, estimate factors. |
+| \-\-scaling-factors | float1,float2, ... | estimate | Scaling factor for each BAM file (not control input-DNA) as comma separated list for each BAM file in config file. If option is not chosen, follow normalization strategy (TMM or HK approach) |
+| \-\-save-input | boolean | False | Save input DNA bigwig (if input was provided). |
+
+### Advanced options
+
+| **Option Name** | **Type** | **Default** | **Description** |
+| :--------------- | :-------- | :----------- | :--------------- |
+| \-\-regions | file (bed format) | none | Define regions (BED format) to restrict the analysis, that is, where to train the HMM and search for DPs. It is faster, but less precise. |
+| \-\-binsize | integer | 100 | Size of bins for creating the signal. |
+| \-\-step | integer | 50 | Stepsize with which the window consecutively slides across the genome to create the signal. |
+| \-\-debug | boolean | False | Output debug information. Warning: space consuming! |
+| \-\-no-gc-content | boolean | False | Do not normalize towards GC content. |
+| \-\-norm-regions | file (bed format) | none | Restrict normalization to particular regions (BED format). |
+| \-\-threshold | float | 95 | Minimum signal support for differential peaks to define training set as percentage (t2, see [paper](http://nar.oxfordjournals.org/content/early/2016/08/01/nar.gkw680.abstract)). |
+| \-\-size | integer | 10000 | Number of bins the HMM's training set constists of. |
+| \-\-par | integer | 1 | Percentile for p-value postprocessing filter. |
+| \-\-poisson | boolean | False | Use binomial distribution as emmission. |
+| \-\-single-strand | boolean | False | Allow single strand BAM file as input. |
+| \-\-m\_threshold | integer | 80 | Define the M threshold of percentile for training TMM. |
+| \-\-a\_threshold | integer | 95 | Define the A threshold of percentile for training TMM. |
+| \-\-rmdup | boolean | False | Remove the duplicate reads. |
+
+### Config File
+
+The config file contains all necessary information for the experiment.
+
+The file contains several groups, each group is indicated by a header symbol (#) followed by the group’s name. The information of each group is then given in the following lines.
+
+The headers are:
+
+* #rep1
+* #rep2
+* #chrom\_sizes
+* #genome (optional)
+* #inputs1 (optional)
+* #inputs2 (optional)
+
+The headers #inputs1, #inputs2 and #genome are optional. All other headers are required. An example for a config file is given [here](https://reg-gen.readthedocs.io/en/latest/thor/introduction.html).
+
+**\#rep **
+
+Headers #rep1 and #rep2 give all BAM files to be analysed. All files are listed line by line for the first (#rep1) and second (#rep2) biological condition.
+
+**\#chrom\_sizes **
+
+The chromosome sizes is a tab limited file with the chromosome name and the size of the chromosome. To download the chromosome size file of an organism, follow these [instructions](http://wiki.bits.vib.be/index.php/FetchChromSizes).
+
+**\#inputs **
+
+For each BAM file, we can also provide control input-DNA BAM files by using headers #inputs1 and #inputs2. Input-DNA helps to handle bias in ChIP-seq profiles and can therefore improve the differential peak estimation.
+
+**\#genome **
+
+Input-DNA and the genome (in [fasta](http://blast.ncbi.nlm.nih.gov/blastcgihelp.shtml) format) is necessary to correct for GC-content. GC-content correction can lead to more precise DP estimates. For instance, to download the human genome hg19.fa, run
+
+```shell
+wget 'ftp://hgdownload.cse.ucsc.edu/goldenPath/hg19/chromosomes/*'
+gunzip *gz
+cat *fa > hg19.fa
+rm md5sum.txt chr*fa README.txt
+```
+
+## Output
+
+By default, THOR creates four types of files. The file(s)
+
+* \*.bw give the postprocessed ChIP-seq signal (in bigWig format),
+* \*-setup.info gives information about the setting,
+* \*-diffpeaks.bed describes the estimated differential peaks in a proprietary BED format, and
+* \*-diffpeaks.narrowPeak describes the estimated differential peaks in narrowPeak format.
+
+### setup.info file
+
+The setup file contains information about the experiment and the configuration of THOR used to predict DPs.
+
+This includes the extension sizes used for the reads, the normalization scaling factors for each single ChIP-seq signal and the overdispersion factor c describing the quadratic relationship between the mean and the variance in the data.
+
+Moreover, the setup.info file gives information about the HMM configuration, that is, the transition matrix as well as the emission distribution. The emission distribution is a Negative Binomial distribution described by mu and alpha.
+
+### BigWig files
+
+THOR creates a BigWig file for each ChIP-seq signal (\*s[1/2]rep[0/1/2/3/…].bw), that is, for each biological condition (s[1/2]) and each replicate (rep[0/1/2/3/…]).
+
+### narrowPeak and BED files
+
+As the narrowPeak format does not give the possibility to store the counts of condition 1 and condition 2 for each peak, we output an BED format with additional information in the 11th column.
+
+The 11th column in the BED file gives a semicolon separated list for each differential peak. The first (second) element of the semicolon separated list contains a comma separated list of the counts of each replicate of the first (second) biological conditions. The third element of the list gives the calculated p-value.
+
+For downstream analysis of this BED file, we provide two [tools](http://www.regulatory-genomics.org/wp-content/uploads/2015/07/THOR-tools.tar.gz). The first tool separates the BED file by differential peaks that gain peaks in condition 1 and that gain peaks in condition 2. The second tool filters the BED file by p-value.
+
+**Please note:** In both files, \*-diffpeaks.bed and \*-diffpeaks.narrowPeak, the strand column (6th column) indicates whether the peak gains condition 1 (positive strand) or gains condition 2 (negative strand).
+
+To output valid narrowPeak and BED files, we add several dummy columns. Only then it is possible to use them in other tools, like e.g. the IGV browser.
+
+For the BED file, columns 7 and 8, give the same information as columns 2 and 3; column 9 gives a colour code for the peaks (red for a differential peak in signal 1, and green for a differential peak in signal 2). Column 10 is always 0.
+
+For the narrowPeak file, column 5, 7, 9 and 10 do not give any information.
+
+## Advanced Example
+
+In this example, we study the Dendritic Cell development. Dendritic cells (DC) are professional antigen presenting cells that develop from hematopoietic stem cells in bone marrow. This system allows to differentiate ex-vivo mul-
+tipotent progenitors (MPP) to DC progenitors (CDP). We want to call differential peaks in our provided ChIP-seq data for histone mark H3K27ac of MPP and CDP (see [here](http://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE73143) the GEO entry for the fastq files).
+
+We create the config file DC\_THOR.config:
+
+```shell
+#rep1
+MPP_WT_H3K27ac_1.bam
+MPP_WT_H3K27ac_2.bam
+#rep2
+CDP_WT_H3K27ac_1.bam
+CDP_WT_H3K27ac_2.bam
+#genome
+mm9.fa
+#chrom_sizes
+mm9.chrom.sizes
+#inputs1
+Input_MPP_H3K27ac.bam
+Input_MPP_H3K27ac.bam
+#inputs2
+Input_CDP_H3K27ac.bam
+Input_CDP_H3K27ac.bam
+```
+
+Here, we do not give only the required (#rep1, #rep2 and #chrom\_sizes), but also the optional header (#genome, #inputs1, #inputs2) to take advantage of input-DNA and GC-content normalization.
+
+The following command finds differential peaks in the two conditions MPP and CDP.
+
+```shell
+rgt-THOR THOR_DC.config --report --housekeeping-genes hk_genes_new_promotor_mm9.bed --no-correction --output-dir ~/my_MPP_CDP_exp -n THOR_DC
+```
+
+THOR\_DC.config contains all information necessary to run THOR. With *\-\-output-dir* we determine the working directory ~/my\_MPP\_CDP\_exp. THOR stops, if the folder already exists. All files created in ~/my\_MPP\_CDP\_exp have the prefix *THOR\_DC* (specific by the option *-n*). The p-value is not corrected (option *\-\-no-correction*). Instead of using a TMM based approach to normalize the ChIP-seq profiles (default), we use a housekeeping genes approach. For that, please [download](http://www.regulatory-genomics.org/wp-content/uploads/2015/07/hk_THOR.tar.gz) the housekeeping genes and use the option *\-\-house-keeping-genes*.
+
+Moreover, we create a HTML report of our experiment (option *-\-report*). The HTML report contains useful information about our experiment, such as, the the experimental configuration (also described by _\*-setup.info_), the function THOR is using to model the mean-variance relationship, the function used to estimate the fragmentation size, and the quality check of the house-keeping-gene normalization approach.
diff --git a/requirements.txt b/requirements.txt
index 887303252..444f8c277 100644
--- a/requirements.txt
+++ b/requirements.txt
@@ -5,7 +5,6 @@ numpy
scipy
pysam
pyBigWig
-pyVCF
fisher
HTSeq
hmmlearn
diff --git a/setup.py b/setup.py
index 05264a684..a030891a7 100644
--- a/setup.py
+++ b/setup.py
@@ -74,8 +74,7 @@ def find_version(*file_paths):
"numpy>=1.4.0",
"scipy>=1.0.0",
"pysam>=0.12.0",
- "pyBigWig",
- "PyVCF"]
+ "pyBigWig"]
if platform.startswith("darwin"):
bin_dir = "mac"
@@ -252,6 +251,13 @@ def recursive_chown_chmod(path_to_walk, uid, gid, file_permission, path_permissi
data_config_file.write("# gene_regions: " + path.join(genome, "genes_RefSeq_mm10.bed # alternative to Gencode\n"))
data_config_file.write("annotation: " + path.join(genome, "gencode.vM25.annotation.gtf\n"))
data_config_file.write("gene_alias: " + path.join(genome, "alias_mouse.txt\n\n"))
+genome = "mm39"
+data_config_file.write("[" + genome + "]\n")
+data_config_file.write("genome: " + path.join(genome, "genome_mm39.fa\n"))
+data_config_file.write("chromosome_sizes: " + path.join(genome, "chrom.sizes.mm39\n"))
+data_config_file.write("gene_regions: " + path.join(genome, "genes_Gencode_mm39.bed\n"))
+data_config_file.write("annotation: " + path.join(genome, "gencode.vM25.annotation.gtf\n"))
+data_config_file.write("gene_alias: " + path.join(genome, "alias_mouse.txt\n\n"))
genome = "hg19"
data_config_file.write("[" + genome + "]\n")
data_config_file.write("genome: " + path.join(genome, "genome_hg19.fa\n"))
@@ -284,6 +290,13 @@ def recursive_chown_chmod(path_to_walk, uid, gid, file_permission, path_permissi
data_config_file.write("gene_regions: " + path.join(genome, "genes_zv10.bed\n"))
data_config_file.write("annotation: " + path.join(genome, "Danio_rerio.GRCz10.84.gtf\n"))
data_config_file.write("gene_alias: " + path.join(genome, "alias_zebrafish.txt\n\n"))
+genome = "tair10"
+data_config_file.write("[" + genome + "]\n")
+data_config_file.write("genome: " + path.join(genome, "genome_tair10_ensembl_release_51.fa\n"))
+data_config_file.write("chromosome_sizes: " + path.join(genome, "chrom.sizes.tair10\n"))
+data_config_file.write("gene_regions: " + path.join(genome, "genes_tair10.bed\n"))
+data_config_file.write("annotation: " + path.join(genome, "Arabidopsis_thaliana.TAIR10.51.gtf\n\n"))
+data_config_file.write("gene_alias: " + path.join(genome, "alias_tair10.txt\n\n"))
data_config_file.write("[MotifData]\n")
data_config_file.write("pwm_dataset: motifs\n")
@@ -354,6 +367,7 @@ def recursive_chown_chmod(path_to_walk, uid, gid, file_permission, path_permissi
"hg38": ["genes_Gencode_hg38.bed", "chrom.sizes.hg38", "alias_human.txt", "genes_RefSeq_hg38.bed"],
"mm9": ["genes_Gencode_mm9.bed", "chrom.sizes.mm9", "alias_mouse.txt", "genes_RefSeq_mm9.bed"],
"mm10": ["genes_Gencode_mm10.bed", "chrom.sizes.mm10", "alias_mouse.txt", "genes_RefSeq_mm10.bed"],
+ "mm39": ["genes_Gencode_mm39.bed", "chrom.sizes.mm39", "alias_mouse.txt"],
"zv9": ["genes_zv9.bed", "chrom.sizes.zv9", "alias_zebrafish.txt"],
"zv10": ["genes_zv10.bed", "chrom.sizes.zv10", "alias_zebrafish.txt"],
"fp_hmms": ["dnase.hmm", "dnase_bc.hmm", "histone.hmm", "dnase_histone.hmm", "dnase_histone_bc.hmm",
diff --git a/tools/update_alias.py b/tools/update_alias.py
index 3837e3d59..e55d25dc3 100644
--- a/tools/update_alias.py
+++ b/tools/update_alias.py
@@ -6,10 +6,11 @@
ALIAS_hg19 = "../data/hg19/alias.txt"
ALIAS_hg38 = "../data/hg38/alias.txt"
ALIAS_mm9 = "../data/mm9/alias.txt"
+ALIAS_mm39 = "../data/mm39/alias_mouse.txt"
# HG19 and HG38
response = urllib.request.urlopen('http://www.genenames.org/cgi-bin/download?col=gd_hgnc_id&col=gd_app_sym&col=gd_aliases&col=gd_pub_acc_ids&col=gd_pub_eg_id&col=gd_mgd_id&col=gd_pubmed_ids&col=gd_pub_refseq_ids&col=gd_ccds_ids&col=gd_vega_ids&col=md_eg_id&col=md_prot_id&col=md_ensembl_id&status=Approved&status=Entry+Withdrawn&status_opt=2&where=&order_by=gd_app_sym_sort&format=text&limit=&hgnc_dbtag=on&submit=submit')
-table = response.read()
+table = response.read().decode("utf-8")
alias = []
for i,line in enumerate(table.split("\n")):
@@ -24,7 +25,11 @@
# for line in alias:
# print(line, file=f)
-with open(ALIAS_hg38,"w") as f:
+# with open(ALIAS_hg38,"w") as f:
+# for line in alias:
+# print(line, file=f)
+
+with open(ALIAS_mm39,"w") as f:
for line in alias:
print(line, file=f)