Commit 76627f6b authored by aditya.bhagwat's avatar aditya.bhagwat
Browse files

Fix double counting of pam mismatches

parent 40100b33
......@@ -161,6 +161,10 @@ read_bowtie_results <- function(outfile){
outfile,
col.names = c( 'readname', 'strand', 'target', 'position',
'readseq', 'quality', 'matches', 'mismatches'))
pattern <- '20:[ACGT][>][ACGT]'
dt %<>% extract(!stri_detect_regex(mismatches, pattern))
dt[ is.na(mismatches), mismatch := 0]
dt[!is.na(mismatches), mismatch := stringi::stri_count_fixed(
mismatches, '>')]
......@@ -233,6 +237,8 @@ match_seqs <- function(seqs, indexdir, norc, mismatches = 2,
readdt <- data.table(readname = names(reads),
readseq = unname(as.character(reads)))
readdt %<>% merge(matches, by='readname', all=TRUE, sort=FALSE)
readdt <- cbind(readdt[, c(1, 2)],
data.table::setnafill(readdt[, -c(1, 2)], fill = 0))
# Return
readdt[, 'readname' := NULL]
......
......@@ -10,22 +10,22 @@ reticulate::use_condaenv('azienv')
bsgenome <- BSgenome.Mmusculus.UCSC.mm10::BSgenome.Mmusculus.UCSC.mm10
bedfile <- system.file('extdata/SRF.bed', package='multicrispr')
targets <- multicrispr::bed_to_granges(bedfile, genome='mm10', plot = FALSE)
png('graphs/srf_karyogram.png')
multicrispr::plot_karyogram(targets, title = NULL)
dev.off()
#png('graphs/srf_karyogram.png')
#multicrispr::plot_karyogram(targets, title = NULL)
#dev.off()
extended <- extend(targets, -22, +22)
targets['T0151'] %>% plot_intervals()
targets['T0151'] %>% extend() %>% plot_intervals()
targets['T0151'] %>% up_flank() %>% plot_intervals()
targets['T0151'] %>% down_flank() %>% plot_intervals()
# targets['T0151'] %>% plot_intervals()
# targets['T0151'] %>% extend() %>% plot_intervals()
# targets['T0151'] %>% up_flank() %>% plot_intervals()
# targets['T0151'] %>% down_flank() %>% plot_intervals()
spacers <- extended %>% find_spacers(bsgenome)
spacers %<>% add_specificity(extended, bsgenome)
spacers %<>% add_efficiency(bsgenome, method = 'Doench2016')
spacers <- extended %>%
find_spacers(bsgenome) %>%
add_specificity(extended, bsgenome) %>%
add_efficiency(bsgenome, method = 'Doench2016')
spacers %>% subset(seqnames == 'chr1') %>%
subset(specific == TRUE) %>%
gr2dt() %>%
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment