对于一些前期没有研究基础的研究团队来说,通过TCGA数据库挖掘候选基因以确定研究方向是一个不错的选择,既可以降低研究风险和成本,又可以把数据挖掘的结果作为论文的一部分,可谓两全其美。
本文主要通过R、Python完成基因的筛选工作,之所以用这些晦涩的脚本,是因为只有从原始数据出发才能保证数据的可靠性。比如我们之前借用别人开发的软件下载TCGA数据会遇到突然使用不了的情况,该软件对临床数据的解析很混乱,搞得我们当时看得一脸懵逼;对官网cases数据和files数据理解不全面,导致数据的滥用,癌与癌旁数据选择不当等等。这些令人抓狂的问题回归到对原始数据的解析可以得到很好的解决。
1
下载Clinicai和RNA-seq数据
➀进入https://portal.gdc.cancer.gov网站➙搜索胃癌数据(TCGA-STAD),RNA-seq数据选择HTSeq-FPKM(Counts是未经处理的原始表达量,而FPKM和FPKM-UQ是两种处理方法得到的数据)➙将文件加入Cart。
➁点击Cart➙页面跳转到如下图所示的界面。点击Download➙选择Manifest即为下载引导文件,由于文件较大,需要使用官方提供的下载工具(使用方法:更改环境变量,在cmd中输入命令gdc-client download -m gdc_manifest.2018-07-02.txt)➙点击Sample Sheet➙下载的文件为RNA-seq文件名与Sample ID对应关系的表格。
Clinical数据下载与此类似,不再赘述。
2
原始数据的解析、合并与整理
1
每个临床样本存储在单一的xml文件中,首先需要将这些xml文件放入同一文件夹下。
2
使用python minidom模块解析所需的临床数据。
3
这些数据包括bcr_patient_barcode、gender、vital_status、days_to_deathdays_to_birth、neoplasm_histologic_grade、pathologic_stage、pathologic_T、pathologic_N、pathologic_M、person_neoplasm_cancer_status等指标。部分结果,如下表所示:
4
和临床样本一样,每一个样本表达量数据存储在单一的txt文件中,将RNA-seq的数据解压到同一文件中。
5
将不同文件的数据合并到一块,将文件名与Sample ID对应上,再将Ensemble ID转化为Gene Symbol。
6
部分结果如下表所示:
3
基因筛选
生存分析,以基因表达量的中位数为分界点,将样本分为高低表达组,进行生存分析。
筛选出p值小于0.05的基因,部分结果如表所示:
不同预后指标分组下对各个基因表达量进行独立t检验。
共有7个指标,每个指标包括各组均值及p值,部分结果如下表所示。
Sample_ID前三个字段代表Case_ID,若相同即为同一Case,第四字段01代表癌组织,11代表正常组织。我们用同一样品的癌与癌旁组织做配对t检验。
部分结果如下表所示:
我们规定,IHC staining为高和中的数量大于等于2且IHC staining为低和未检出的数量也大于等于2的基因即为我们所需要的候选基因。
利用python检测生存曲线的交叉点,因为两条生存曲线开始会混杂在一起,所以实验从第20个像素开始检测,若红线全在蓝线上则为促癌基因,若蓝线均在红线上则为抑癌基因,若两者存在交叉则不予考虑。
4
结果
最后找到了37个基因,当然这个结果并不固定,我们可以适当放宽或提高筛选标准,保证有足够的基因进行后续的文献检索,以便能查找到未在客户所研究的癌症中被报道且在其他癌症中有研究的基因,进行接下的实验。
5
结语
“
通过编程可以大大降低重复性工作,使人从这些繁芜的工作中抽身出来思考更宏观的问题,同时能降低人犯错误的概率,以及人的主观性带来的标准不统一的问题。总之,因为客户的不断鞭策促使我们努力学习新方法、更新现有技术,用科学、严谨、创新的精神回报客户的每一份重托!
”
以上就是利用TCGA Clinical和RNA-Seq数据筛选候选基因的全部过程
ღ如有任何疑问,欢迎给小编留言ღ
♥期待您的来稿♥
领取专属 10元无门槛券
私享最新 技术干货