十一 23

关于multiple mapping我想说的

很多时候,我们都要选取unique mapped的reads,尤其是在RNA-seq和CHIP-seq的时候,但是如何保留,各种教程都不一致,我稍微总结了一下,是因为使用的比对工具不一样导致的!但是主要都反应在sam文件的一系列tag里面~

首先对bwa来说,如果它遇到一个reads可以比对到参考基因在的多个序列,只会随机的选取一个位置来输出到sam文件,但是会加上一个tag是XS:I:<N>来告诉我们第二好的比对情况的比对得分是多少,bowtie也是一样。但是它们都有参数来决定是否只对每个reads输出一条信息,还是输出全部的信息,在bwa是-a的参数,在bowtie里面是-m参数。

但是bowtie2里面取消了这个参数,它们都必须用XS:I:<N>这个tag来挑选unique mapped的reads

但是如果是用hisat来比对的话,决定是否是唯一比对的却是NH这个tag信息。默认情况下一条reads可以输出多条比对结果。

我想起了再补充吧,其实应该找几个例子用IGV看看,就明白了,可是我暂时没有时间了,只是觉得这个很重要,就提一下。

 

十一 23

quantile normalization到底对数据做了什么?

提到normalization很多人都烦了,几十种方法,而对于芯片或者其它表达数据来说,最常见的莫过于quantile normalization啦。那么它到底对我们的表达数据做了什么呢?首先要么要清楚一个概念,表达矩阵的每一列都是一个样本,每一行都是一个基因或者探针,值就是表达量咯。quantile normalization 就是对每列单独进行排序,排好序的矩阵求平均值,得到平均值向量,然后根据原矩阵的排序情况替换对应的平均值,所以normalization之后的值只有平均值了。具体看下面的图: Continue reading

十一 23

用R的bioconductor里面的stringDB包来做PPI分析

PPI本质上是根据一系列感兴趣的蛋白质或者基因(可以是几百个甚至上千个)来去PPI数据库里面找到跟这系列蛋白质或者基因的相互作用关系!

本次的主角是stringDB,顾名思义用得是大名鼎鼎的string数据库,
本来还以为需要自己上传自己的基因给这个数据库去做分析,没想到他们也开发了R包,主页见: http://www.bioconductor.org/packages/release/bioc/html/STRINGdb.html 而我比较喜欢用编程来解决问题,所以就学了一下这个包,非常好用!
它只需要一个3列的data.frame,分别是logFC,p.value,gene ID,就是标准的差异分析的结果。
然后用string_db$map函数给它加上一列是 string 数据库的蛋白ID,然后用string_db$add_diff_exp_color函数给它加上一列是color。
用string_db$plot_network函数画网络图,只需要 string 数据库的蛋白ID,如果需要给蛋白标记不同的颜色,需要用string_db$post_payload来把color对应到每个蛋白,然后再画网络图。
也可以直接用get_interactions函数得到所有的PPI数据,然后写入到本地,再导入到cytoscape进行画图

Continue reading

十一 23

java环境变量的问题

有篇文章提到了cytoscape,想着一直没用过这个神器对不起我生信大神的称号呀,就下载了准备安装,居然报错了,简直不可思议,因为一直以为它是java软件,一般不需要安装,结果是exe的,只是依赖于java,报错是EXE4J_JAVA_HOME, No JVM could be found on your system,这是个很常见的错误,我 简单搜索了解决方案https://wincrunch.com/exe4j-java-home-no-jvm-could-be-found-on-your-system/ 居然无效,但是里面有句话引起了我的注意,通常64位的window电脑的java是安装在Program Files 而不是Program Files (x86),这才是问题所在,我当初图简单,直接用了JDK来安装JRE,所以导致软件安装目录错误。有非常多的生物信息学软件都依赖与java,比如IGV,GSEA,cytoscape,一般来说window电脑安装好了java之后这些软件都挺好用的。那么关于java问题,我整理了3个:

Continue reading

十一 23

【直播】我的基因组(八):原始测序数据质量报告

由于我是分期付款,所以我先拿到了我的测序数据的质控结果和比对情况分析报告,需要补齐全款后才能拿到原始测序数据!(中间还出了个小意外,打款的时候不小心多打了30块钱!(⊙o⊙)…不过多打的30块钱想拿回来估计不太可能了,需要填写书面申请表格并且自费快递到公司,这边跨境快递费都不止这个数了) Continue reading

十一 23

【直播】我的基因组(七):从整体理解全基因组测序数据的变异位点

首先记住一个很重要的知识点,变异是相对的!

简单说一下什么是找变异,变异跟突变有什么区别呢?举个栗子:有国际组织规定了人类的参考基因组(如UCSC,ENSEMBL,NCBI等,前面帖子都有讲),就是 AAAAA(这里简化一下,就5个碱基,其实人类基因组多达30亿个)  。现在通过给自己测序得知,我与之对应的是AGCAA,那么我相比国际基因组来说,就是2个变异位点,位于基因组的坐标2和3,但是它们还不能说就是突变。 Continue reading

十一 23

【直播】我的基因组(六):变异位点注释数据库的准备

通常一个人的全基因组测序数据可以挖掘到四百万个SNVs(跟参考基因组不一样的单碱基位点),还有五十万的indels(insertions or deletions),但是得到的数据通常是以vcf文件格式给出的(自行搜索什么是vcf格式),比如下面:

很明显,正常人是看不懂这些变异位点有啥子一样的,只知道第20条染色体的1230237坐标上面本来是一个T碱基的,但是突变成了G,那么我们必然还想知道,这个位点是在某个基因上面吗?如果是,在基因的外显子还是内含子?它的突变有没有改变该基因的功能呢?有没有影响它的转录和翻译呢?还有世界上有没有其他正常人也是这个位点变异呢?如果有,是哪些人种呢?有没有癌症病人也发现了这个变异呢?如果有,是什么癌症呢?所以我们必须下载一系列的变异位点注释数据库,来全方位的解释我们自己找到那四百万个SNVs和五十万的indels。下面我们一起进行数据库准备。

Continue reading