首先我们看一下传统测序(非链特异性测序)存在的问题:
https://www.researchgate.net/figure/mRNA-library-construction-workflow-for-Illumina-from-David-Corney-2013_fig1_323731847
普通的RNA-seq,一段原始的RNA序列,完成cDNA的合成后,我们无法得知哪一条链是原始的了,再加上Y型接头,进行PCR扩增后,有的原始片段可能5‘连接的是P5,有的可能连的P7,所以最后测到的序列一半和原始序列的方向相同一半相反,无法得知序列来自DNA的那一条链。
https://bmcgenomics.biomedcentral.com/articles/10.1186/s12864-015-1876-7
https://htseq.readthedocs.io/en/master/htseqcount.html
基因ICAM4在正链上,基因CTD-2369P2.8在负链上。基因ICAM4和基因CTD-2369P2.8是在正负链上相互重叠的基因。在non-stranded RNA-seq中,所有reads被同时比对到基因ICAM4和基因CTD-2369P2.8,使用htseq-count计数时重叠区域的ambiguous reads计数被排除在外,因此结果显示基因ICAM4没有表达。(所有基因、转录物和序列读数如果在“+”链中,则显示为蓝色,如果在“−”链中则显示为绿色)
non-stranded RNA-seq
htseq-count
ambiguous reads
基因ICAM4的序列在正链,当reads也被映射到正链时,htseq-count不计数,因为reads的序列不是与基因ICAM4的序列反向互补,也就是说reads来源于负链CTD-2369P2.8。
首先看一下HISAT2的命令参数--rna-strandness,其作用的用于指定strand-specific的信息。默认情况下是unstranded。
--rna-strandness
这里很容易理解,以单端测序为例,只要确保reads的方向与转录本的方向相同或相反就行。这里假设reads的方向与转录本的方向相同,HISAT2在将reads回贴到基因组时,如果reads的反向互补贴到
根据gtf文件从fastq文件中提取转录本序列
grep ENSG00000100271 ../genomic/chr22_with_ERCC92.gtf > TTLL1.gtf grep ENSG00000100304 ../genomic/chr22_with_ERCC92.gtf > TTLL12.gtf gffread -w TTLL12_transcripts.fa -g ../genomic/chr22_with_ERCC92.fa TTLL12.gtf # hisat2 --new-summary -p 30 \ # -x ../hisat2_index/chr22_with_ERCC92 \ # --dta --rna-strandness RF \ # -U out.fq \ # | samtools view -bS - | samtools sort - -o out.bam hisat2 --new-summary -p 30 \ -x ../hisat2_index/chr22_with_ERCC92 \ -U out.fq \ -S out.sam 2> out.txt samtools view -bS out.sam -o out.bam samtools sort out.bam -o out.sort.bam samtools index out.sort.bam
图片alt
https://www.biostars.org/p/70833/