基于Salmon的转录组批量定量流程和差异分析
循环定量多个样品的表达量
整理样本信息表，命名为sampleFile，内容如下：
Samp conditions individual
untrt_N61311 untrt N61311
untrt_N052611 untrt N052611
untrt_N080611 untrt N080611
untrt_N061011 untrt N061011
trt_N61311 trt N61311
trt_N052611 trt N052611
trt_N080611 trt N080611
trt_N061011 trt N061011采用for循环进行批量定量 (参考这个为生信学习打造的开源Bash教程真香！！，理解更多)：
for samp in `tail -n +2 sampleFile | cut -f 1`; do salmon quant --gcBias -l A -1 ${samp}_1.fq.gz -2 ${samp}_2.fq.gz -i genome/GRCh38.salmon_sa_index -o ${samp}/${samp}.salmon.count -p 4 >${samp}.salmon.log 2>&1; done &整理Salmon定量文件用于DESeq2差异基因鉴定
找到Salmon的输出文件并压缩起来，用于下载到本地进行差异分析。
# 列出salmon的输出文件
find . -name quant.sf
# 这个压缩包下载解压到本地
zip quant.sf.zip `find . -name quant.sf`./trt_N080611/trt_N080611.salmon.count/quant.sf
./trt_N061011/trt_N061011.salmon.count/quant.sf
./untrt_N61311/untrt_N61311.salmon.count/quant.sf
生成辅助文件，指出每个样品对应的自己的quant.sf文件，便于导入tximport包。
# 生成一个两列文件方便R导入
# xargs接收上一步的输出，按批次提供给下游程序作为输入
# -i: 用{}表示传递的值
cut -f 1 sampleFile | xargs -i echo -e "{}\t{}/{}.salmon.count/quant.sf" >salmon.output
head salmon.output两列文件
# Samp Samp/Samp.salmon.count/quant.sf
# untrt_N61311 untrt_N61311/untrt_N61311.salmon.count/quant.sf
# untrt_N052611 untrt_N052611/untrt_N052611.salmon.count/quant.sf获得基因和转录本的对应关系，获取基因的表达量
# 如果没有GTF文件，可以用其他文件，只需获取转录本和基因名字对应关系就可以
# 如果不知道对应关系，也可以把每个转录本当做一个基因进行分析
# Trinity拼装时会生成这个文件
# 注意修改$14, $10为对应的信息列，
# tx2gene为一个两列文件，第一列是转录本没名字，第二列是基因名字。
sed 's/"/\t/g' genome/GRCh38.gtf | awk 'BEGIN{OFS=FS="\t"}{if(FNR==1) print "TXname\tGene"; if($3=="transcript") print $14, $10}' >GRCh38.tx2gene
head GRCh38.tx2gene转录本->基因
# TXname Gene
# ENST00000608838 ENSG00000178591
# ENST00000382410 ENSG00000178591
# ENST00000382398 ENSG00000125788
# ENST00000542572 ENSG00000125788至此就完成了基于Salmon的所有样本基因和转录本的定量。然后下载sampleFile、GRCh38.tx2gene、salmon.output、quant.sf.zip文件到本地进行下游分析。
具体差异基因鉴定可参考高通量数据中批次效应的鉴定和处理 - 系列总结和更新。
往期精品(点击图片直达文字对应教程)
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
后台回复“生信宝典福利第一波”或点击阅读原文获取教程合集
（请备注姓名-学校/企业-职务等）
评论