前提已經(jīng)獲得了trinity生成的.fasta文件(trinity_out_dir_all.Trinity.fasta)
在此，我們采用的是Trinity自帶的腳本進行定量計算。
一：首先生成如下的一個**samples.txt**文檔，包含了所有的樣本：
cond_J??? cond_J1??? J1_1.fq.gz??? J1_2.fq.gz
cond_J??? cond_J2??? J2_1.fq.gz??? J2_2.fq.gz
cond_J??? cond_J3??? J3_1.fq.gz??? J3_2.fq.gz
cond_SF??? cond_SF1??? SF1_1.fq.gz??? SF1_2.fq.gz
cond_SF??? cond_SF2??? SF2_1.fq.gz??? SF2_2.fq.gz
cond_SF??? cond_SF3??? SF3_1.fq.gz??? SF3_2.fq.gz
cond_SM??? cond_SM1??? SM1_1.fq.gz??? SM1_2.fq.gz
cond_SM??? cond_SM2??? SM2_1.fq.gz??? SM2_2.fq.gz
cond_SM??? cond_SM3??? SM3_1.fq.gz??? SM3_2.fq.gz
二：用trinity自帶的腳本align_and_estimate_abundance.pl進行轉(zhuǎn)錄本表達定量。
$align_and_estimate_abundance.pl? #指明調(diào)用的腳本
--transcripts? /data1/spider/ytbiosoft/data/trinity.all/trinity_out_dir_all.Trinity.fasta? #讀取之前Trinity拼接的轉(zhuǎn)錄本
--seqType fq --samples_file /data/spider_data/tamu/sample.txt? #讀取第一步創(chuàng)建的樣品信息
--est_method? RSEM? #進行表達量計算的軟件是RSEM
--aln_method bowtie2? #由于RSEM是通過比對進行的表達量計算，因此會采用的bowtie2進行比對
--trinity_mode? #這個加上會采用Trinitymode以調(diào)用前期assembly過程中的一個gene_trans_map文件
--prep_reference? #會根據(jù)拼接的fasta文件構(gòu)建index?
--output_dir /data1/spider/ytbiosoft/data/RSEMResult? #輸出文件夾，如果不指定路徑的話，由于會采用讀取樣品信息，因此會輸入到樣品信息的文件夾
--thread_count 6? #線程數(shù)
三：最后會生成6個文件夾（我有6個樣），每個文件夾內(nèi)容如下，最主要的是RSEM.genes.results與RSEM.isoforms.results。
├── bowtie2.bam #bowtie2 生成的 bam文件
├── bowtie2.bam.for_rsem.bam #用于RSEM計算的 bam文件
├── bowtie2.bam.ok
├── RSEM.genes.results #基于基因的EM Reads Count
├── RSEM.isoforms.results #基于轉(zhuǎn)錄本的 EM Reads Count
├── RSEM.isoforms.results.ok
└── RSEM.stat
? ? ? ? ├── RSEM.cnt
? ? ? ? ├── RSEM.model
? ? ? ?└── RSEM.theta
四：由于我們在每一個文件夾中的Reads count 沒有經(jīng)樣本間的均一化，因此需要做一個樣本均一化，構(gòu)建轉(zhuǎn)錄本-基因表達矩陣并得到不同樣本中的均一化表達數(shù)據(jù)TMM是后期要做的一個工作。我們采用以下的矩陣的到了三個結(jié)果
* 首先手動整理到一個文件夾中:
$find * -name '*.isoforms.results'> quant.file? #這個地方，我們采用了find命令將子文件夾中的isoforms基因表達量結(jié)果全部查找出來然后路徑放到一個文件中（quant.file）后期要使用
可以$cat quant.file文件看一下內(nèi)容
* 再調(diào)用trinity自帶腳本abundance_estimates_to_matrix.pl構(gòu)建轉(zhuǎn)錄本-基因表達矩陣
abundance_estimates_to_matrix.pl? #采用的腳本
--est_method RSEM? #由于是對RSEM的結(jié)果進行矩陣構(gòu)建，因此需要說明這個
--gene_trans_map trinity_out_dir_all.Trinity.fasta.gene_trans_map? #通過這個map構(gòu)建基因的表達量矩陣
--name_sample_by_basedir? #**這個必須要選**，不然會導(dǎo)致程序沒辦法合并之前的結(jié)果進行計算
--quant_files quant.file? #這個指明需要的上游文件的位置
* 經(jīng)過計算后得到的結(jié)果如下：
同理做*.genes.results的表達矩陣。
最后獲得*.genes.results以及*.isoforms.result的表達矩陣：
* 以上結(jié)果中分為基因的表達矩陣和轉(zhuǎn)錄本的表達矩陣。
其一： '.counts.matrix' 文件用于后期的差異表達分析
其二： '.TMM.EXPR.matrix'文件可以用于其他基因表達的分析
到此，就把相關(guān)的基因表達均一化這邊的工作做完了，應(yīng)該說經(jīng)過這個工作后所有的基因的表達量已經(jīng)完成計算，后期就是通過差異表達的軟件進行一個差異分析就可以進入可視化階段了。
歡迎聯(lián)系互相學(xué)習(xí)：909474045@qq.com
