Interpret IGV output inversion
0
0
Entering edit mode
15 months ago
pablo ▴ 310

Hi,

I have these IGV outputs (I do not filter supplementary and secondary alignments) :igv

I suspect a "duplicated inversion" in an other region on the same chromosome :

igv

All the soft clipped reads in the first region are found in the other one. Also, when a soft-clipped read is aligned on the (+) strand in the first region, it is aligned on the reverse strand in the other region. And vice versa.

Can we consider this variation as a "duplicated inversion" ?

Best

variant-calling IGV alignment • 1.9k views
ADD COMMENT
0
Entering edit mode

if these are paired end, make sure to look at the orientation https://software.broadinstitute.org/software/igv/interpreting_pair_orientations

ADD REPLY
0
Entering edit mode

There are Pacbio reads (~8kb mean), these are not paired end then.

ADD REPLY
0
Entering edit mode

gotcha that is good to know. I have actually been working on some ways in jbrowse 2 (disclaimer: i'm a jb2 dev) to try to help with this. for example here is a screenshot showing an inversion with some long reads on hg002. in this case the nanopore reads span the whole inversion but these visualizations may help in your case too

enter image description here

i can try to make a tutorial about how to set it up since i know it's a bit tricky to use jbrowse 2 but you can get jbrowse desktop from jbrowse.org or setup the web version. then the basic idea is after you have added your bam or cram track, you can navigate to location of interest, then perform "replace lower panel with..."->"arc display" for example, and you can make multiple copies of your track to show different forms of the track at once

enter image description here

in igv if you are able to load the whole region in view, you can try link supplementary alignments perhaps also

ADD REPLY
0
Entering edit mode

Thanks for you reply. I am able to import and visualize my alignment in jb2, which looks like incredible btw! I'm trying to get used to it.

First, I am not able to load the whole region (neither with IGV nor with JB2) : chr9:104,601,207 - chr9:23,773,494 = 80kb.

I can spot supplementary alignments : the corresponding primary alignment is located to the other region, on the other strand. I mean, if I spot a sec. alignment at 23,773,494 pb on the (+) strand , I'll find "its" primary alignment at 104,601,207 pb on the (-) strand. And vice versa.

Nothing to report when I "link supplementary alignments" with IGV, the same when I perfom the "arc display" mode on JB2.

I did a variant calling run with pbsv. It found a translocation (called "BND") : pbsv.BND.CM000671.2:23773494-CM000671.2:104601207 ; pbsv.BND.CM000671.2:104601207-CM000671.2:23773494 which corresponds to my regions.

ADD REPLY
1
Entering edit mode

oh ya that is 80 megabase...that is tricky. I would maybe try creating a BigWig file showing coverage to see any large scale patterns in the data.

you can also try loading the vcf from pbsv into the "sv inspector" in jbrowse 2, this is an entirely different workflow in jbrowse 2 that can show read evidence for a breakpoint, see https://jbrowse.org/jb2/docs/user_guides/sv_inspector_view/

another thing you can do in jbrowse 2 is put "discontinuous regions" in the search bar, so you can put "just the breakpoints" in a view, so you can enter a space separated list of locstrings like "chr9:104,600,000-104,602,000 chr9:23,772,000-23,774,000"

http://splitthreader.com/ and https://genomeribbon.com are also pretty interesting, might be able to help

ADD REPLY
1
Entering edit mode

Thanks for the tips! I created a BigWig file but I do not notice any drop coverage along these regions.

That's why I get when I run the SV inspector tool (amazing tool)

jb2

I have a look at Ribbon which also seems a great tool.

ADD REPLY
1
Entering edit mode

nice...that does sort of put together a picture like this which is something like a large inversion. could plot mappability bigwig and repeat track just to see that it looks sort of confident but it does look pretty confident to me!

example diagram if the lines were "straightened out" into linear segments in powerpoint

enter image description here

ADD REPLY
1
Entering edit mode

Is there a way to "straighten it out" in JB2 ?

Btw, I had a look at the bw coverage along the chromosome. There is a 0X coverage between 40mb-60mb (it is the centromeric region, but I could expect a bit of coverage) :

coverage

Indeed, I think I can consider your opinion as the fact it is a large inversion. But I cannot explain yet why there is a coverage between 23-40mb and after 60mb-104mb. I mean, if it a inversion, should the reads map to the reference in these intervals?

ADD REPLY
0
Entering edit mode

I mean, if it a inversion, should the reads map to the reference in these intervals?

yes i think this is expected, reads will map inside the inversion ("in the opposite direction" perhaps, but the aligner doesnt care about this, and it's hard to tell bioinformatically what the "right" direction is unless you have stranded sequencing which is somewhat uncommon)

Is there a way to "straighten it out" in JB2 ?

there could be a couple ways i could try to imagine but it is worth brainstorming more... there is an "arc" track type, but it only works on e.g. plain BED data with the start and end of features connected.

you can also enter locstrings like

"chr9:1-23,000,000 chr9:23,000,001-104,000,000[rev] chr9:104,000,001-200,000,000"

and it will then create these three "displayed regions" side by side, with that middle segment reversed, and any tracks that you open will obey that. the locstrings can be discontiguous also if you wanted to zoom in just on the breakpoints (as mentioned above it could be like "chr9:104,600,000-104,602,000 chr9:23,772,000-23,774,000" to just be the breakpoints)

ADD REPLY
1
Entering edit mode
  • Well! You agree, I've already seen in some papers reads mapped on a inverted region (large or not). Also, any idea why there's no coverage between 43959908 - 60518558 positions (I checked with bedtools). It still remains weird for me.
  • I'm gonna have a look at your suggestion.

Thanks for your time! I appreciate.

ADD REPLY
0
Entering edit mode

i'm not very knowledgeable about centromere stuff but I would imagine that it's sort of expected to go to zero coverage in the centromere areas. sometimes these are fully repeatmasked by aligners, a long string of N's for example in the fasta file. the T2T assemblies might be better and not have N gaps but i dont have experience with aligning much and none with aligning to T2T

here are some random bigwig files coverage on my jbrowse with hg19, the coordinates are different on hg38 though

enter image description here

ADD REPLY

Login before adding your answer.

Traffic: 1800 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6