现在纯生信SCI文章真是太难发了，需结合基础实验。但如果生信文章有思想的创新，发文章还是有办法的。今天小伍和大家来聊下一个新的通过TCGA数据挖掘双向基因对的表征数据，给大家快速毕业及职称评定提供一条新的通路。
“双向基因对”概念：表示特定的头对头基因，其中两个基因的转录起始位点位于一个kb区域内的基因组DNA的相反链上。尽管双向基因对已被充分表征，但对其在肿瘤发生中的表达谱和调节特征知之甚少。
我们使用来自癌症基因组图谱（TCGA）数据集的RNA-seq数据来系统分析13个癌症数据集中双向基因对的表达谱。在相对链上的基因对与转录末端位点距离在一个kb内或在相同链上，两个基因在1-10kb之间的距离和包含两个随机选择的基因的基因对被用作对照基因对（CG1，CG2和随机）。我们通过比较13个TCGA数据集中肿瘤与邻近正常组织之间的表达水平来鉴定和表征上调/下调基因。在大多数TCGA数据集中，双向和对照/随机基因之间上调/下调基因的百分比没有一致的显着差异。
然而，TCGA数据集和所有13个TCGA数据集中的随机基因对中，包含两个上调或下调基因的双向基因对的百分比显着高于来自CG1 / 2的基因对。然后我们确定了甲基化相关的双向基因，以探索双向基因的调控机制。像差异表达的基因对一样，与CG2 /随机基因对比较，一对中的双向基因在12/13 TCGA数据集中显着倾向于是低甲基化或超甲基化相关基因，尽管低/高甲基化百分比之间没有显着差异大多数TCGA数据集中双向和CG2 /随机基因的基因。
最后，我们探讨了双向基因与患者生存之间的相关性，确定了每个TCGA数据集中的预后双向基因和预后双向基因对。值得注意的是，我们发现了一组预后双向基因对，其中具有不同表达水平的两种蛋白质编码基因的组合与OS的存活分析中的不同存活预后相关。
这就是一篇生信文章的思路，相信花时间实现这个过程不会是很难。
步骤：
从GENCODE（https://gdc.cancer.gov/about-data/data-harmonization-and-generation/gdc-reference-files）检索基因注释GDC.h38 v.22 。每个基因的基因类型的定义包含在注释文件中。
为了鉴定双向基因对，我们使用了先前工作中提出的方法（Trinklein等，2004）。简而言之，对于每个基因，我们将TSS和转录终点（TES）的位置分别定义为分析基因的5'最边界和3'最多边界转录物（头对头结构）。我们将双向基因定义为相反链上的两个基因，TSS距离在1kb内。我们排除了基因完全位于另一个基因中的基因对，这种情况通常被定义为嵌套基因。
为了鉴定两个对照基因对（CG1和CG2），我们采用相同的方式，除了在相反的链上搜索TES距离在1kb内作为对照基因对1（CG1）和在同一条链上的距离为作为对照基因对2（CG2）的1至10kb之间的两个基因（Trinklein等，2004）。最后，我们排除了含有双向基因的对照基因对。
由于上述双向基因对和对照基因对的基因类型的数量和组成是不同的，我们还收集了一组基因对，其包含在基因组中随机选择的两个基因作为对照。每组随机基因对的基因类型的数量和组成与双向基因对相同。我们重复了100次选择并获得了100组随机基因对。
来自RNA-seq的所有TCGA表达数据集均从TCGA网站（https://portal.gdc.cancer.gov/）下载。原始读数计数是从具有后缀“htseq.counts”的文件中提取的。我们仅将≥200个癌症样本和≥15个正常样本的癌症数据集合并，以准确鉴定差异表达的基因。将正常样品定义为标记为“固体组织正常”的样品，并且没有血液样品作为正常样品包括在分析中。共有13个癌症数据集符合这些标准。除甲基化数据分析外，对这13个癌症数据集进行所有后续分析。我们在甲基化分析中排除了STAD数据集中的甲基化数据，该数据仅适用于该数据集中的两个正常样本。
我们在“RUVSeq”包（v3.8）中使用“RUVg”功能来校正先前研究后的RNA-seq数据集中的批次效应（Aran等，2017 ; Risso等，2014）。差异表达基因的鉴定通过“edgeR”包（v3.24.0）中的“glmTreat”功能进行（Robinson，McCarthy＆Smyth，2010））。我们在每个TCGA数据集中保留基因，当它们的百万分之数在至少是正常样本的大小时≥0.1时，并将这些基因定义为表达基因。然后我们通过比较癌症和正常样本组之间特定基因的表达谱来鉴定差异表达的基因。将FDR≤0.05且表达倍数变化≥1.5（上调基因）或≤0.67（下调基因）的基因定义为差异表达基因（Xiao等，2018）。修正的M值（TMM）归一化表达值的平均值也由“edgeR”包生成。
基因本体论（GO）分析在DAVID GO生物信息学平台（V6.8）（https://david.ncifcrf.gov/）上进行（Dennis等，2003）。每个基因的ensembl基因ID用作输入。未经DAVID GO网站识别的基因被排除在GO分析之外。
来自Human Methylation 450 BeadChip平台的DNA甲基化数据从TCGA网站（https://portal.gdc.cancer.gov/）下载。我们提取β值以评估每个探针的DNA甲基化水平。探针对特定基因的注释被定义为这些探针位于基因的启动子区域。启动子区域定义为基因TSS下游2kb上游至500bp下游的区域（Tang＆Epstein，2007）。我们排除了STAD数据集中仅有两个正常样本可用的甲基化数据。
我们在R中的“ChAMP”包中使用“champ.DMP”功能来鉴定差异甲基化探针（Morris等，2014）。我们将调整后的p值≤0.05的探针定义为差异甲基化探针。
然后我们将该基因定义为与甲基化相关的基因，当该基因在癌症中上调时，至少有一个低甲基化探针被注释到该基因，或者当该基因在癌症中被至少一个超甲基化探针下调时注释了这个基因。
生存分析基本上如先前报道的那样进行（Anaya等，2015）。简而言之，对于每个癌症数据集，从“TCGA-CDR”数据集中检索存活信息，其为每个患者保留最新的随访信息（Liu等，2018）。由于疾病特异性存活的信息在大多数TCGA数据集中是近似的，我们将其排除在我们的分析之外并保留其他三个终点（总生存期（OS），无病间期（DFI）和无进展期间隔（PFI）） ）。
Cox模型使用R中“生存”包的“coxph”函数运行，双向基因的方程式为“coxph（Surv（时间，检验员）~exrsrs”），其中时间是OS或疾病的生存时间 - DFI / PFI的免费/无进展时间，即使对于每个患者的DFI / PFI的OS或间隔事件，检查员也存活，并且exprs是每个基因的TMM标准化表达值。
对于双向基因对，该等式是“coxph（Surv（时间，检验员）~exprs1 + exprs2），”其中exprs1和exprs2是一对中gene1和gene2的TMM标准化表达值。我们用Logrank p值≤0.05 标记基因或基因对作为预后基因或基因对。
对于每个预后双向基因对，我们基于一对中两个基因的表达水平的组合将样品分成四组（high_vs_high，low_vs_low，high_vs_low和low_vs_high）。对组进行分类的截止值是一对中基因1和基因2的TMM标准化表达值的中值。然后使用“coxph”函数运行Cox模型，双向基因对的方程是“coxph（Surv（时间，检验员）〜组）。”我们选择Logrank p值≤0.05 的双向基因对并考虑这些双向基因配对作为对，其中具有不同表达水平的两个基因的组合与不同的存活结果相关。
我们没有对生存分析进行多项测试调整，因为之前的一项研究表明，16个TCGA数据集中有9个预测基因的数量不超过50个，调整后16个TCGA数据集中有5个中有0个（Anaya等。，2015）。
结果：
（A）双向基因对（BG），对照基因对（CG1和CG2）和随机基因对组织的示意图。TSS / TES，转录开始/结束站点。基因上的箭头表示转录方向。（B）双向和对照基因（CG1）中基因类型的分布。（C）双向和对照基因对（CG1）中基因类型组合的分布。（D）基因本体论（GO）显示双向基因的功能富集。红色虚线表示显着截止。
（A）TCGA数据集中表达和未表达的双向基因的百分比。X轴，TCGA数据集。（B）TCGA数据集中表达和未表达的对照基因（CG1）的百分比。X轴，TCGA数据集。（C）TCGA数据集中上调的双向和对照基因（CG1）的百分比。X轴，TCGA数据集。的p -值通过在R.“NS”，“chisq.test”功能指示不显著计算。“**”表示p值<0.01。（D）TCGA数据集中下调的双向和对照基因（CG1）的百分比。X轴，TCGA数据集。该p-values由R.“ns”中的“chisq.test”函数计算，表示不重要。“**”表示p值<0.01。“***”表示p值<0.001。
（A）不同TCGA数据集中表达和未表达的双向基因对的百分比。X轴，TCGA数据集。（B）不同TCGA数据集中表达和未表达的对照基因对（CG1）的百分比。X轴，TCGA数据集。（C）成对中差分表达的双向组合的模式的百分比。X轴，TCGA数据集。（D）成对（CG1）中差异表达的对照组合的模式的百分比。X轴，TCGA数据集。
（A）在特定数量的TCGA数据集中上调的双向，对照和随机基因的百分比。X轴上的每个数字“n” 表示在13个TCGA数据集的“n”中Y轴上的上调基因的相应百分比。随机基因的百分比代表100个随机组的平均值，误差条代表标准偏差。的p -值通过在R.“chisq.test”函数来计算在“chisq.test”使用TCGA数据集的特定数目的上调或下调的基因的随机数是100个随机组的平均值。（B）在特定数量的TCGA数据集中下调的双向，对照和随机基因的百分比。X上的每个数字“n”-axis表示在13个TCGA数据集的“n”中Y轴上下调基因的相应百分比。随机基因的百分比代表100个随机组的平均值，误差条代表标准偏差。的p -值通过在R.“chisq.test”函数来计算在“chisq.test”使用TCGA数据集的特定数目的上调或下调的基因的随机数是100个随机组的平均值。（C）在特定数量的TCGA数据集中具有上调或下调基因的特定组合的双向，对照和随机基因对的百分比。X轴上的每个数字“n” 表示Y上的相应百分比 - 13个TCGA数据集的“n”中上调或下调基因对的特定组合的轴。随机基因的百分比代表100个随机组的平均值，误差条代表标准偏差。的p -值通过在R.“chisq.test”函数来计算在“chisq.test”使用TCGA数据集的特定数量的特定组合的随机基因对的数量是100个随机组的平均值。（D）在特定数量的TCGA数据集中具有上调或下调基因的特定组合的双向，对照和随机基因对的百分比。X轴上的每个数字“n” 表示Y上的相应百分比 - 13个TCGA数据集的“n”中上调或下调基因对的特定组合的轴。随机基因的百分比代表100个随机组的平均值，误差条代表标准偏差。的p -值通过在R.“chisq.test”函数来计算在“chisq.test”使用TCGA数据集的特定数量的特定组合的随机基因对的数量是100个随机组的平均值。
（A）TCGA数据集中低甲基化相关双向和对照基因（CG1）的百分比。X轴，TCGA数据集。的p -值通过在R.“NS”，“chisq.test”功能指示不显著计算。“**”表示p值<0.01。“***”表示p值<0.001。（B）TCGA数据集中超甲基化相关双向和对照基因（CG1）的百分比。X轴，TCGA数据集。的p -值通过在R.“NS”，“chisq.test”功能指示不显著计算。“*”表示p值<0.05。“**”表示p值<0.01。“***”表示p - 值<0.001。（C）TCGA数据集中低或高甲基化相关双向基因组合的百分比。X轴，TCGA数据集。（D）TCGA数据集中低或超甲基化相关对照基因（CG2）组合的百分比。X轴，TCGA数据集。
（A）在OS的存活分析中包含两种蛋白质编码基因的所有预后双向和对照基因对（CG1）中包含两种蛋白质编码基因的候选预后基因对的百分比。候选预后基因对被鉴定为该对中的两个基因，其不同的表达水平与不同的存活结果相关。的p -值通过在R.“NS”，“chisq.test”功能指示不显著计算。“*”表示p - 值<0.05。（B）在OS的存活分析中包含两种蛋白质编码基因的所有预后双向和对照基因对（CG2）中包含两种蛋白质编码基因的候选预后基因对的百分比。候选预后基因对被鉴定为该对中的两个基因，其不同的表达水平与不同的存活结果相关。的p -值通过在R.“NS”，“chisq.test”功能指示不显著计算。“*”表示p值<0.05。“**”表示p - 值<0.01。（C）双向基因的存活图，ENSG00000040531.13。使用TMM标准化表达值的中值作为截止值将样品分成四组。（D）双向基因的存活图，ENSG00000262304.1。使用TMM标准化表达值的中值作为截止值将样品分成四组。（E）双向基因对的存活图。使用TMM标准化表达值的中值作为截止值将样品分成四组。
如果您觉得有价值，请把此文放到您朋友圈，大家都会感谢你
联系客服
