Question: Is there a way to save a specific Gene in IGV of a whole exome sample?
0
gravatar for mhmtgenc85
14 months ago by
mhmtgenc8530
Turkey
mhmtgenc8530 wrote:

Dear all, Is there a way to save a specific Gene in IGV of a whole exome sample? For example I have the whole exome sequencing of a sample but I would like to view and save only one gene out of that whole BAM as only one gene specific BAM?

Thanks in advance

igv bam sequence • 397 views
ADD COMMENTlink modified 14 months ago by sacha1.8k • written 14 months ago by mhmtgenc8530
0
gravatar for cpad0112
14 months ago by
cpad011211k
India
cpad011211k wrote:

Make a bam for gene of interest (using gene coordinates), index the bam and load the indexed bam in IGV.

samtools view -h test.bam gene_coordinates
samtools view -h aln.sorted.bam chr2:20,100,000-20,200,000

Extract Alignment From Very Large Bam File: Extract Alignment From Very Large Bam File

or you can write a session file in IGV, so that IGV always zooms in to that gene of interest.

ADD COMMENTlink modified 14 months ago • written 14 months ago by cpad011211k

Thank you for the tip. I tried it and get the BAM file. But now I need to get the index file of the bam. When I use

samtools index test.bam > test.bai I get an error like this

EOF marker is absent. invalid BAM binary header (this is not a BAM file) Invalid BAM header. Fail to index..

I have tried this command and it worked.

samtools view -b test.bam > output.bam

Changed -h to -b ...

Thank you.

ADD REPLYlink modified 14 months ago • written 14 months ago by mhmtgenc8530

Are you using latest version of samtools? If so, command to index is samtools index test.bam. You should also sort the test.bam file to be safe by doing samtools sort -o test_sorted.bam test.bam.

ADD REPLYlink written 14 months ago by genomax70k

Probably bam headers are not extracted.

@OP: Try this instead:

samtools view -b aln.sorted.bam chr2:20,100,000-20,200,000 > out.bam
samtools index out.bam
ADD REPLYlink modified 14 months ago • written 14 months ago by cpad011211k

@OP: Btw, you do not have to output bai file like samtools index test.bam > test.bai. Run samtools index test.bam. That would create .bai file automatically.

ADD REPLYlink written 14 months ago by cpad011211k

if you are on ubuntu, could you please post the output from test.bam:

$ od -x test.bam | head -1

It should be 0000000 8b1f 0408 0000 0000 ff00 0006 4342 0002

ADD REPLYlink written 14 months ago by cpad011211k

Command for bam files is -b (implies -h).

ADD REPLYlink written 14 months ago by ATpoint19k

@ ATpoint. Thanks and updated the code.

ADD REPLYlink written 14 months ago by cpad011211k
0
gravatar for sacha
14 months ago by
sacha1.8k
France
sacha1.8k wrote:

Something like this ?

bedtools intersect -a exome.bam -b mygene.bed  > mygene.bam
ADD COMMENTlink written 14 months ago by sacha1.8k
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.3.0
Traffic: 1314 users visited in the last hour