Hello, recently I've found the mismatch issue with my RNAseq data.
Basically I have both normal & tumor RNAseq data. I did the variant calling following the pipeline below:
- STAR (index using GENCODE GRCh38.primary_assembly.genome.fa & gencode.v40.annotation.gtf)
- GATK (GRCh38.primary_assembly.genome.fa provided as fasta)
2a. MarkDuplicatesSpark
2b. SplitNCigarReads
2c. BQSR
2d. Mutect2
- VEP (version 106)
vep --input_file vcf_files/${sample}_RNA_match_RNA_somatic.filtered.vcf.gz \
--format vcf \
--output_file vep/${sample}_RNA_match_RNA_somatic.filtered.vep.vcf \
--vcf --symbol --terms SO --tsl --hgvs \
--fasta reference/vep/homo_sapiens/106_GRCh38/Homo_sapiens.GRCh38.dna.toplevel.fa.gz \
--offline --cache --dir_cache reference/vep/ \
--plugin Frameshift --plugin Wildtype --dir_plugins reference/vep/plugins/ \
--synonyms reference/vep/homo_sapiens/106_GRCh38/chr_synonyms.txt
- pVACtools (version 3.0.1)
pvacseq run --iedb-install-directory ~/.conda/miniconda/4.9/envs/neo/ \
--binding-threshold 500 --n-threads 24 \
--net-chop-method "cterm" --netmhc-stab --exclude-NAs \
--normal-sample-name ${sample}_N \
vep/${sample}_RNA_match_RNA_somatic.filtered.vep.vcf \
${sample} ${hla1},${hla2} \
MHCflurry MHCnuggetsI NetMHC NetMHCpan \
pvacseq/${sample}_RNA_match_RNA
The first issue is when running VEP, I got the follow message for a few lines, but it still generated a VEP annotated vcf.
Use of uninitialized value in numeric lt (<) at /hpctmp/renyi04/reference/vep/plugins/Frameshift.pm line 129, <$fh> line 5840.
Use of uninitialized value in numeric lt (<) at /hpctmp/renyi04/reference/vep/plugins/Frameshift.pm line 129, <$fh> line 5840.
Use of uninitialized value in addition (+) at /hpctmp/renyi04/reference/vep/plugins/Frameshift.pm line 129, <$fh> line 5840.
The second error is when I run pvacseq, I got warnings and errors below:
KI270733.1 135027 C CT
Variant doesn't have protein position information. Skipping.
KI270733.1 135027 C CCG
Variant doesn't have protein position information. Skipping.
KI270733.1 135028 CGCGCGCGCGGGAGGGCGCGTGCCCCG C
Variant doesn't have protein position information. Skipping.
KI270733.1 168056 T G ENST00000615130
Variant doesn't have protein position information. Skipping.
KI270733.1 177433 G GGGC ENST00000618998
Completed
Generating Variant Peptide FASTA and Key File
Warning: Amino acid change is not sane - contains multiple stops. Skipping entry 2208.EIF4E2.ENST00000258416.inframe_del.176-177QVSHP*ARLVSCVAFALLAAPDLL**FLLFLW*E/QE
Warning: Amino acid change is not sane - contains multiple stops. Skipping entry 2209.EIF4E2.ENST00000409098.inframe_del.176-177QVSHP*ARLVSCVAFALLAAPDLL**FLLFLW*E/QE
Warning: Amino acid change is not sane - contains multiple stops. Skipping entry 2210.EIF4E2.ENST00000409167.inframe_del.131-132QVSHP*ARLVSCVAFALLAAPDLL**FLLFLW*E/QE
Warning: Amino acid change is not sane - contains multiple stops. Skipping entry 2211.EIF4E2.ENST00000409322.inframe_del.131-132QVSHP*ARLVSCVAFALLAAPDLL**FLLFLW*E/QE
Warning: Amino acid change is not sane - contains multiple stops. Skipping entry 2212.EIF4E2.ENST00000409394.inframe_del.131-132QVSHP*ARLVSCVAFALLAAPDLL**FLLFLW*E/QE
Warning: Amino acid change is not sane - contains multiple stops. Skipping entry 2213.EIF4E2.ENST00000409495.inframe_del.176-177QVSHP*ARLVSCVAFALLAAPDLL**FLLFLW*E/QE
Warning: Amino acid change is not sane - contains multiple stops. Skipping entry 2214.EIF4E2.ENST00000409514.inframe_del.176-177QVSHP*ARLVSCVAFALLAAPDLL**FLLFLW*E/QE
Warning: Amino acid change is not sane - contains multiple stops. Skipping entry 2215.EIF4E2.ENST00000454501.inframe_del.171-172QVSHP*ARLVSCVAFALLAAPDLL**FLLFLW*E/QE
Warning: Amino acid change is not sane - contains multiple stops. Skipping entry 2216.EIF4E2.ENST00000687222.inframe_del.131-132QVSHP*ARLVSCVAFALLAAPDLL**FLLFLW*E/QE
Warning: Amino acid change is not sane - contains multiple stops. Skipping entry 2217.EIF4E2.ENST00000690794.inframe_del.76-77QVSHP*ARLVSCVAFALLAAPDLL**FLLFLW*E/QE
ERROR: There was a mismatch between the actual wildtype amino acid sequence (VHILGEPRPHLFGQMFVRLQLLRAVREVLHTGLAMLGL) and the expected amino acid sequence (VHILGVSTQQNGVGRAERYRVKAKEAEDETSRPFLFQE). Did you use the same reference build version for VEP that you used for creating the VCF?
{'chromosome_name': 'chr3', 'start': '49015706', 'stop': '49015817', 'reference': 'TCCTGAAAGAGAAAGGGCCTGCTGGTCTCATCCTCTGCTTCCTTTGCCTTTACCCTATACCTCTCTGCACGTCCCACCCCATTTTGCTGTGTGCTCACCCCCAGGATGTGTA', 'variant': 'T', 'gene_name': 'DALRD3', 'transcript_name': 'ENST00000313778', 'transcript_support_level': '2', 'amino_acid_change': 'VHILGVSTQQNGVGRAERYRVKAKEAEDETSRPFLFQE/E', 'codon_change': 'gTACACATCCTGGGGGTGAGCACACAGCAAAATGGGGTGGGACGTGCAGAGAGGTATAGGGTAAAGGCAAAGGAAGCAGAGGATGAGACCAGCAGGCCCTTTCTCTTTCAGGag/gag', 'ensembl_gene_id': 'ENSG00000178149', 'hgvsc': 'ENST00000313778.9:c.998_1012del', 'hgvsp': 'ENSP00000323265.5:p.Val333_Gln337del', 'wildtype_amino_acid_sequence': 'MLTFLQQLRVDWPAASERASSHTLRSHALEELTSANDGRTLSPGILGRLCLKELVEEQGRTAGYDPNLDNCLVTEDLLSVLAELQEALWHWPEDSHPGLAGASDTGTGGCLVVHVVSCEEEFQQQKLDLLWQKLVDKAPLRQKHLICGPVKVAGAPGTLMTAPEYYEFRHTQVCKASALKHGGDLAQDPAWTEIFGVLSVATIKFEMLSTAPQSQLFLALADSSISTKGTKSGTFVMYNCARLATLFESYKCSMEQGLYPTFPPVSSLDFSLLHDEGEWLLLFNSILPFPDLLSRTAVLDCTAPGLHIAVRTEMICKFLVQLSMDFSSYYNRVHILGEPRPHLFGQMFVRLQLLRAVREVLHTGLAMLGLPPLSHI', 'frameshift_amino_acid_sequence': '', 'fusion_amino_acid_sequence': '', 'variant_type': 'inframe_del', 'protein_position': '333-338', 'transcript_expression': 'NA', 'gene_expression': 'NA', 'normal_depth': 'NA', 'normal_vaf': 'NA', 'tdna_depth': '473', 'tdna_vaf': '0.016', 'trna_depth': 'NA', 'tr
[CV023_RNA_match_RNA_somatic.filtered.vep.vcf.zip](https://github.com/griffithlab/pVACtools/files/9018016/CV023_RNA_match_RNA_somatic.filtered.vep.vcf.zip)
na_vaf': 'NA', 'index': '2823.DALRD3.ENST00000313778.inframe_del.333-338VHILGVSTQQNGVGRAERYRVKAKEAEDETSRPFLFQE/E', 'protein_length_change': ''}
I tried to use VAtools (5.0.1) ref-transcript-mismatch-reporter to exclude the mismatch. But got the error after it generated some entries
Traceback (most recent call last):
File "/home/svu/renyi04/.conda/miniconda/4.9/envs/neo/bin/ref-transcript-mismatch-reporter", line 8, in <module>
sys.exit(main())
File "/home/svu/renyi04/.conda/miniconda/4.9/envs/neo/lib/python3.8/site-packages/vatools/ref_transcript_mismatch_reporter.py", line 162, in main
position = int(protein_position.split('-', 1)[0]) - 1
ValueError: invalid literal for int() with base 10: ''
Last, after I removed the mismatched entries mannually (about 4-7 in different samples), pVACtools worked.
I saw the same issue about VAtools but it should be solved in the lastest version. Any suggestion how I can solve the issue?
My vep.vcf is attached
CV023_RNA_match_RNA_somatic.filtered.vep.vcf.zip