Identification of genetic variation with bulk tissue RNA-seq

RNA was collected using Qiagen RNeasy plus kit (Catalog ID: 74034, Qiagen, Hilden, Germany). RNA-seq libraries were prepared following the Takara Pico v3 low input protocol and sequenced on Illumina NextSeq500.

Raw RNA-seq data were trimmed for polyA sequences and low-quality bases using cutadapt (v2.4)⁵¹. Reads were aligned to human genome hg38 using STAR(v 2.7.2b)⁵². Duplicate read identification and other quality control checks for read alignment were performed using CollectRNASeqMetrics and MarkDuplicates in Picard Tools.⁵³ Reads containing N were split using SplitNCigarReads function in the Genome Analysis Toolkit (GATK)^54,55. Bases quality scores were recalibrated using known variants from the GATK resource bundle and with the BaseRecalibrator and ApplyBQSR functions in GATK^54,55. Somatic SNV and indels were called with Mutect2 in tumor-only mode^54,55. Only variants with at least read depth of 10, 5% allele frequency, read depth of 5 for the alternate allele were kept for analysis. The variants were then filtered for variants in sex or mitochondrial chromosomes, RNA editing sites, repeat masker regions, and variants in Panel of Normal (from GATK) references. Variants were then annotated using the Funcotator function in GATK^54,55.