Scripts for spike-in control statistics

Spike-in RNA-seq controls provides an important QC information. A simple QC use case is to calculate the ratio of the reads from spike-in sequences to that of the total reads in one experiment. The following code uses the idxstat function of samtools to do this.

echo "bam_file Spike-in_reads total_reads spike-in/total_reads"
for i in `ls mapping/*.bam`;
do 
  echo -n $i "";
  samtools idxstats $i|awk '{s+=$3; total+=$3 + $4} END {print s, total, s/total}';
done
Advertisements

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s