filtergRNAs {CRISPRseek} | R Documentation |
Filter gRNAs containing restriction enzyme cut site
filtergRNAs(all.gRNAs, pairOutputFile = "", findgRNAsWithREcutOnly = FALSE, REpatternFile = system.file("extdata", "NEBenzymes.fa", package = "CRISPRseek"), format = "fasta", minREpatternSize = 4, overlap.gRNA.positions = c(17, 18),overlap.allpos = TRUE)
all.gRNAs |
gRNAs as DNAStringSet, such as the output from findgRNAs |
pairOutputFile |
File path with paired gRNAs |
findgRNAsWithREcutOnly |
Indicate whether to find gRNAs overlap with restriction enzyme recognition pattern |
REpatternFile |
File path containing restriction enzyme cut patters |
format |
Format of the REpatternFile, default as fasta |
minREpatternSize |
Minimum restriction enzyme recognition pattern length required for the enzyme pattern to be searched for, default 4 |
overlap.gRNA.positions |
The required overlap positions of gRNA and restriction enzyme cut site, default 17 and 18 |
overlap.allpos |
Default TRUE, meaning that only gRNAs overlap with all the positions are retained FALSE, meaning that gRNAs overlap with one or both of the positions are retained |
gRNAs.withRE |
gRNAs as DNAStringSet that passed the filter criteria |
gRNAREcutDetails |
a data frame that contains a set of gRNAs annotated with restriction enzyme cut details |
Lihua Julie Zhu
offTargetAnalysis
all.gRNAs <- findgRNAs( inputFilePath = system.file("extdata", "inputseq.fa", package = "CRISPRseek"), pairOutputFile = "testpairedgRNAs.xls", findPairedgRNAOnly = TRUE) gRNAs.RE <- filtergRNAs(all.gRNAs = all.gRNAs, pairOutputFile = "testpairedgRNAs.xls",minREpatternSize = 6, REpatternFile = system.file("extdata", "NEBenzymes.fa", package = "CRISPRseek"), overlap.allpos = TRUE) gRNAs <- gRNAs.RE$gRNAs.withRE restriction.enzyme.cut.sites <- gRNAs.RE$gRNAREcutDetails