Unconfigured Ad

Collapse
X
 
  • Time
  • Show
Clear All
new posts
  • jans
    Junior Member
    • Mar 2014
    • 1

    #1

    help with counting reads (only 0 counts were returned)

    Hi,
    This is my very first experience analysing RNAseq data. My goal is to do differential analysis between two strains of a bacteria. So far, i managed to align and produce SAM and BAM files. I'm having problems to annotate and count my reads. Here are the commands that I used.
    My reads are from SOLID and hence in colourspace


    Code:
    $ nohup solid2fastq.pl 291_01_01 291_01_01-bwa  #Convert .csfasta and .qual to .fastq
    
    $ nohup bwa index -c TbruceiTreu927Genomic_TriTrypDB-4.0.fasta 
    
    $ nohup bwa aln -c TbruceiTreu927Genomic_TriTrypDB-4.0.fasta 291_01_01-bwa.singleF3.fastq 291_01_01-bwa.sai 
    
    $ perl -ne 'if($_ !~ m/^\S+?\t4\t/){print $_}' 291_01_01-bwa.sam > 291_01_01-bwa.sam.filtered #Convert to SAM file
    
    $ samtools sort 291_01_01-bwa.bam 291_01_01-bwa.bam.sorted
    
    $ samtools index 291_01_01-bwa.bam.sorted.bam
    #to produce .rpkm file
    
    $ java -jar ~/bin/bam2rpkm-0.06/bam2rpkm-0.06.jar  -i 291_01_01-bwa.bam.sorted.bam -f Tbrucei427_TriTrypDB-4.0.gff > 291_01_01-bwa.RPKM2.out  # i get an error here
    $ERROR: Problem encountered whilst reading gtf file. Could not interpret line 'GeneDB|Tb427_01_v4 EuPathDB supercontig 1
    
    so i tried different method to count
    
    $ htseq-count -i ID 291_01_01-bwa.sam Tbrucei427_TriTrypDB-4.0.gff > 291_01_01-bwa.sam_htseq-count #still error 
    $Error occured when processing GFF file (line 37060 of file Tbrucei427_TriTrypDB-4.0.gff):
    need more than 1 value to unpack
    
    and different method
    
    $ python make_bed_from_fasta.py ~/Downloads/reference/TbruceiTreu927AnnotatedCDS_TriTrypDB-4.0.fasta > 927_reference.bed #this python script converts .fasta into .bed file since the .gff file cannot be processed
    $multiBamCov -q 30 -p -bams 291_01_01-bwa.bam.sorted.bam -bed 927_reference.bed > test_counts.txt
    now I only get 0 counts for all genes. Does this mean that there is something wrong with my alignment files or something wrong with the counting method . And it seems like my .gff (version 3) file was unable to be read by htseq-count and also the java script . I downloaded the gff file from GeneDB and it seems like in many tutorials .gtf files are used instead. So I'm stuck at counting the read part and I really need some help . Help please .

Latest Articles

Collapse

  • SEQadmin2
    Multiomic Liquid Biopsy: What DNA/RNA Panels Reveal That DNA Alone Can’t
    by SEQadmin2



    Tumors contain complex mixtures of cell populations that evolve in response to treatment and host-specific factors. Tissue biopsies provide only static snapshots of a tumor’s state. But what if we could follow that evolution over time using a minimally invasive technique?

    Liquid biopsies can detect tumor-derived material circulating in the bloodstream, providing a way to study treatment response, detect emerging resistance, and discover new biomarkers. Circulating
    ...
    10-14-2026, 09:55 AM
  • SEQadmin2
    New Genomics Technologies Take Aim at Long-Standing Limits
    by SEQadmin2


    Researchers using sequencing and genomics tools often have to make trade-offs. They can choose between speed or scale, short reads or long-range information, or targeted panels or a view of the whole transcriptome. New technologies that have been released this year are built to address those tough choices.

    We asked six companies the same four questions to learn about their latest products. The new technologies bring a lot to the table, including rethinking sequencing
    ...
    09-28-2026, 10:25 AM

ad_right_rmr

Collapse

News

Collapse

Topics Statistics Last Post
Started by SEQadmin2, 10-22-2026, 11:13 AM
0 responses
40 views
0 reactions
Last Post SEQadmin2  
Started by SEQadmin2, 10-22-2026, 10:48 AM
0 responses
29 views
0 reactions
Last Post SEQadmin2  
Started by SEQadmin2, 10-15-2026, 10:52 AM
0 responses
42 views
0 reactions
Last Post SEQadmin2  
Started by SEQadmin2, 10-14-2026, 05:39 AM
0 responses
45 views
0 reactions
Last Post SEQadmin2  
Working...