????????起因是在閱讀一篇文獻(xiàn)時(shí)發(fā)現(xiàn)了這么一段話： In the initial step,we spiked 2 × 106 mRNA-GFP and 2 × 106 mRNA-RFP, transcribed in vitro, into each lysed sample for RNA isolation.

關(guān)于什么是spike-in，可以先移步從spike-in到DESeq2：文庫(kù)normalization這篇。

概括一下就是在建庫(kù)前往不同樣本的total RNA中摻入已知序列的定量外源mRNA，然后通過比較它們的reads數(shù)，能比較不同樣本間mRNA的總量，對(duì)基因表達(dá)數(shù)據(jù)進(jìn)行絕對(duì)值的矯正。

目前使用最多的spike-in是ERCC( The External RNA Control Consortium)的spike-in MIX，有mix1和mix2兩只試劑，混合有92種不同的mRNA。ERCC的試劑盒售價(jià)高達(dá)2萬(wàn)元，若非必要一般也不會(huì)去買。這篇文獻(xiàn)的作者沒有使用ERCC，所以就簡(jiǎn)單的拿GFP和RFP來(lái)當(dāng)spike in。

首先，我從EBI上下載了原文獻(xiàn)的原始數(shù)據(jù)，然后質(zhì)量過濾，得到clean fastq。由于文中未提到具體的GFP和RFP序列，我便在SnapGene里挑了常用的GFP，EGFP,mChery，DsRed1試試能不能map出來(lái)。

我用的軟件是HISAT2，首先要建立index。有兩種選擇，一是將四個(gè)熒光蛋白序列加入到小鼠基因組的fasta文件里，這樣可以一次性全讀出來(lái)；二是就用這四條序列建個(gè)index。我選擇了第二種，將四條序列加到一個(gè)fasta文件中，大致就是下圖這樣。

然后就hisat2 build ，index建好后就比對(duì)，具體指令可以看相應(yīng)的文章。

比對(duì)得到的sam文件用samtools 轉(zhuǎn)換為bam格式，然后以位置排序。我使用的計(jì)reads數(shù)的軟件是stringtie，用其他計(jì)數(shù)軟件也沒問題。

重點(diǎn)是用來(lái)計(jì)數(shù)的annotation文件該怎么處理？

GFP?havana??exon????1???717?.???+???.???gene_id?"GFP";?transcript_id?"GFP";?gene_name?"GFP";

EGFP????havana??exon????1???720?.???+???.???gene_id?"EGFP";?transcript_id?"EGFP";?gene_name?"EGFP";

mCherry?havana??exon????1???711?.???+???.???gene_id?"mCherry";?transcript_id?"mCherry";?gene_name?"mCherry";

DsRed1??havana??exon????1???681?.???+???.???gene_id?"DsRed1";?transcript_id?"DsRed1";?gene_name?"DsRed1";

創(chuàng)建一個(gè)gtf文件，將上面的信息復(fù)制進(jìn)去就好了，注意windows和linux編碼的區(qū)別。

計(jì)數(shù)結(jié)果如上，可見作者是使用EGFP和mCherry作為spike in，DsRed1和mCherry有部分序列接近，所以也顯示有少量reads。

文獻(xiàn)鏈接：BTG4 is a meiotic cell cycle–coupled maternal-zygotic-transition licensing factor in oocytes