自学miRNA-seq分析第三讲~公共测序数据下载

Posted on by

前面已经讲到了该文章的数据已经上传到NCBI的SRA数据中心，所以直接根据索引号下载，然后用SRAtoolkit转出我们想要的fastq测序数据即可。下载的数据一般要进行质量控制，可视化展现一下质量如何，然后根据大题测序质量进行简单过滤。所以需要提前安装一些软件来完成这些任务，包括： sratoolkit /fastx_toolkit /fastqc/bowtie2/hg19/miRBase/SHRiMP

下面是我用新服务器下载安装软件的一些代码记录，因为fastx_toolkit /fastqc我已经安装过，就不列代码了，还有miRBase的下载，我在前面第二讲里面提到过，传送门：自学miRNA-seq分析第二讲~学习资料的搜集

## pre-step: download sratoolkit /fastx_toolkit_0.0.13/fastqc/bowtie2/hg19/miRBase/SHRiMP
 ## 我这里特意挑选的二进制版本程序下载的，这样直接解压就可以用，但是需要挑选适合自己的操作系统的程序。
cd ~/biosoft
mkdir sratoolkit &&  cd sratoolkit
##
##  Length: 63453761 (61M) [application/x-gzip]
##  Saving to: "sratoolkit.2.6.3-centos_linux64.tar.gz"
tar zxvf sratoolkit.2.6.3-centos_linux64.tar.gz
cd ~/biosoft
mkdir bowtie &&  cd bowtie
#Length: 27073243 (26M) [application/octet-stream]
#Saving to: "download"
 mv download  bowtie2-2.2.9-linux-x86_64.zip
 unzip bowtie2-2.2.9-linux-x86_64.zip
mkdir SHRiMP &&  cd SHRiMP
tar zxvf SHRiMP_2_2_3.lx26.x86_64.tar.gz 
cd SHRiMP_2_2_3
export SHRIMP_FOLDER=$PWD  ## 这个软件使用的时候比较奇葩，需要设置到环境变量，不能简单的调用全路径
SHRiMP这个软件比较小众，我也是第一次听说过，本来我计划是能用bowtie搞定，就不麻烦了，但是第一次比对出了一个bug，就是下载的miRNA序列里面的U没有转换成T，所以导致比对率非常之低，所以我不得不根据文章里面记录的软件SHRiMP 来做比对，最后发现比对率完全没有改善，搞得我都在怀疑是不是作者乱来了。
下面是下载数据，质量控制的代码，希望大家可以照着运行一下：
## step1 : download raw data
mkdir miRNA_test && cd miRNA_test
echo {14..19} |sed 's/ /\n/g' |while read id; \
done
## step2 :  change sra data to fastq files.
## 主要是用shell脚本来批量下载
ls *sra |while read id; do ~/biosoft/sratoolkit/sratoolkit.2.6.3-centos_linux64/bin/fastq-dump $id;done
rm *sra
##  33M --> 247M
#Read 1866654 spots for SRR1542714.sra
#Written 1866654 spots for SRR1542714.sra
## step3 : download the results from paper
mkdir paper_results && cd paper_results
## tar xvf GSE60292_RAW.tar
ls *gz |while read id ; do (echo $id;zcat $id | cut -f 2 |perl -alne '{$t+=$_;}END{print $t}');done
ls *gz |xargs gunzip
## step4 : quality assessment
ls *fastq | while read id ; do ~/biosoft/fastqc/FastQC/fastqc $id;done
## Sequence length 8-109
## %GC 52
## Adapter Content passed
## write a script : :: cat >filter.sh
ls *fastq |while read id
do
echo $id
~/biosoft/fastx_toolkit_0.0.13/bin/fastq_quality_filter -v -q 20 -p 80 -Q33  -i $id -o tmp ;
~/biosoft/fastx_toolkit_0.0.13/bin/fastx_trimmer -v -f 1 -l 27 -i tmp  -Q33 -z -o ${id%%.*}_clean.fq.gz ;
done
rm tmp
## discarded 12%~~49%%
ls *_clean.fq.gz | while read id ; do ~/biosoft/fastqc/FastQC/fastqc $id;done
mkdir QC_results
mv *zip *html QC_results
这个代码是我自己根据文章的理解写出的，因为我本身不擅长miRNA数据分析，所以在进行QC的时候参数选择可能并不是那么友好，如果有高手能指正就最好了，可以直接打我电话告诉我，或者发邮箱给我，邮箱用户名是jmzeng1314，是163邮箱。
~/biosoft/fastx_toolkit_0.0.13/bin/fastq_quality_filter -v -q 20 -p 80 -Q33  -i $id -o tmp ;
~/biosoft/fastx_toolkit_0.0.13/bin/fastx_trimmer -v -f 1 -l 27 -i tmp  -Q33 -z -o ${id%%.*}_clean.fq.gz ;
最后得到的clean.fq.gz系列文件，就是我需要进行比对的序列啦。

Comments are closed.