Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

kpPlotRainfall error when given a path to VCF file that contains other variants than SNVs #153

Open
fawaz-dabbaghieh opened this issue Jun 14, 2024 · 0 comments

Comments

@fawaz-dabbaghieh
Copy link

fawaz-dabbaghieh commented Jun 14, 2024

Hi,

Thanks for the great library!

I am trying to generate some rainfall plots from somatic variants over chromosomes. I tested my code on a table generated from somatic SNVs and that works fine. However, I later noticed that kpPlotRainfall also takes a path to VCF as data, I tried that with VCF which contains more complex variants like indels, duplications, and BNDs, and I got the following error:

Error in `[[<-`(`*tmp*`, name, value = c(`NA` = "#888888")) : 
  1 elements in value to replace 0 elements
Calls: kpPlotRainfall -> $<- -> $<- -> [[<- -> [[<-
Execution halted

I am using custom genome and tracks (chm13), and generating the plots for each chromosome on a PDF page. This works fine when the data given to kpPlotRainfall is a GRanges object, but not when giving a VCF

# genome size and bands for chm13 t2t
sample = "sample1"
custom.genome <- toGRanges(genome_size)
custom.cytobands <- toGRanges(genome_bands)
chromosomes = c("chr1")
i <- 1

pdf("test.pdf", height=10, width=16)
# this is supposed to be part of a loop that plots for each chromosome
kp <- plotKaryotype(plot.type=1, genome = custom.genome[seqnames(custom.genome)==chromosomes[i]]
                    ,cytobands = custom.cytobands[seqnames(custom.cytobands)==chromosomes[i]],
                    chromosome = chromosomes[i])
kpAddCytobandLabels(kp, force.all = TRUE, srt = 90, col = 'orange')
title = paste0("[vcf] Somatic Mutations of ", sample, " and chromosome ", chromosomes[i])
kpAddMainTitle(kp, main=title, cex=1.5)
kpAddBaseNumbers(kp)
# only chromosome one VCF
kpPlotRainfall(kp, data = "sniffles_svs_chr1.vcf")
kpAxis(kp, ymax = 8, r0=0.01, r1=0.9)
kpAddLabels(kp, labels = c("Distance between mutations (log10)"), srt=90, pos=1, label.margin = 0.04,
            r0=0.01, r1=0.9)
dev.off()

To debug this, I kept removing different SVs (e.g. DEL, INS, INV, DUP, BND), until it worked with only SNVs left.

Is this the intended usage? Or am I doing something wrong?

Thanks for the help,
Fawaz

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

1 participant