有小伙伴们建议把Cytoscape 讲一讲，那我们今天我们先来学习共表达网络的第一部分：如何筛选共表达基因对！

照例回顾

MEV 聚类—助力转录组分析，成就大paper!

这篇算是转录组后期分析的一种思路！

好像木有其他相关文章，好吧，今天的就算是第一篇！
本人之前一直写基因组的文章，貌似近几年关注度不是特别高啊，反而转录组的，群体的关注度相当高，那我也分享些相关内容，来！

共表达网络定义

往往整个转录组中许多基因通常是同时表达的，表达趋势一致。那么问题来了，什么基因倾向于共表达呢？为什么要共表达呢？如何找到趋势一致的基因呢？（看回顾中的聚类文章）。

共表达的基因往往共同调控一类代谢过程，或者被某些转录因子共调控，所以有些共表达的基因启动子区往往有共同的启动子。或者某些共表达的基因可以形成某个调控网络，来行使一定的功能！

所以，有意思了，怎么做共表达网络呢。很简单，其实就是看表达相关性！

如何计算相关性

方案一

来看数据（test.data）：

第一列是基因名称，从第二列开始都是不同的样品，这是测试数据

来看脚本，我这里恰巧有一个脚本来做这个数据：

a=read.table("test.data",head=T)
library(Hmisc)
mat=matrix(ncol=4,nrow=sum(1:(ncol(a)-1)))
m=1
for(i in 1:(ncol(a)-1)){
     for(j in (i+1):ncol(a)){
           mat[m,2]=names(a)[j]
               mat[m,1]=names(a)[i]
               mat[m,3]=cor(a[,i],a[,j],method="spearman")
               w=rcorr(a[,i],a[,j],type="spearman")
               mat[m,4]=w$P[1,2]
          m=m+1
  }    
}

colnames(mat)

这是R脚本，可以将输入文件的名字修改成自己的文件名字，输出文件呢，就是p-value.xls，也可以自己修改！

来看输出文件：

这个结果呢，就将所有的基因对遍历了一遍，并且求得了他们的相关性和P 值，当然呢，这个P值有点偏大，因为这些数据都是我自己造的，大家可以用自己的数据做一下。

没错！
这个脚本用到了  Hmisc， 什么！还需要下载！ 懒人看下面的方案：

方案二

a=read.table("FPKM_log2.txt",head=T)
library(Hmisc)
mat=matrix(ncol=4,nrow=sum(1:(ncol(a)-2)))
m=1
for(i in 2:(ncol(a)-1)){
     for(j in (i+1):ncol(a)){
           mat[m,2]=names(a)[j]
               mat[m,1]=names(a)[i]
               mat[m,3]=cor(a[,i],a[,j],method="spearman")
               w=cor.test(a[,i],a[,j],type="spearman")
               mat[m,4]=w$p.value
          m=m+1
  }    
}

colnames(mat)

这个脚本呢，没有用到安装包，直接用自带的cor 函数计算的。

请到「今天看啥」查看全文