0

我有以下格式的 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 包来查找跨多个启动子（具有重复）的差异表达基因，但由于我是 r 新手并且无法完全理解它，所以我总是遇到错误。

这是我使用的代码：

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)

其中 modified 是以上述格式输入的 rma 归一化数据矩阵。我收到以下错误：

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

lmfit（设计，t（M））中的错误：0（非NA）案例

4

1 回答 1

0

应该只是：

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"
于 2017-10-27T10:58:54.827 回答