Obtaining contig stats from Velvet+Oases

[Blaze9]

Hey all,

Here's a laydown of my information.

We sequenced a transcriptome. The company returned to us 4 different reads of that same transcriptome.

We ran FastQC and trimmed all 4 paired reads individually. We had some contamination issues, so we ran DeconSeq on all 4 reads.

We ran through the velvet+oases pipeline using Kmergenie for our k-mer count.

We concatenated the final oases outputs into one fasta file. This contained 4 repeats of our transcripts.

We ran through the Velvet+Oases + kmergenie pipeline again, this time utilizing the "merged" option.

We took the single output from Oases and ran it against CD-HITs to merge the similar sequences together.

Now, we have a final fasta file containing around 37,000 transcripts.

How can I obtain the contig stats for this final file? I've read countless papers which outline the following core-information:

Contig Number
Maximum Contig Length
Minimum Contig length
Average Contig Length
N50 Length
Number of Reads per contig

I looked at the "stats.txt" file from Oases, but nothing is given in this format.

How would I go about generating that info?

Thank you.

[Brian Bushnell]

If you download the BBTools package, it contains a program for generating assembly statistics:

stats.sh in=contigs.fa

For the number of reads per contig, though, you will need to use mapping, and the answer will be subjective, particularly for a transcriptome. That said, you can determine it like this:

bbmap.sh ref=contigs.fa in=reads.fq covstats=covstats.txt

The columns of "covstats" labeled "Plus_reads" and "Minus_reads" will tell you how many reads map to each contig. It's probably easiest if you first concatenate all of your read input files together (all the read1 files into one file and all the read2 files into another) and feed them with "in1" and "in2".

-Brian