Question: How can I create a GmapGenome that doesn't start at index 1?
0
gravatar for Kyle Johnsen
4 months ago by
Kyle Johnsen20
United States/Brigham Young University
Kyle Johnsen20 wrote:

Let me explain: I want to use a small "reference genome" of only a gene, for use in examples. But I need it to still have the coordinates that correspond to the actual genome coordinates. But when I create a GmapGenome using this truncated fasta, it begins at index 1, which prevents me from using it as wanted.

The header of the fasta file: >18 dna:chromosome chromosome:GRCz11:18:5213338:5227420:1

This is how I build the genome: refGenome <- gmapR::GmapGenome('slc24a5.fa', create=TRUE)

I want to be able to reference the genome using coordinates between 5213338 and 5227420, instead of 1 and 14083 (the length of the gene). Is this possible?

To clarify, what I need to do in the end is call variants (with VariantTools) using BAM files (aligned to genomic coordinates).

gmapr • 97 views
ADD COMMENTlink modified 4 months ago • written 4 months ago by Kyle Johnsen20
Answer: How can I create a GmapGenome that doesn't start at index 1?
1
gravatar for Michael Lawrence
4 months ago by
United States
Michael Lawrence11k wrote:

You'll have to map the resulting coordinates using e.g. mapToTranscripts(). See the GenomicAlignments package.

ADD COMMENTlink written 4 months ago by Michael Lawrence11k

Ok, in my case, I want to call variants in this region, but the BAM files are in genomic coordinates. If I understand right, mapToTranscripts() would help if I had the BAM files in transcript coordinates, called variants, and then converted to genomic coordinates. Is that right, and if so, is that the only way?

ADD REPLYlink written 4 months ago by Kyle Johnsen20

Sorry, I meant mapFromTranscripts. You could call variants in your restricted genome space and then map to the true genome using that function.

ADD REPLYlink written 4 months ago by Michael Lawrence11k

Ok. That won't exactly work for what I wanted to do, since I wanted to have the BAM files and genome all ready to plug into a pipeline of which variant calling is only a part, but thank you for the response!

ADD REPLYlink modified 4 months ago • written 4 months ago by Kyle Johnsen20
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 16.09
Traffic: 176 users visited in the last hour