Hello everyone,
when I try the last step of purecn using cnvkit outputs:
$RSCRIPT_PATH $PURECN/PureCN.R \
--out $CNV_T_RESULTS/${TUMOR_PFX} \
--sampleid ${TUMOR_PFX}\
--tumor $CNV_T_RESULTS/${TUMOR_PFX}_T.cnr \
--segfile $CNV_T_RESULTS/${TUMOR_PFX}.seg \
--vcf $CNV_T_RESULTS/${TUMOR_PFX}.vcf.gz \
--genome hg19 \
--force --postoptimize --seed 123
I get this error:
> log2-ratio in tumor coverage data. WARN [2020-11-10 23:04:51] Allosome
> coverage missing, cannot determine sex. WARN [2020-11-10 23:04:51]
> Allosome coverage missing, cannot determine sex. INFO [2020-11-10
> 23:04:51] Using 8016 intervals (8016 on-target, 0 off-target). INFO
> [2020-11-10 23:04:51] No off-target intervals. If this is
> hybrid-capture data, consider adding them. INFO [2020-11-10 23:04:51]
> Loading VCF... INFO [2020-11-10 23:04:51] Found 837 variants in VCF
> file. INFO [2020-11-10 23:04:51] Removing 150 triallelic sites. WARN
> [2020-11-10 23:04:52] vcf.file has no DB info field for membership in
> germline databases. Found and used valid population allele frequency >
> 0.001000 instead. INFO [2020-11-10 23:04:52] 447 (65.1%) variants annotated as likely germline (DB INFO flag). INFO [2020-11-10
> 23:04:53] AR2T.060A is tumor in VCF file. INFO [2020-11-10 23:04:53] 0
> homozygous and 4 heterozygous variants on chrX. INFO [2020-11-10
> 23:04:53] Sex from VCF: F (Fisher's p-value: 0.32, odds-ratio: 0.00).
> INFO [2020-11-10 23:04:53] Detected MuTect2 VCF. INFO [2020-11-10
> 23:04:53] Removing 0 MuTect2 calls due to blacklisted failure reasons.
> INFO [2020-11-10 23:04:54] Initial testing for significant sample
> cross-contamination: unlikely INFO [2020-11-10 23:04:55] Removing 332
> variants with AF < 0.030 or AF >= 0.970 or less than 4 supporting
> reads or depth < 15. Error in which(as.numeric(geno(vcf)$MBQ[,
> tumor.id.in.vcf]) >= min.base.quality) : (list) object cannot be
> coerced to type 'double' Calls: runAbsoluteCN ... filterVcfMuTect2 ->
> filterVcfBasic -> .filterVcfByBQ -> which Execution halted
I used Mutect2, GATK4 and purecn 1.12.2. Can anyone advice me on what to do? Thank you
Hi Markus,
Thank you for your reply. I do not have the option of upgrading my current R version which is 3.5.1 and I can only install pureCN 1.12.2 on this version of R. Is it possible to work around this? I am also not clear on whether my workflow is correct. I would like to run PureCN with Mutect2, GATK4 and integrate cnvkit results with PureCN.
Currently what I do is I create a pooled reference by: 1) running each germline (unmatched but process matched) sample in Mutect2 tumor mode, and then run SomaticPanelOFNormals on this.
2) using bcftools merge to merge all the Mutect2 vcf files for the germlines to create one merged germline.vcf. I then compress using bgzip and index using tabix.
3) I run this step:
4) I create a mapping bias file using an R script
contents of create_mapping_bias_file.r:
5) I run this command for every tumor sample:
5) I run mutect2 on every single tumor sample
6) I run pureCN
In the last step, somehow I get the error that the
and if I remove that, I get the error above ((list) object cannot be
I have searched extensively for documentation on how to integrate Mutect2 and GATK4 with PureCN while using cnvkit output and I could not get much information and so I thought documenting it here and asking for advice might help others in future. Thank you.
You can install the latest version from GitHub. If BiocManager is not available for your R version, try devtools::install_github. PureCN 1.20.0 is the only version that works with GATK4. I would also recommend updating GATK4.
Hi Markus,
Thank you for your reply. I am now running PureCN 1.21.0 but still see the same error in the last step:
Could you advice on whether I might be making a mistake anywhere in the workflow I am following? Thank you!
Hi, I don't know if I'm correct here.
One suggestion: At step3 you have mentioned above NormalDB.R seems to be internally calling the calculateMappingBiasVcf() to create a mapping bias file. So, I think you don't need to perform the step4 using create_mapping_bias_file.r again.
Also, according to the vignette, --normaldb and --mappingbiasfile are two different argument and needs different inputs. But as per the step6 you have showed above, you seems to be using the wrong argument and input combination like,
--normaldb mapping_bias_agilent_v6_hg19.rds
, but I guess it supposed to be--mappingbiasfile mapping_bias_agilent_v6_hg19.rds
Actually, I just started working on this tool, not sure if I'm correct. Hope this helps.