表达谱数据中相同基因如何处理

2022-09-21 16:26:56 浏览数 (1)

在分析表达谱芯片的时候,我们经常会遇到多个探针对应同一个基因的情况。一般遇到这种情况,最常见的两种处理方法是

1)取平均

2)取表达值高的那个探针

那么今天我们就用R来实现这两种处理方式。至于,如何将探针转换成相应的基因名字,相对来说还是比较容易的。一般的芯片数据都会有一个相应的注释文件,从中可以找到探针对应的基因名字。对于一些Agilent的商用芯片和一些比较特殊的芯片平台,可能找不到探针的注释文件。前面我们也简单介绍过

☞探针注释文件中没有基因名字怎么办?

☞探针注释文件中没有基因名字怎么办?(二)

首先我们先来随便造一个基因名有重复的表达谱数据。

代码语言:javascript复制
#设置随机过程的seed,保证结果可重复
set.seed(123)
#随机生成一个30行10列的矩阵
expr=matrix(runif(300,5,10),ncol=10)
#列名字为sample1-10
colnames(expr)=paste0("sample",1:10)
#行名从26个大写字母里面有放回的抽取30个字母,作为基因名
genes=sample(LETTERS,30,replace=T)
#合并得到基因名有重复的表达谱矩阵
expr=data.frame(genes,expr)
expr

接下来我们先用第一种方法

1)取平均

代码语言:javascript复制
#利用aggregate函数,对相同的基因名按列取平均
expr_mean=aggregate(.~genes,mean,data=expr)
expr_mean

会得到如下结果,感兴趣的小伙伴可以随便挑几个check一下

2)对于重复的基因名字,取表达值最大的哪一行

其实aggregate也可以对相同的基因使用max函数取最大值,但是这样处理是有问题的。例如同一个基因出现了三次,那么会有三行数据。如果使用aggregate max,对于每一个样本,他会从三个值中挑选最大的那个值最为这个样本的表达值,这样做是不科学的。我们先来看看效果

代码语言:javascript复制
#利用aggregate函数,对相同的基因名按列取取最大值
expr_max=aggregate(.~genes,max,data=expr)
expr_max

原始数据

处理之后的数据

所以这个做法不可取。

对于相同的基因,我们应该挑选行平均值大的那一整行,而不应该打乱

代码语言:javascript复制
#计算行平均值,按降序排列
index=order(rowMeans(expr[,-1]),decreasing = T)
#调整表达谱的基因顺序
expr_ordered=expr[index,]
#对于有重复的基因,保留第一次出现的那个,即行平均值大的那个
keep=!duplicated(expr_ordered$genes)
#得到最后处理之后的表达谱矩阵
expr_max=expr_ordered[keep,]
expr_max

最后结果是这样的

0 人点赞