Question: RNAseq consensus transcript sequence
0
gravatar for Neuls
3 months ago by
Neuls0
Neuls0 wrote:

Dear all,

I have a RNA seq pair-end data from a bacteria. I mapped the data against a well-known reference genome using bowtie2. Then, I converted and sorted the resulting .sam file to .bam using samtools and I assembled the transcripts using StringTie:

samtools view -u bowtie2.sam | samtools sort -o bowtie2.bam
stringtie bowtie2.bam -G ../GCF_000203855.3_ASM20385v3_genomic.gff -o output.gff

Finally, I would like to get the .fasta sequence for each transcript. Searching in internet I came across gffread. However, this sequence is the one contained in the reference genome and I am interested in the consensus sequence of all overlapping reads mapped in this region.

I found in internet you can a consensus sequence out of .sam files and reference genome using the following command:

samtools mpileup -uf ../GCF_000203855.3_ASM20385v3_genomic.fna BEE1ST_sorted.bam | bcftools view -cg - | vcfutils.pl vcf2fq > cns.fq

However I got this error:

[warning] samtools mpileup option `u` is functional, but deprecated. Please switch to using bcftools mpileup in future.
[mpileup] 1 samples in 1 input files
Error: Could not parse --min-ac g
Use of uninitialized value $l in numeric lt (<) at /usr/bin/vcfutils.pl line 566.
Use of uninitialized value $l in numeric lt (<) at /usr/bin/vcfutils.pl line 566

I also found that you can get a consensus sequence using IGV software. After dealing a bit with it I managed to get a consensus sequence of a region containing my gene of interest (I was having trouble with selecting regions!). I blasted the consensus against refseq_genomes and it gives me the reference genome I used. The was good, a query coverage of 100% and 99% identity.

I would like to know if there is any command line-based tool where you can extract the consensus sequences for the assembled transcripts after mapping against a reference genome. Or what is you opinion on IGV sequence consensus method.

Thank you,

rna-seq alignment assembly • 177 views
ADD COMMENTlink modified 3 months ago • written 3 months ago by Neuls0
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: 2305 users visited in the last hour