match_spacers.Rd
Count matches to indexed target/genome and add to GRanges
match_spacers( spacers, indexdir, norc, mismatches = 2, outdir = OUTDIR, pam = "NGG", verbose = TRUE )
spacers | spacer |
---|---|
indexdir | string: dir containing indexed reference.
This can be an indexed genome( |
norc | TRUE or FALSE: whether to run bowtie also with revcompls Generally TRUE for genome and FALSE for target matches, because target ranges generally include both strands. |
mismatches | number (default 2): max number of mismatches to consider |
outdir | string: file where to output bowtie results |
pam | string (default 'NGG') pam pattern to expand |
verbose | TRUE (default) or FALSE |
data.table
Expands iupac amgiguities in the pam sequence. Matches all resulting sequences against (indexes) target and genome. Adds match counts to GRanges object, and then returns it.
# PE example #----------- require(magrittr) bsgenome <- BSgenome.Hsapiens.UCSC.hg38::BSgenome.Hsapiens.UCSC.hg38 gr <- char_to_granges(c(PRNP = 'chr20:4699600:+', # snp HBB = 'chr11:5227002:-', # snp HEXA = 'chr15:72346580-72346583:-', # del CFTR = 'chr7:117559593-117559595:+'), # ins bsgenome) spacers <- find_pe_spacers(gr, bsgenome)# indexdir <- genome_dir(indexedgenomesdir = INDEXEDGENOMESDIR, bsgenome) # match_spacers(spacers, indexdir, norc=TRUE, mismatches = 1) # TFBS example #------------- bsgenome <- BSgenome.Mmusculus.UCSC.mm10::BSgenome.Mmusculus.UCSC.mm10 bedfile <- system.file('extdata/SRF.bed', package = 'multicrispr') targets <- extend(bed_to_granges(bedfile, genome = 'mm10'))#>#>#>#>#>#>#>#>#>match_spacers(spacers, indexdir, norc=FALSE, mismatches = 1)#>#>#>#>#>#> crisprspacer MM0 MM1 #> 1: ATATAAGGGCATTGGAAGAA 2 1946 #> 2: AATATAAGGGCATTGGAAGA 2 1920 #> 3: TGGAGACAATATAAGGGCAT 2 1950 #> 4: TTCTGCTGGAGACAATATAA 4 1654 #> 5: CTTCTGCTGGAGACAATATA 1650 484 #> --- #> 3013: AGATGAGGAATATGCAAATA 2 0 #> 3014: TTGCATATTCCTCATCTGAT 2 0 #> 3015: AGTGTGCTTATAAGGGGGGA 2 0 #> 3016: AGAGAGTGTGCTTATAAGGG 2 4 #> 3017: AGCACACTCTCTTAGTAAAT 2 6