5软件安装在上一篇资料搜集的链接中大神们已经详细介绍
以下是我参照链接自己跑的流程命令
比对 hisat2
$ hisat2 -p 16 -x /path/index -1 R1.fq.gz -2 R2.fq.gz -S R1-2.sam
hisat2比对时报错：
File "/home/panda/miniconda3/envs/py3/bin/hisat2_read_statistics.py", line 210, in <module>
reads_stat(args.read_file, args.read_count)
File "/home/panda/miniconda3/envs/py3/bin/hisat2_read_statistics.py", line 168, in reads_stat
for id, seq in fstream:
File "/home/panda/miniconda3/envs/py3/bin/hisat2_read_statistics.py", line 44, in parser_FQ
if line[0] == '@':
IndexError: index out of range
可参考：Question: hisat2 error index out of range ； Question: HISAT2: an "IndexError: index out of range"
samtools对SAM文件处理
sam转bam
$ samtools view -bS R1-2.SAM > R1-2.BAM
bam文件排序（-n是按照read名称，默认是-p染色体位置。这里就按照默认的染色体位置进行排序）
$ samtools sort R1-2.bam -o R1-2_sorted.bam
排序后的bam文件建立索引（按照read名称排序建立索引时报错，按照默认染色体位置时不报错，原因没有探究）
$ samtools index R1-2_sorted.bam
产生一个R1-2_sorted.bam.bai索引文件，在将R1-2_sorted.bam载入IGV可视化时可以用到！
RSeQC质控（使用方法官方网站这里，还可以参考这里和这里）
统计分析，对bam文件进行统计分析
$ bam_stat.py -i R1-2_sorted.bam
查看基因覆盖率，需要用到RefSeq.bed文件，前往RSeQC的网站这里下载，或者可以用gtf转换（还不晓得怎样转换）
$ read_distrbution.py -i R1-2_sorted.bam -r mm10_RefSeq.bed
reads计数htseq-count
这里需要下载小鼠注释文件，前往这里下载（此处GTF文件不用排序，若是传入到IGV则需要排序，可以使用IGV的tools中run igvtools的sort功能排序）
$ htseq-count -s no -r pos -f bam R1.sorted.bam gencode.vM23.chr_patch_hapl_scaff.annotation.gtf > R1.sorted_matrix.count 2> R1.sorted_matrix.log
网友评论