1. fastqc 检查下数据质量

结果显示已去除接头

2. flash 双端测序reads 进行拼接

flash -t 6 MTO1_IP.1.fq MTO1_IP.2.fq -p 33 -r 150 -s 100 -o MTO1_IP_merge

生成了6个文件

image.png

FLASH拼接默认输出6个结果文件：

extendeFrags.fastq 为拼接后的扩增片段序列文件；

output.flash.log 为日志文件，详细记录了拼接过程中的参数和拼接统计的数据；

output.hist 为拼接后的reads长度的统计信息文件；

image.png

output.histogram 为拼接后的reads长度直方图文件；

image.png

output.notCombined_1.fastq 为拼接不上的reads1序列文件；

output.notCombined_2.fastq 为拼接不上的reads2序列文件；
没匹配上的文件根据需要再进行下步分析

3 fastq_quality_filter 去除低质量的reads

fastq_quality_filter -q 20 -p 80 -i MTO1_IP_merge.extendedFrags.fastq -o MTO1_IP_mer_clean.fq -Q 33

image.png

4 bowtie2 对比

bowtie2-build 建立索引, ensemble ftp下载人线粒体序列fasta文件和gff3注释文件

bowtie2-build ./index/Homo_sapiens.GRCh38.dna.chromosome.MT.fa ./index/human_mt_index

开始比对

bowtie2 -p 10 -x ./index/human_mt_index  MTO1_IP_mer_clean.fq \
| samtools sort -O bam -@ 10 -o - > MTO1_IP_mer_clean.bam

得到bam文件

image.png

5 计算counts

samtools idxstats适合标注的序列进行的counts计数，如小RNA或者tRNA的序列进行索引和匹对，该结果得出来得结果如下

samtools idxstats MTO1_IP_mer_clean.bam | cut -f 1,3 > MTO1_IP_mer_clean.bam.counts.txt

image.png

对线粒体16596bp 序列没有进行注释。

6 stringtie 计算count.txt (包括FPKM和RPM)

sudo apt install stringtie

stringtie MTO1_IP_mer_clean.bam -p 16 -G ./gtf/Danio_rerio_Ensemble_97.gtf -B -o ./A/TRANHOM.gtf

这里生成了结果gtf文件和ballgown需要的.ctab文件，还有基因的表达量文件gene_abund.tab，该文件包括基因的表达量FPKM以及TPM等。当然如果你想要转录本的表达量，直接打开t_data.ctab这个文件，这里面有转录本的FPKM值。

image.png

image.png

#################################################
多个转录本 stringtie --merge 合并多个gtf文件

image.png

$ stringtie --merge -p 8 -G Danio_rerio_Ensemble_97.gtf -o \
./A/stringtie_merged.gtf ./A/mergelist.txt.txt
或
stringtie --merge -p 8 -G ./GTF/Danio_rerio_Ensemble_97.gtf -o \
./A/stringtie_merged_MCK_HOM.gtf ./a/mergelist.txt

###########################################

7. featureCounts也可以对bam文件计算counts数

featureCounts -T 6 -t exon -g exon_id -a ./gtf/Homo_sapiens.GRCh38.107.chromosome.MT.gff3/Homo_sapiens.GRCh38.107.chromosome.M
T.gff3 -o MTO1_IP.ele.txt MTO1_IP_mer_clean.bam

-g exon_id 需要根据实际的gff文件进行选择