R-基因表达的Limma包装



i具有以下格式的RMA归一化矩阵:

ID_REF GSM362180    GSM362181  GSM362188    GSM362189  GSM362192
244901 5.094871713 4.626623079 4.554272515 4.748604391 4.759221647
244902 5.194528083 4.985930299 4.817426064 5.151654407 4.838741605
244903 5.412329253 5.352970877 5.06250609  5.305709079 8.365082403
244904 5.529220594 5.28134657  5.467445095 5.62968933  5.458388909
244905 5.024052699 4.714631878 4.792865831 4.843975286 4.657188246
244906 5.786557533 5.242403911 5.060605782 5.458148567 5.890061836

不同的列对应于四种不同类型的启动子,而四个启动子中的每一个都具有生物学重复,因此完全有8列。

我尝试使用limma软件包在几个启动子(带有重复)中找到差异表达的基因,并且我总是会遇到一个错误,因为我是新手,无法完全理解它。

这是我使用的代码:

Group <- factor(c("p1", "p1", "p2", "p2", "p3","p3","p3","p4","p4"), levels = c("GSM362180","GSM362181","GSM362188","GSM362189","GSM362192","GSM362193","GSM362194","GSM362197","GSM362198"))
design <- model.matrix(~0 + Group)
colnames(design) <- c("GSM362180","GSM362181","GSM362188","GSM362189","GSM362192","GSM362193","GSM362194","GSM362197","GSM362198")
fit <- lmFit(modified, design)

在上述格式中,修改为RMA归一化数据矩阵。我收到以下错误:

Coefficients not estimable: GSM362180 GSM362181 GSM362188 GSM362189 GSM362192 GSM362193 GSM362194 GSM362197 GSM362198 

lmfit中的错误(设计,t(m)):0(non-na)案例

应该只是:

Group <- factor(c("p1", "p1", "p2", "p2", "p3", "p3", "p3", "p4", "p4"))

组类别顺序应对应于colnames(modified)

中样本的顺序

在原始代码中,您将不同的唯一levels分配给Group中的每个因素,因此从本质上讲,您没有分组类别,没有重复,也无法估计系数。

创建您的设计矩阵:

design <- model.matrix(~0 + Group)

然后,您想将样本的名称分配给设计矩阵的rownames,而不是colnames,因为列代表分组。

rownames(design) <- c("GSM362180","GSM362181","GSM362188","GSM362189","GSM362192","GSM362193","GSM362194","GSM362197","GSM362198")
# or simply
rownames(design) <- colnames(modified)

您的设计矩阵应该看起来像这样,带有一个 1,其中一个样本属于特定组类别。

design
          Groupp1 Groupp2 Groupp3 Groupp4
GSM362180       1       0       0       0
GSM362181       1       0       0       0
GSM362188       0       1       0       0
GSM362189       0       1       0       0
GSM362192       0       0       1       0
GSM362193       0       0       1       0
GSM362194       0       0       1       0
GSM362197       0       0       0       1
GSM362198       0       0       0       1
attr(,"assign")
[1] 1 1 1 1
attr(,"contrasts")
attr(,"contrasts")$Group
[1] "contr.treatment"

最新更新