salmon
salmon copied to clipboard
COMBINE-lab
Feature Request: Genome Mappings on --writeMappings
Is it possible to implement a utility function within the Salmon binary that converts the transcript mappings to genome coordinates of an input GTF file?
Hi @AndrewSkelton ,
Yes, this is definitely possible. I'll put it on the short-list of feature requests. We're also happy to accept any PRs if you know someone who might be interested in implementing this.
This seems like a fairly simple job for an external tool as well. What did you have in mind? Just transcript --> genomic coordinates in a table or are you interested in visualizing the pseudobam output mapped to genomic coordinates?
Hey @rob-p @mdshw5 ,
First, thanks a lot for implementing salmon, it is super helpful.
Has this been achieved? I would like to visualise the salmon quantoutput in IGV.
I have generate the .bam file with salmon quant -lISR -r Parent_NGSC3_DI_HodgkinsLymphoma_S1_L002_R2_001.fastq.gz -i $SALMON_INDEX -p 26 -o alignment_file --validateMappings --no-version-check --writeMappings | samtools view -Sb - | samtools sort -T sort.tmp -o - > salmon.bam
But if I explore the bam file with samtools viewwhat I see is:
A00519:603:H7KMNDSXY:2:1307:5900:21699 16 ENST00000390436.2 204 255 91M * 0 0 AAAAGACCTCAGCTTATTATAGACATTCGTTCAAATGTGGGCGAAAAGAAAGACCAACGAATTGCTGTTACATTGAACAAGACAGCCAAAC * NH:i:1 HI:i:1 XT:A:T AS:i:176
A00519:603:H7KMNDSXY:2:1307:5801:2018 0 ENST00000535880.2 326 255 91M * 0 0 GACGGTTTTCTGTGAAACACATTCTGACCCAGAAAGCCTTTCACTTGGTGATCTCTCCAGTAAGGACTGAAGACAGTGCCACTTACTACTG * NH:i:1 HI:i:1 XT:A:T AS:i:182
Where in the 3rd column I am seeing the transcript ID. I would like to have the chromosome coordinates of that read, like here:
A00519:603:H7KMNDSXY:2:1401:1814:22357 0 chr1 11261 0 90M * 0 0 CTATTGCTTAGACTGGTGGCCAGCGCCCCCTGCTGGCGCCGGGGCACTGCAGGGCCCTCTTGCTTACTGTATAGTGGTGGCACGCCGCCT FF,FFFFFFFFF,FFFFFFFFFFFFF:FFFF,FFFFFFFFFFFF,FFFFFFF:,FFFFFFFFFFFF,:,:FFFFF,FFFF::F,F,F:FF NH:i:5 HI:i:1 AS:i:84 nM:i:2 RG:Z:Parent_NGSC3_DI_HodgkinsLymphoma:0:1:H7KMNDSXY:2 RE:A:I xf:i:0 CR:Z:GATTCGACAACGCATT CY:Z:FFFFFFFFFF:F,F:F CB:Z:GATTCGACAACGCATT-1 UR:Z:CGATCACGGAAC UY:Z:FFFFFFF,:FFF UB:Z:CGATCACGGAAC
A00519:603:H7KMNDSXY:2:1250:26503:19272 16 chr1 11570 1 90M * 0 0 GATTACCATCAGAATTGTACTGTTCTGTATCCCACCAGCAATGTCTAGGAATGCCTGTTTCTCCACAAAGTGTTTACTTTTGGATTTTTG FFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFF NH:i:3 HI:i:1 AS:i:88 nM:i:0 RG:Z:Parent_NGSC3_DI_HodgkinsLymphoma:0:1:H7KMNDSXY:2 RE:A:I xf:i:0 CR:Z:CCGGACAGTACAAACA CY:Z:FFFFFFFFFFFFFFFF CB:Z:CCGGACAGTACAAACA-1 UR:Z:AAGATTCACTAT UY:Z:FFFFFFFFFFFF UB:Z:AAGATTCACTAT
Is that possible? It would be awesome to have the genomic coordinates.
Thanks in advance and sorry for the inconveniences,
Best,
Kike
@kikegoni this is something I started to work on but then quickly abandoned. I've uploaded the utility (transcoorder) that I started writing (https://github.com/mdshw5/transcoorder) for this purpose. If you find it useful and want to help finish the work please do. Currently the transcoord command will take a SAM/BAM, a GTF and matching genomic FASTA, and will convert reads from transcript to genomic coordinates, with appropriate reference names and offsets. However, it is not tested, and will only properly handle reads that fall entirely within an exon. Spliced reads shouldn't be too hard to add, but I just don't have the time right now.
$ transcoord -h
usage: transcoord [-h] [-o OUT] [-t TAG_NAME] [--debug] [--version] gtf bam fasta
positional arguments:
gtf GTF file containing transcripts
bam SAM or BAM files aligned to transcriptome
fasta FASTA format assembly coresponding to GTF
optional arguments:
-h, --help show this help message and exit
-o OUT, --out OUT output file for genomic SAM (default: stdout)
-t TAG_NAME, --tag-name TAG_NAME
SAM tag name for storing transcript identifier. default: ZT
--debug enable debugging
--version display version number
The command is veeeeeery slow, but should process a typical RNAseq library in a few hours.
The example output for SRR5024127 is:
Transcript coordinate SAM:
s/docsource/analysis/SRR5024127/salmon/gencode_basic_exonic/mapped.bam | head -n10
@HD VN:1.0 SO:coordinate
@SQ SN:ENST00000456328.2 LN:1657 DS:T
@SQ SN:ENST00000450305.2 LN:632 DS:T
@SQ SN:ENST00000488147.1 LN:1351 DS:T
@SQ SN:ENST00000619216.1 LN:68 DS:T
@SQ SN:ENST00000473358.1 LN:712 DS:T
@SQ SN:ENST00000469289.1 LN:535 DS:T
@SQ SN:ENST00000607096.1 LN:138 DS:T
@SQ SN:ENST00000417324.1 LN:1187 DS:T
@SQ SN:ENST00000606857.1 LN:840 DS:T
...
SRR5024127.14884291.CTGAAGCT 163 ENST00000456328.2 99 1 75M = 242 213 TTGAGAGGTCACAGGGTCTTGATGCTGTGGTCTTGATCTGCAGGTGTCTGACTTCCAGCAACTGCTGGCCTGTGC * NH:i:4 HI:i:1 XT:A:T AS:i:144
SRR5024127.17622738.CTGAAGCT 163 ENST00000456328.2 120 1 75M = 162 103 ATGCTGTGGTCTTGATCTGCAGGTGTCTGACTTCCAGCAACTGCTGGCCTGTGCCAGGGTGCAAGCTGAGCACTG * NH:i:4 HI:i:1 XT:A:T AS:i:144
SRR5024127.17649643.CTGAAGCT 163 ENST00000456328.2 120 1 75M = 162 103 ATGCTGTGGTCTTGATCTGCAGGTGTCTGACTTCCAGCAACTGCTGGCCTGTGCCAGGGTGCAAGCTGAGCACTG * NH:i:4 HI:i:1 XT:A:T AS:i:144
SRR5024127.19781953.CTGAAGCT 99 ENST00000456328.2 125 1 75M = 183 133 GTGGTCTTCATCTGCAGGTGTCTGACTTCCAGCAACTGCTGGCCTGTGCCAGGGTGCAAGCTGAGCACTGGAGTG * NH:i:4 HI:i:1 XT:A:T AS:i:150
SRR5024127.17622738.CTGAAGCT 83 ENST00000456328.2 162 1 75M = 120 -103 GCTGGCCTGTGCCAGGGTGCAAGCTGAGCACTGGAGTGGAGTTTTCCTGTGGAGAGGAGCCATGCCTAGAGTGGG * NH:i:4 HI:i:1 XT:A:T AS:i:150
SRR5024127.17649643.CTGAAGCT 83 ENST00000456328.2 162 1 75M = 120 -103 GCTGGCCTGTGCCAGGGTGCAAGCTGAGCACTGGAGTGGAGTTTTCCTGTGGAGAGGAGCCATGCCTAGAGTGGG * NH:i:4 HI:i:1 XT:A:T AS:i:150
SRR5024127.19781953.CTGAAGCT 147 ENST00000456328.2 183 1 75M = 125 -133 AGCTGAGCACTGGAGTGGAGTTTTCCTGTGGAGAGGAGCCATGCCTAGAGTGGGATGGGCCATTGTTCATCTTCT * NH:i:4 HI:i:1 XT:A:T AS:i:150
SRR5024127.14884291.CTGAAGCT 83 ENST00000456328.2 242 1 75M = 99 -213 CCATTGTTCATCTTCTGGCCCCTGTTGTCTGCATGTAACTTAATACCACAACCAGGCATAGGGGAAAGATCGGAG * NH:i:4 HI:i:1 XT:A:T AS:i:144
SRR5024127.12905971.CTGAAGCT 163 ENST00000456328.2 264 1 75M = 286 82 TGTTGTCTGCATGTAACTTAATACCACAACCAGGCATAGGGGAAAGATTGGAGGAAAGATGAGTGAGAGCATCAA * NH:i:4 HI:i:1 XT:A:T AS:i:150
SRR5024127.10923795.CTGAAGCT 137 ENST00000456328.2 286 1 75M = 286 0 ACCACAACCAGGCATAGGGGAAAGATTGGAGGAAAGATGAGTGAGAGCATCAACTTCTCTGACAACCTAGGCCAG * NH:i:4 HI:i:1 XT:A:T AS:i:144
Genomic coordinate SAM:
@HD VN:1.0 SO:unknown
@SQ SN:chr1 LN:248956422
@SQ SN:chr2 LN:242193529
@SQ SN:chr3 LN:198295559
@SQ SN:chr4 LN:190214555
@SQ SN:chr5 LN:181538259
@SQ SN:chr6 LN:170805979
@SQ SN:chr7 LN:159345973
@SQ SN:chr8 LN:145138636
@SQ SN:chr9 LN:138394717
@SQ SN:chr10 LN:133797422
...
SRR5024127.14884291.CTGAAGCT 163 chr1 11968 1 75M = 242 213 TTGAGAGGTCACAGGGTCTTGATGCTGTGGTCTTGATCTGCAGGTGTCTGACTTCCAGCAACTGCTGGCCTGTGC * AS:i:144 HI:i:1 NH:i:4 XT:Z:T ZT:Z:ENST00000456328.2
SRR5024127.17622738.CTGAAGCT 163 chr1 11989 1 75M = 162 103 ATGCTGTGGTCTTGATCTGCAGGTGTCTGACTTCCAGCAACTGCTGGCCTGTGCCAGGGTGCAAGCTGAGCACTG * AS:i:144 HI:i:1 NH:i:4 XT:Z:T ZT:Z:ENST00000456328.2
SRR5024127.17649643.CTGAAGCT 163 chr1 11989 1 75M = 162 103 ATGCTGTGGTCTTGATCTGCAGGTGTCTGACTTCCAGCAACTGCTGGCCTGTGCCAGGGTGCAAGCTGAGCACTG * AS:i:144 HI:i:1 NH:i:4 XT:Z:T ZT:Z:ENST00000456328.2
SRR5024127.19781953.CTGAAGCT 99 chr1 11994 1 75M = 183 133 GTGGTCTTCATCTGCAGGTGTCTGACTTCCAGCAACTGCTGGCCTGTGCCAGGGTGCAAGCTGAGCACTGGAGTG * AS:i:150 HI:i:1 NH:i:4 XT:Z:T ZT:Z:ENST00000456328.2
SRR5024127.17622738.CTGAAGCT 83 chr1 12031 1 75M = 120 -103 GCTGGCCTGTGCCAGGGTGCAAGCTGAGCACTGGAGTGGAGTTTTCCTGTGGAGAGGAGCCATGCCTAGAGTGGG * AS:i:150 HI:i:1 NH:i:4 XT:Z:T ZT:Z:ENST00000456328.2
SRR5024127.17649643.CTGAAGCT 83 chr1 12031 1 75M = 120 -103 GCTGGCCTGTGCCAGGGTGCAAGCTGAGCACTGGAGTGGAGTTTTCCTGTGGAGAGGAGCCATGCCTAGAGTGGG * AS:i:150 HI:i:1 NH:i:4 XT:Z:T ZT:Z:ENST00000456328.2
SRR5024127.19781953.CTGAAGCT 147 chr1 12052 1 75M = 125 -133 AGCTGAGCACTGGAGTGGAGTTTTCCTGTGGAGAGGAGCCATGCCTAGAGTGGGATGGGCCATTGTTCATCTTCT * AS:i:150 HI:i:1 NH:i:4 XT:Z:T ZT:Z:ENST00000456328.2
SRR5024127.14884291.CTGAAGCT 83 chr1 12111 1 75M = 99 -213 CCATTGTTCATCTTCTGGCCCCTGTTGTCTGCATGTAACTTAATACCACAACCAGGCATAGGGGAAAGATCGGAG * AS:i:144 HI:i:1 NH:i:4 XT:Z:T ZT:Z:ENST00000456328.2
SRR5024127.12905971.CTGAAGCT 163 chr1 12133 1 75M = 286 82 TGTTGTCTGCATGTAACTTAATACCACAACCAGGCATAGGGGAAAGATTGGAGGAAAGATGAGTGAGAGCATCAA * AS:i:150 HI:i:1 NH:i:4 XT:Z:T ZT:Z:ENST00000456328.2
SRR5024127.17645097.CTGAAGCT 83 chr1 14583 1 75M = 1302 -135 TACCCCAACACCAGCAATTGTGCCAAGGGCCATTAGGCTCTCAGCATGACTATTTTTAGAGACCCCGTGTCTGTC * AS:i:150 HI:i:1 NH:i:5 XT:Z:T ZT:Z:ENST00000456328.2
SRR5024127.20859782.CTGAAGCT 147 chr1 14613 1 75M = 1349 -118 CATTAGGCTCTCAGCATGACTATTTTTAGAGACCCCGTGTCTGTCACTGAAACCTTTTTTGTGGGAGACTATTCC * AS:i:150 HI:i:1 NH:i:5 XT:Z:T ZT:Z:ENST00000456328.2