加帽端测序(CAGE-seq)

一、CAGE-seq数据分析方案流程

二、CAGE-seq测序数据质量分析

Clean reads数据质量分析

方法描述：碱基质量值是衡量测序质量的重要指标，碱基质量（Q）与测序错误率（P）密切相关，受测序仪状态，测序试剂质量，样本特性等的影响。质量值计算公式如下：

结果展示

GC含量分析

方法描述：对测序reads中四种碱基的分布比例进行评估，检查是否存在AT、CG分离现象，理论上A与T、C与G的含量在整个测序反应中分别相同，且维持 在稳定水平。

有效长度统计

方法描述：去掉index序列、建库平衡用随机碱基及截取掉后面低质量的碱基后，我们用获得的clean reads进行有效长度分析。

reads冗余度统计

在cDNA文库构建的过程中对捕获的mRNA/ncRNA 进行随机片段化，随后加接头并进行RT-PCR。一个多样性的文库中大多数序列应该只出现一次，低水平的序列冗余度往往表明高水平的靶标序列覆盖度，而高水平的序列冗余度则意味着一定程度上的偏好富集性，如文库构建过程中PCR过度扩增。通常测序深度越高，越容易产生一定程度的重复reads，属于正常的现象。实际操作中，由于数据量较大，为了降低计算中对内存的要求，仅选取了每个文件的前200,000条reads进行分析，认为其可以代表全部序列的冗余度。

PCR duplication level计算方法为：从测序数据中随机挑选20万reads作为Total Reads，按照如下公式进行计算：PCR duplication level=Duplication Reads/Total Reads

三、RNA-seq结果展示

全基因组定位分析

Reads比对到参考基因组结果

方法描述：分析中对真核生物的转录组测序数据采用tophat2，原核生物的转录组测序数据采用bowtie2将reads比对到参考基因组上，比对中根据需要会设定一定的容错率。对于转录组数据，在参考基因组选择合适并且实验不存在污染的情况下，测序序列的定位到参考基因组的百分比通常会高于70%(Total Mapped Reads)，uniquely mapped reads（单位置比对）绝大多数来自成熟的mRNA，而多位置比对（Multiple mapped）的reads以rRNA和tRNA为主，其在clean reads中所占比例通常不到10%。下游所有分析仅使用了uniquely mapped reads，以保证分析结果的可信度。

Reads在基因组不同区域的分布情况

方法描述：将在参考基因组上有唯一定位的reads（uniquely mapped Reads)在基因组上各区域的分布情况进行统计，转录组测序中reads通常会在CDS区域富集。定位到Intron (内含子) 区域的测序序列可能是由于非成熟的mRNA的污染或者基因组注释不完全导致的，而定位到Intergenic（基因间隔区域）的测序序列可能为基因组注释不完全以及背景噪音。

Reads在染色体上的分布情况

方法描述：把基因组平均分成100000个bin，根据比在基因组上有唯一定位的（uniquely mapped）reads数，统计落在每个bin中的reads的平均depth，然后取log2，使用circos作图。

基因组的覆盖度和特征分析

方法描述：reads随转录单元长度的覆盖强度分析，以距转录起始位点和转录终止位点为标准，把cDNA平均分成100份，每一份称为一个bin，求落在每个bin中的reads平均数之和，从而得到每个bin上整体的reads覆盖度。

方法描述：把基因平均分成100份，每一份称为一个bin，求落在每个bin中的reads平均数之和，从而得到每个bin上整体的reads覆盖度。

饱和曲线检查

分别对10%，20%，30% … 90%的测序量各自进行基因定量分析，并将以完整测序量分析得到的基因表达水平作为最终表达水平。用各个百分比的数据量得到的基因表达水平和和最终表达水平进行比较，如果差异小于15%，则认为该基因在该数据量条件下被准确定量。

定量饱和曲线检查反映了基因表达水平定量对数据量的要求。表达量高的基因，容易被准确定量，而表达量低的基因，则需要较大的测序量才能被准确定量。

reads在转录起始位点，转录终止位点，起始密码子和终止密码子附近的分布

方法描述：分别以转录起始位点（TSS）和转录终止位点（TTS）为原点，统计其上下游1kb范围内reads的分布情况，结果如下：

方法描述：分别以起始密码子（start codon）和终止密码子（stop codon）为原点，统计其上下游1kb范围内reads的分布情况，结果如下：

基因表达情况

方法描述：一个基因表达水平的直接体现是其转录本的丰度，转录本丰度程度越高，则基因表达水平越高。在RNA-seq分析中，我们通过定位到基因组区域或基因外显子区的测序序列（reads）的计数来估计基因的表达水平。Reads计数除了与基因的真实表达水平成正比外，还与基因的长度、测序深度成正相关。为了使不同基因、不同实验间估计的基因表达水平具有可比性，引入了RPKM的概念，RPKM（reads per kilobase per millionmapped reads）是每百万reads中来自某一基因每千碱基长度的reads数目。RPKM同时考虑了测序深度和基因长度对reads计数的影响，是目前最为常用的基因表达水平估算方法。通过对mRNA长度和测序深度进行均一化(RPKM)，使不同测序样本之间的表达丰度具有可比性，消除了因mRNA长度和不同样本之间测序丰度差异带来的偏差。计算公式：

RPKM=exon reads/(unique mapped reads*exon length)，在此公式中，unique mapped reads为每个样本中在参考基因组上有唯一定位的reads数目（以百万位单位），exon length为每个基因的外显子长度之和（以kb为单位），而分子为每个基因上的reads数目。此处图为基因的RPKM值累计分布图，横坐标表示RPKM值，纵坐标表示基因百分比，蓝色线上的每个点的横纵坐标分别表示RPKM值小于等于该值的基因所占的百分比。

四、CAGE位点分析

CAGE位点检出

方法描述：TC在基因组上的分布情况统计如下：

CAGE位点在基因组上的分布

方法描述：TC在基因组上的分布情况统计如下：

CAGE位点距离统计

方法描述：每个样本中TC之间的距离统计结果如下：

样本相关性分析

方法描述：我们通过比较同一CAGE位点在不同样本间的TPM，进行任意两个样本之间的相关性分析，检查不同样品之间相关性。

TC差异分析

方法描述：使用CAGE进行样本间比较，根据位置差异计算出CAGE位点的shitf score，并用柯尔莫诺夫-斯米尔诺夫检验(Kolmogorov-Smirnov test)进行显著性分析，获得不同样本间位置和丰度发生显著变化的CAGE位点。

TC发生差异的基因的Gene Ontology富集分析

方法描述：

• 利用blast将参考基因组的基因序列比对到Gene Ontology数据库，进行GO注释；
• 提取比对结果，作为背景即background；
• 根据shifting promoter分析结果，统计每个promoter所属基因所在的GO Term，根据每个Term的基因数目，以及背景中此Term的基因数目，用Fisher Exact Test分析每个Term的显著性；
• 选取排名前10的GO Term及其校正p-value和百分比作图进行展示。

TC发生差异的基因的KEGG Pathway富集分析

方法描述：

1. 利用blast将参考基因组的基因序列比对到KEGG数据库，进行KEGG注释；
2. 提取比对结果，作为背景即background；
3. 根据shifting promoter分析结果，统计每个promoter所属基因所在的KEGG pathway，根据每个pathway的基因数目，以及背景中此pathway的基因数目，用Fisher Exact Test分析每个pathway的显著性；
4. 选取排名前10的pathway及其校正p-value和百分比作图进行展示。