我以为我写完了R包终极解决方案！ 之后，应该不会再有任何关于R包安装的问题产生了，但仔细回过头来看才发现，我介绍的都是如何从CRAN或者bioconductor里面安装正规发布的包，但是有很多人开发的是自己私人的包，而我们有的确非常需要用怎么办？？这个时候就需要下载别人开发的包来安装了。比如我R包地址见github：https://github.com/jmzeng1314/humanid Continue reading
十一
29
Monthly Archives: 11月 2016
十一
29
如何开发自己的R包
随着R语言的流行度的提高，开发一个R包已经不再是专业程序猿才有的技能了。我这里讲的不是如何写一个包含了复杂统计公式或者发表一篇SCI文章的包，而是简简单单的用Rstudio自带的创建包的功能把自己的几个函数和数据打包！！！我R包地址见github：https://github.com/jmzeng1314/humanid Continue reading
十一
29
R来完成表达芯片分析全流程
包括如何从GEO下载数据，如何分组，两组直接如何找差异，差异基因如何去注释，包括GO/KEGG注释，还有特殊数据库，自定义数据库的注释，比如oncogene或者tumor suppress genes，TF的gene注释，还有GSEA软件的分析。
然后是对选择好的差异基因去string等PPI数据库拿到网络数据，在R或者cytoscape里面画网络图，然后是用MCODE插件和bioNet包来对网络找sub-network或者module，和hub genes。
就拿GSE42872 这个数据来做例子吧，希望听众具有基础R知识，了解什么是bioconductor，然后具有基础生物学知识，知道什么是基因，什么是表达，什么是通路，什么是富集，什么是注释。
总共10讲，每次半小时，每周3,4,6的晚上十一点开讲！
讲义的草稿如下，如果你能看懂草稿，能自己学会，就不用参加本次课程啦。
如果需要问我问题，就赶快找我申请加入交流群，提供本次培训的全套视频和代码！！！
十一
28
R语言画网络图三部曲之igraph
经过热心的小伙伴的提醒，我才知道我以前写的R语言画网络图三部曲竟然漏掉了最基础的一个包，就是igraph，不了解这个，后面的两个也是无源之水。
R语言画网络图三部曲之networkD3
R语言画网络图三部曲之sna
其实包括了3个包：igraph/RBGL/Rgraphviz
用到了一个测试数据，是构建好的PPI网络对象：We will first analyse a curated data set of protein-protein interactions in the yeast Saccharomyces cerevisiae extracted from published papers. This data set comes from with an R package called “yeastExpData”, which calls the data set “litG”. This data was first described in a paper by Ge et al (2001) in Nature Genetics (http://www.nature.com/ng/journal/v29/n4/full/ng776.html).
十一
25
hisat2+stringtie+ballgown
早在去年九月，我就写个博文说 RNA-seq流程需要进化啦！ http://www.bio-info-trainee.com/1022.html ，主要就是进化成hisat2+stringtie+ballgown的流程，但是我一直没有系统性的讲这个流程，因为我觉真心木有用。我只用了里面的hisat来做比对而已！但是群里的小伙伴问得特别多，我还是勉为其难的写一个教程吧，你们之间拷贝我的代码就可以安装这些软件的！然后自己找一个测试数据，我的脚本很容易用的！ Continue reading
十一
25
用BioNet这个bioconductor包来找 maximal-scoring subgraph
## 此包是为了解决一个难题： maximal-scoring subgraph (MSS) problem ，在一个巨大的复杂网络里面找到significantly differentially expressed subnetworks，就是说，得到了几百个差异基因，去PPI数据库做网络图的时候，发现还是巨大无比，所以需要用这个包来精简我们的网络图。
heuristically的中文意思：启发性地
## 而这个R包可以整合多种数据结果来给一个网络打分，
它整合了PPI网络分析和寻找功能模块的需求。
重点就是根据一个"igraph" or "graphNEL"对象和打分来找最大的MSS
subnet <- subNetwork(dataLym$label, interactome)
module <- runFastHeinz(subnet, scores)
plotModule(module, scores=scores, diff.expr=logFC) #这个就是精简后的我们的网络图。
其实另外一个函数也有类似的功能，dNetFind https://rdrr.io/cran/dnet/man/dNetFind.html
十一
25
最终还是把博客的全半角中英文标点符合的bug解决了
已经有非常多的小伙伴跟我反映了直接拷贝我的代码无法运行的问题，其实报错的原因很简单，就是中英文标点的bug而已。所以我给大家的理由是不用那么懒，拷贝我的代码，我就是希望你们能手动敲打每一个命令，来熟练记忆使用。
其实，我没那么好心，我就是懒而已。因为这个博客是host在阿里云的免费服务器上面的，各种IP密码我懒得去记忆，就差不多忘记了。当初弄好了我就懒得管了，正好现在博客免费时期快到了，也就想把这个问题解决掉。
我很简单搜索了一下，http://shiyun1013.blog.163.com/blog/static/10774036201301824446708/ 需要连接我博客的ftp，去修改博客里面的文件，
<?php remove_filter('the_content','wptexturize'); ?> 正好看看这个标点符号被改变了吗？
好像还不错，以后大家就可以直接拷贝我的代码去运行啦！
下面是我登陆了ftp，发现以前用rmarkdown写的几个教程，感兴趣的小伙伴可以随便看看！
十一
25
kegg在线链接图的颜色设置
一般来说， 有了kegg的ID，就可以直接去官网查看具体的通路图片，但是需要把差异基因给标注上去，就有点麻烦了，我以前做过类似的工作，结果没有做笔记，这次相当于重新造了个轮子，好惨！
简单的KEGG图片，看下面的url：
可以看出来只需要变化KEGG的ID即可。
如果要做下面的这个，上调基因用红色表示，下调基因用绿色表示：
十一
24
cytoscape五步曲之二：在cytoscape里面生成网络图
通过上一讲大家应该明白了，网络图是为了展现分子之间的连接关系的，并不是一定要用cytoscape来做，只需要根据连接关系给我们的所有点安排一个坐标，然后把相应的线连接起来即可！那么既然我们要学习cytoscape，肯定是要用cytoscape做好第一步，就是根据输入数据来做网络图。
可以先了解一下cytoscape定义好的输入数据，
http://wiki.cytoscape.org/Cytoscape_User_Manual/Network_Formats 当然，其实木有意义！因为我们不可能拿到cytoscape的输入文件（cys格式的），除非是你朋友传给你的。我们肯定是根据txt.csv等分割的文本文件来做网络图。
十一
24
cytoscape五步曲之一：明白什么是网络图
想了想还是写一个系列教程吧，问的朋友也太多了，主要是因为cytoscape跟python一样，经历了从2到3的进化阵痛过程，而且进化的面目全非了！！！很多人拿着2.x的说明书教程，视频，然后下载的却是3.x版本的cytoscape，真可怕！！！
已经从两万个芯片探测到的基因里面找到了近千个差异基因了，对它们做了GO/KEGG分析还是抓不住重点，看到文献说可以用PPI数据库做network analysis之后找hub基因，也也许可以说明一些问题！
提到 network analysis ，我想起来我以前总结过 R语言画网络图的三部曲，里面讲到过网络分析的基本原理！
十一
24
cytoscape五步曲之三：安装各种插件
软件安装我就不多说了，直接去官网下载即可，请务必下载3.x版本，我讲的是 最新版教程！
本次讲解如何给cytoscape安装插件，cytoscape本身是一个平台，学者可以在上面开发各种各样功能的插件实现不同的分析需求，类似于R语言这个平台，人们在上面安装包一样。R里面如何安装包我博客讲了4次，基本上看完的人都会懂。而cytoscape不一样，它的插件安装非常简单！非常简单！非常简单！
你只需要去cytoscape的APP中心找到包，如果你打开了cytoscape的界面，那么网页就会有install的字样，非常显眼，点击就自动安装了，这个时候会安装到
C:\Users\jimmy1314\CytoscapeConfiguration\3\apps\installed 这个目录！！~ 在你的电脑里面 jimmy1314 不一样
如果你这个时候并没有打开cytoscape的界面，那么网页就会有download的字样，也是非常显眼，点击就可以下载， 下载之后你需要自己把下载的jar文件放到cytoscape的安装路径，一般默认是
C:\Program Files\Cytoscape_v3.3.0\apps
最后，cytoscape提供了APP中心，就跟苹果手机安卓手机安卓软件一样，直接在cytoscape软件的菜单栏app中心就可以点击安装！
一般来说每一个插件都会对应的文章发表，而且一定在cytoscape的APP中心有展示，如果没有，可能是冒牌的哦！比如下面的插件
我要说的就是这么多了，我安装了十几个插件了，都没有什么问题，如果大家有遇到安装不了的，随时报告我，我来更新教程！联系jmzeng1314@163.com
下面的链接选择性观看：
十一
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进行画图
十一
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个：
十一
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格式)，比如下面：
十一
15
htseq-counts跟bedtools的区别
我以前写过bedtools和htseq-counts的教程，它们都可以用来对比对好的bam文件进行计数，正好群里有小伙伴问我它们的区别，我就简单做了一个比较，大家可以先看看我以前写的软件教程。写的有的挫：
言归正传，我这里没精力去探究它们的具体原理，只是看看它们数一个read是否属于某个基因的时候，区别在哪里，大家看下图： Continue reading
十一
14
TPM值就是RPKM的百分比嘛！
很久以前就有人问过这个问题啦，虽然目前主流还是用RPKM/FPKM来形容一个基因的表达量。但是既然大家都说TPM更好，我也来探究一下吧！
我不喜欢看公式，直接说事情，我有一个基因A，它在这个样本的转录组数据中被测序而且mapping到基因组了 5000个的reads，而这个基因A长度是10K，我们总测序文库是50M，所以这个基因A的RPKM值是 5000除以10，再除以50，为10. 就是把基因的reads数量根据基因长度和样本测序文库来normalization 。 Continue reading