CRISPR package demo
PGFE Umass Medical School
September 10 2014
========================================================
Please login to http://rstudio.umassmed.edu with your email username and password 

First load the required packages and specify the input file path.
We are going to use a sequence from human as input, which has been included as 
as fasta file in the CRISPRseek package.
To perform off-target analysis, we need to load Human BSgenome package
To annotate the target and off-targets, we need to load Human Transcript and gene identifier mapping packages.
In addition, you need to specify the output directory which will be the
directory to look for all the output files.
For the current release, you also need to specify the file containing all restriction enzyme (RE) cut patterns. You have the option to use the RE pattern file in the CRISPR package or specify your own RE pattern file.
In the future releases, you do not have to specify the RE pattern file unless you want to. For your convenience, the development version has been installed in rstudio.umassmed.edu. 


```{r}
library(CRISPRseek)
library(BSgenome.Hsapiens.UCSC.hg19)
library(TxDb.Hsapiens.UCSC.hg19.knownGene)
library(org.Hs.eg.db)

outputDir <- file.path(getwd(),"CRISPRseekDemo")
inputFilePath <- system.file('extdata', 'inputseq.fa', package = 'CRISPRseek')
REpatternFile <- system.file('extdata', 'NEBenzymes.fa', package = 'CRISPRseek')

```

Here is the command to learn more about CRISPR package and how to use the two workflow functions to find gRNAs and perform off-target analysis with different use cases.
========================================================

```{r}
args(offTargetAnalysis)
args(compare2Sequences)
?offTargetAnalysis
?compare2Sequences
?CRISPRseek
browseVignettes('CRISPRseek')
```

Scenario 1: Find all gRNAs and perform off-target analysis for each gRNA
========================================================
Please note that chromToSearch is set to chrX here for speed purpose, usually you do not need to set it, by default it is set to all.
```{r}
offTargetAnalysis(inputFilePath,  
    REpatternFile = REpatternFile, 
    BSgenomeName = Hsapiens, 
    chromToSearch ="chrX", 
    txdb = TxDb.Hsapiens.UCSC.hg19.knownGene,  
    orgAnn = org.Hs.egSYMBOL,
    outputDir = outputDir,  
    overwrite = TRUE) 
```

CRISPR with paired nickases decreases off-target cleavage by requiring the 
independent binding of two separate gRNAs around a genomic region.

Scenario 2: Find gRNAs in paired configuration and perform off-target analysis for each gRNA
========================================================
```{r}

offTargetAnalysis(inputFilePath, 
    findPairedgRNAOnly = TRUE, 
    min.gap = 0, 
    max.gap = 20,
    BSgenomeName = Hsapiens, 
    chromToSearch ="chrX",  
    txdb = TxDb.Hsapiens.UCSC.hg19.knownGene, 
    orgAnn = org.Hs.egSYMBOL, 
    max.mismatch = 0, 
    outputDir = outputDir,  
    overwrite = TRUE) 


```

You can specify the criteria of off-target search by specifying max.mismatch, e.g., allowing up to 2 mismatches to be considered as potential off-targets, by default it is set to 4.

```{r}
offTargetAnalysis(inputFilePath,  
    findPairedgRNAOnly = TRUE,  min.gap = 0, max.gap = 20, 
    BSgenomeName = Hsapiens, 
    chromToSearch ="chrX", 
    txdb = TxDb.Hsapiens.UCSC.hg19.knownGene,  
    orgAnn = org.Hs.egSYMBOL, 
    max.mismatch = 2, 
    outputDir = outputDir,  
    overwrite = TRUE) 

```

Paired gRNAs in proper spacing and orientation give more specificity and gRNAs 
overlap with restriction enzyme cut sites facilitates cleavage monitoring. 
Calling the function offTargetAnalysis with findPairedgRNAOnly = TRUE and 
findgRNAsWithREcutOnly = TRUE results in searching, scoring and annotating gRNAs 
that are in paired configuration and at least one of the pairs overlap a 
restriction enzyme cut site. To be considered as a pair, gap between forward 
gRNA and the corresponding reverse gRNA needs to be (min.gap, max.gap) inclusive
and the reverse gRNA must sit before the forward gRNA. The default 
(min.gap, max.gap) is (0,20). In order for a gRNA to be considered overlap with restriction enzyme cut site, the enzyme cut pattern must overlap with one of the gRNA positions specified in overlap.gRNA.positions, default position 17 and 18. 

Scenario 3: Find paired gRNAs with restriction enzyme cut site(s) and perform
off-target analysis
========================================================
```{r}

offTargetAnalysis(inputFilePath, 
	findgRNAsWithREcutOnly = TRUE, 
        minREpatternSize = 6,
 	overlap.gRNA.positions = c(17, 18),
	findPairedgRNAOnly = TRUE, min.gap = 0, max.gap = 20, 
    	BSgenomeName = Hsapiens, 
	chromToSearch ="chrX",  
    	txdb = TxDb.Hsapiens.UCSC.hg19.knownGene,  
	orgAnn = org.Hs.egSYMBOL,
	max.mismatch = 0, 
    	outputDir = outputDir,  
	overwrite = TRUE)

```



Scenario 4: Find gRNAs with restriction enzyme cut site(s) and  perform off-target 
analysis
========================================================
```{r}
offTargetAnalysis(inputFilePath, 
    findgRNAsWithREcutOnly = TRUE, 
    minREpatternSize = 6, 
    overlap.gRNA.positions = c(17, 18), 
    findPairedgRNAOnly = FALSE,  
    BSgenomeName = Hsapiens, 
    chromToSearch ="chrX",
    max.mismatch = 0, 
    txdb = TxDb.Hsapiens.UCSC.hg19.knownGene,  
    orgAnn = org.Hs.egSYMBOL, 
    outputDir = outputDir,  
    overwrite = TRUE) 
```

Scenario 5: Target and off-target analysis for user specified gRNAs
========================================================
Calling the function offTargetAnalysis with findgRNAs = FALSE will skip the gRNA search step and go directly to off-target search, scoring and annotation for the input gRNAs. The input gRNAs will be annotated with restriction enzyme cut sites for users to review later. 
However, paired information will not be available.

```{r}

gRNAFilePath <- system.file('extdata', 'testHsap_GATA1_ex2_gRNA1.fa', 
    package = 'CRISPRseek')

offTargetAnalysis(inputFilePath = gRNAFilePath,
    findPairedgRNAOnly = FALSE, 
    findgRNAs = FALSE,
    BSgenomeName = Hsapiens, 
    chromToSearch = 'chrX',
    txdb = TxDb.Hsapiens.UCSC.hg19.knownGene,  
    orgAnn = org.Hs.egSYMBOL,
    max.mismatch = 2, 
    outputDir = outputDir, 
    overwrite = TRUE)

```

Scenario 6. Quick gRNA finding without off-target analysis
========================================================
Calling the function offTargetAnalysis with chromToSearch = "" results in quick 
gRNA search without performing off-target analysis. Parameters 
findgRNAsWithREcutOnly and findPairedgRNAOnly can be tuned to indicate 
whether searching for gRNAs overlap restriction enzyme cut sites, and
whether searching for gRNAs in paired configuration.

```{r}
offTargetAnalysis(inputFilePath,
    chromToSearch = "", 
    outputDir = outputDir, 
    overwrite = TRUE)
```

Scenario 7. Find potential gRNAs preferentially targeting one of two alleles
without running time-consuming off-target analysis on all possible gRNAs.
========================================================
Below is an example to search for gRNAs that target at least one of the 
alleles. Two files are provided containing sequences that differ by a single
nucleotide polymorphism (SNP). The results are saved in file 
scoresFor2InputSequences.xls in outputDir directory.

Hungtinton disease is caused by mutations in the HTT gene. Expansion of CAG 
repeats in one copy of HTT can result in adult onset neurodegeneration. 
Because HTT is an essential gene, nucleases cannot be used that inactivate both
alleles. Therefore, to identify nuclease target sites that are allele-specific,
we will try to search for sites that overlap a single nucleotide polymorphism 
(SNP), RS362331 is located in a coding exon of HTT. Two sequences that differ only
at the polymorphism site will be used as inputs for compare2sequences. 


```{r}
inputFile1Path <- system.file("extdata", "rs362331C.fa", package = "CRISPRseek")
inputFile2Path <- system.file("extdata", "rs362331T.fa", package = "CRISPRseek")
seqs <- compare2Sequences(inputFile1Path, inputFile2Path,
    outputDir = outputDir, 
    REpatternFile = REpatternFile,
    overwrite = TRUE)
```

Exercise 1
========================================================

To preferentially target one allele, select gRNAs that have the lowest 
score for the other allele. Selected gRNAs can then be examined for potential off-target cleavage as described in Scenario 5.

Exercise 2
========================================================
Identify gRNAs that target the following two input sequences equally well with 
minimized off-target cleavage

>MfSerpAEx2
GACGATGGCATCCTCCGTTCCCTGGGGCCTCCTGCTGCTGGCGGGGCTGTGCTGCCTGGCCCCCCGCTCCCTGGCCTCGAGTCCCCTGGGAGCCGCTGTCCAGGACACAGGTGCACCCCACCACGACCATGAGCACCATGAGGAGCCAGCCTGCCACAAGATTGCCCCGAACCTGGCCGACTTCGCCTTCAGCATGTACCGCCAGGTGGCGCATGGGTCCAACACCACCAACATCTTCTTCTCCCCCGTGAGCATCGCGACCGCCTTTGCGTTGCTTTCTCTGGGGGCCAAGGGTGACACTCACTCCGAGATCATGAAGGGCCTTAGGTTCAACCTCACTGAGAGAGCCGAGGGTGAGGTCCACCAAGGCTTCCAGCAACTTCTCCGCACCCTCAACCACCCAGACAACCAGCTGCAGCTGACCACTGGCAATGGTCTCTTCATCGCTGAGGGCATGAAGCTACTGGATAAGTTTTTGGAGGATGTCAAGAACCTGTACCACTCAGAAGCCTTCTCCACCAATTTCGGGGACACCGAAGCAGCCAAGAAACAGATCAACGATTATGTTGAGAAGGGAACCCAAGGGAAAATTGTGGATTTGGTCAAAGACCTTGACAAAGACACAGCTTTCGCTCTGGTGAATTACATTTTCTTTAAAG

>HsSerpAEx2
GACAATGCCGTCTTCTGTCTCGTGGGGCATCCTCCTGCTGGCAGGCCTGTGCTGCCTGGTCCCTGTCTCCCTGGCTGAGGATCCCCAGGGAGATGCTGCCCAGAAGACAGATACATCCCACCATGATCAGGATCACCCAACCTTCAACAAGATCACCCCCAACCTGGCTGAGTTCGCCTTCAGCCTATACCGCCAGCTGGCACACCAGTCCAACAGCACCAATATCTTCTTCTCCCCAGTGAGCATCGCTACAGCCTTTGCAATGCTCTCCCTGGGGACCAAGGCTGACACTCACGATGAAATCCTGGAGGGCCTGAATTTCAACCTCACGGAGATTCCGGAGGCTCAGATCCATGAAGGCTTCCAGGAACTCCTCCGTACCCTCAACCAGCCAGACAGCCAGCTCCAGCTGACCACCGGCAATGGCCTGTTCCTCAGCGAGGGCCTGAAGCTAGTGGATAAGTTTTTGGAGGATGTTAAAAAGTTGTACCACTCAGAAGCCTTCACTGTCAACTTCGGGGACACCGAAGAGGCCAAGAAACAGATCAACGATTACGTGGAGAAGGGTACTCAAGGGAAAATTGTGGATTTGGTCAAGGAGCTTGACAGAGACACAGTTTTTGCTCTGGTGAATTACATCTTCTTTAAAG

Exercise 3
========================================================
Constraint gRNA Sequence by setting gRNA.pattern to require or exclude specific
features within the target site.

3a. Synthesis of gRNAs in vivo from host U6 promoters is more efficient if the 
first base is guanine. To maximize the efficiency, what can we set gRNA.pattern? 

3b. Synthesis of gRNAs in vitro using T7 promoters is most efficient when the 
first two bases are GG. To maximize the efficiency, what can we set gRNA.pattern? 


Exercise 4
========================================================
In the examples we went through, we deliberately restricted to search off-targets
in chromosome X. If we are interested in genome-wide search, what needs to be changed and how?

Exercise 5
========================================================
Find gRNAs in a paired configeration with distance apart between 5 and 15 without
performing off-target analysis


Exercise 6
========================================================
Create a txdb object

Exercise 7
========================================================
It is known that different CRISPR-cas system uses different PAM sequence, what
parameter needs to be reset for PAM = 'NNNNGGGT'?

Exercise 8
========================================================
It is known that different CRISPR-cas system has different gRNA length, what
parameter needs to be reset?

Exercise 9
========================================================
Which parameter needs to be reset to what if we are interested in finding gRANs with
restriction enzyme pattern of size 8 or above?

Exercise 10
========================================================
New penalty matrix has been recently derived, which parameter needs to be reset
accordingly?

Exercise 11
========================================================
It has been shown that although PAM sequence NGG is preferred, a variant NAG is
also recognized with less effeciency. The researcher is interested in performing
off-target searching to include both NGG and NAG variants. What parameter(s) need to be set correctly to carry such a search?

Exercise 12
========================================================
Which parameter to reset if you would like to skip off-target analysis but still want the summary output file? How to reset the parameter?


Exercise 13
========================================================
Could you think of any other use cases?
