Update README.md

Author: e-sollier

README.md

@@ -1,15 +1,21 @@
 # pyjacker
 
-This is a tool to detect enhancer hijacking events in a cohort profiled with WGS and RNA-seq. It does not require matched normals and can detect enhancer hijacking events occurring in a single sample. Briefly, it looks for outlier high and monoallelic expression of a gene in a sample which has a breakpoint close to the gene.
+This is a tool to detect enhancer hijacking events in a cohort of at least 10 samples profiled with WGS and RNA-seq. It does not require matched normals and can detect enhancer hijacking events occurring in only a single sample. Briefly, it looks for outlier high and monoallelic expression of a gene in a sample which has a breakpoint close to the gene.
 
 ## Usage
 In an environment with python>=3.7: 
 ```
 pip install pyjacker
-python pyjacker.py config.yaml
+pyjacker config.yaml
 ```
 The config file indicates all the parameters, including paths to the input files (see config_AML.yaml as an example). Alternatively, we provide a nextflow workflow that generates pyjacker's inputs from bam files, and run pyjacker: https://github.com/CompEpigen/wf_WGS.
 
+### Docker
+Alternatively, pyjacker can be run using a docker image.
+```
+docker run -t -v `pwd`:`pwd` esollier/pyjacker:latest pyjacker config.yaml
+```
+
 ## Inputs
 ### Gene expression table (required)
 Rows are genes (ensembl IDs) and columns are samples. The expression data must be provided in TPM. See data/TPM_ckAML.tsv for an example.
@@ -18,21 +24,35 @@ Rows are genes (ensembl IDs) and columns are samples. The expression data must b
 tsv file with columns: sample, chr1, pos1, chr2, pos2. 
 The fields chr2 and pos2 are optional (for example if you only have copy number data). See data/breakpoints.tsv for an example.
 
+### Reference files: gtf and cytobands (required)
+gtf file containing gene coordinates for your reference genome (see [data/Homo_sapiens.GRCh37.75.gtf.gz](https://github.com/CompEpigen/pyjacker/blob/main/data/Homo_sapiens.GRCh37.75.gtf.gz) or [data/Homo_sapiens.GRCh38.113.gtf.gz](https://github.com/CompEpigen/pyjacker/blob/main/data/Homo_sapiens.GRCh38.113.gtf.gz), and tsv file containing cytobands (see [data/cytobands_hg19](https://github.com/CompEpigen/pyjacker/blob/main/data/cytobands_hg19.tsv) or [data/cytobands_hg38](https://github.com/CompEpigen/pyjacker/blob/main/data/cytobands_hg38.tsv)).
+
 ### TADs (optional)
-A bed file of topologically-associating domains can be used, in which case only the breakpoints in the same TAD as a gene are considered in the search for enhancer hijacking events. One TAD file is provided in the data directory. If not TAD file is provided, pyjacker will instead look for breakpoints within a fixed distance to the gene (1.5Mb by default). See data/HSPC_TADs.bed for an example.
+A bed file of topologically-associating domains can be used, in which case only the breakpoints in the same TAD as a gene are considered in the search for enhancer hijacking events. See  [data/TADs_Dixon_IMR90_hg19.bed](https://github.com/CompEpigen/pyjacker/blob/main/data/TADs_Dixon_IMR90_hg19.bed) or [data/TADs_Dixon_IMR90_hg38.bed](https://github.com/CompEpigen/pyjacker/blob/main/data/TADs_Dixon_IMR90_hg38.bed) for TADs derived from the data of [Dixon et al.](https://www.nature.com/articles/nature14222) or [data/TADs_HSPC_hg19.bed](https://github.com/CompEpigen/pyjacker/blob/main/data/TADs_HSPC_hg19.bed) for TADs derived from [HSPCs](https://doi.org/10.1182/bloodadvances.2023012161). If not TAD file is provided, pyjacker will instead look for breakpoints within a fixed distance to the gene (1.5Mb by default).
 
 ### Allelic read counts at SNPs in RNAseq (optional)
-This is used to detect monoallelic expression. This requires files generated by [fast_ase](https://github.com/e-sollier/fast_ase) or [GATK ASEReadCounter](https://gatk.broadinstitute.org/hc/en-us/articles/360037428291-ASEReadCounter).
+This is used to detect monoallelic expression. This requires files generated by [fast_ase](https://github.com/e-sollier/fast_ase) or [GATK ASEReadCounter](https://gatk.broadinstitute.org/hc/en-us/articles/360037428291-ASEReadCounter). See [data/ASE_ckAML](https://github.com/CompEpigen/pyjacker/tree/main/data/ASE_ckAML) for example files.
 
 ### Copy number alterations (optional)
-tsv file with the following columns: sample, chr, start, end, cn. See data/CNAs.tsv for an example.
+tsv file with the following columns: sample, chr, start, end, cn. See [data/CNAs_ckAML.tsv](https://github.com/CompEpigen/pyjacker/blob/main/data/CNAs_ckAML.tsv) for an example.
 If provided, this will be used to:
 - correct gene expression based on copy number (so high expression because of amplification will not be reported)
 - filter out SNPs within deletions from the monoallelic expression detection
 
+### Enhancers (optional)
+A file of scored enhancers, generated by ROSE. See [data/enhancers_myeloid_hg19.tsv](https://github.com/CompEpigen/pyjacker/blob/main/data/enhancers_myeloid_hg19.tsv) for an example.
+
+### Fusion transcripts (optional)
+
+Fusion transcripts can also lead to aberrant high and monoallelic expression of a gene. If a list of fusion transcripts detected from RNAseq is provided, they will be used to annotate candidate enhancer hijacking events which are actually due to a fusion. See [data/fusions_ckAML.tsv](https://github.com/CompEpigen/pyjacker/blob/main/data/fusions_ckAML.tsv) for an example file.
+
+## Runtime
+Pyjacker takes approximately 5h to run on the ckAML dataset (39 samples) with default settings and 6 cores. The runtime is essentially proportional to the number of samples in the dataset and to the number of iterations used when estimating the null distribution of scores (used to compute the false discovery rate). This number of iterations is 50 by default, which ensures that accurate p-values are computed, but this can easily be reduced to 5-10 to reduce the runtime, without drastically altering the results.
+
+
+## Citation
+If you use pyjacker in your research, please consider citing:
+
+Sollier E, Riedel A, Toprak UH, Wierzbinska JA, Weichenhan D, Schmid JP, Hakobyan M, Touzart A, Jahn E, Vick B, Brown-Burke F. Pyjacker identifies enhancer hijacking events in acute myeloid leukemia including MNX1 activation via deletion 7q. bioRxiv. 2024:2024-09. https://doi.org/10.1101/2024.09.11.611224
 
-### Fusion transcripts
-Fusion transcripts can also lead to aberrant high and monoallelic expression of a gene. If a list of fusion transcripts detected from RNAseq is provided, they will be used to annotate candidate enhancer hijacking events which are actually due to a fusion. See data/fusions.tsv for an example file.
 
-### Enhancers
-A file of scored enhancers, generated by ROSE. See data/enhancers_myeloid.tsv for an example.