我正在用一种方法来找到甲基化调控的基因。该方法首先计算样品中每个基因甲基化与表达数据之间的相关系数。我应用这个第一个过滤器,我得到了418个基因的子集,其相关系数为显著(小于-0,3)。
下一步是找出甲基化和每个基因表达的b样条三次回归系数。我的数据包含在一个包含两个矩阵的列表中,如下所示:
> met.spl[1:5,1:5]
AAAS ABCA7 ABCC12 ABCF1 ACN9
paciente6 0.15013343 0.01489557 0.7987748 0.02826255 0.02169665
paciente7 0.10827520 0.01215497 0.8515188 0.03378193 0.02452192
paciente8 0.12423923 0.01682172 0.4182180 0.03288906 0.02046130
paciente9 0.11779075 0.02198105 0.6101996 0.06389504 0.04574667
paciente10 0.09234602 0.01526621 0.8366319 0.02868425 0.02095470
> exp.spl[1:5,1:5]
AAAS ABCA7 ABCC12 ABCF1 ACN9
paciente1 -0.82350 -1.20725 0.6000000 -0.6783 0.64500
paciente2 -1.14075 -0.59575 0.2173333 -0.2644 0.54100
paciente3 -1.43000 -1.72750 1.0015000 -1.1413 0.98625
paciente4 -1.16650 -0.76250 0.4378333 -0.6804 -0.58650
paciente5 -0.51125 -1.10325 -0.1231667 -0.1521 0.02750为了计算第一个基因的模型系数,我做了以下步骤:
> lm(exp.spl[,1] ~ bs(met.spl[,1]), data = Lshape.spline)$coef
(Intercept) bs(met.spl[, 1])1 bs(met.spl[, 1])2 bs(met.spl[, 1])3
-1.00616163 -0.44292576 0.08767607 -0.61237162 我的目标是为418个选定的基因(行)生成一个包含系数(4列)的新对象。
发布于 2015-06-24 15:37:14
似乎一个简单的sapply会在这里工作
library(splines)
sapply(1:ncol(exp.spl), function(i) {
lm(exp.spl[,i] ~ bs(met.spl[,i]))$coef
})测试用
met.spl<-read.table(text="AAAS ABCA7 ABCC12 ABCF1 ACN9
paciente6 0.15013343 0.01489557 0.7987748 0.02826255 0.02169665
paciente7 0.10827520 0.01215497 0.8515188 0.03378193 0.02452192
paciente8 0.12423923 0.01682172 0.4182180 0.03288906 0.02046130
paciente9 0.11779075 0.02198105 0.6101996 0.06389504 0.04574667
paciente10 0.09234602 0.01526621 0.8366319 0.02868425 0.02095470", header=T)
exp.spl<-read.table(text="AAAS ABCA7 ABCC12 ABCF1 ACN9
paciente1 -0.82350 -1.20725 0.6000000 -0.6783 0.64500
paciente2 -1.14075 -0.59575 0.2173333 -0.2644 0.54100
paciente3 -1.43000 -1.72750 1.0015000 -1.1413 0.98625
paciente4 -1.16650 -0.76250 0.4378333 -0.6804 -0.58650
paciente5 -0.51125 -1.10325 -0.1231667 -0.1521 0.02750", header=T)https://stackoverflow.com/questions/31030481
复制相似问题