2

所以我使用 r 和包 Bioconductor (oligo), (limma) 来分析一些微阵列数据。

我在配对分析中遇到了麻烦。

所以这是我的表型数据

ph@data ph@data index filename group WT1 WT WT1 WT WT2 WT WT2 WT WT3 WT WT3 WT WT4 WT WT4 WT LT1 LT LT1 LT LT2 LT LT2 LT LT3 LT LT3 LT LT4 LT LT4 LT TG1 TG TG1 TG TG2 TG TG2 TG TG3 TG TG3 TG TG4 TG TG4 TG

所以为了分析我做了这个代码:

colnames(ph@data)[2]="source"
ph@data
groups = ph@data$source
f = factor(groups,levels = c("WT","LT","TG"))
design = model.matrix(~ 0 + f) 
colnames(design)=c("WT","LT","TG")
design
data.fit = lmFit(normData,design)

> design
   WT LT TG
1   1  0  0
2   1  0  0
3   1  0  0
4   1  0  0
5   0  1  0
6   0  1  0
7   0  1  0
8   0  1  0
9   0  0  1
10  0  0  1
11  0  0  1
12  0  0  1
attr(,"assign")
[1] 1 1 1
attr(,"contrasts")
attr(,"contrasts")$f
[1] "contr.treatment"

然后我起诉这个来比较我的组:

contrast.matrix = makeContrasts(LT-WT,TG-WT,LT-TG,levels=design) data.fit.con = contrasts.fit(data.fit,contrast.matrix) data.fit.eb = eBayes(data.fit.con)

所以在这一切之后,我想比较我的组:

ph@data[ ,2] = c("control","control","control","control","littermate","littermate","littermate","littermate","transgenic","transgenic","transgenic","transgenic")
colnames(ph@data)[2]="Treatment"
ph@data[ ,3] = c("B1","B2","B3","B4","B1","B2","B3","B4","B1","B2","B3","B4")
colnames(ph@data)[3]="BioRep"
ph@data```

groupsB = ph@data$BioRep 
groupsT = ph@data$Treatment
fb = factor(groupsB,levels=c("B1","B2","B3","B4"))
ft = factor(groupsT,levels=c("control","littermate","transgenic"))```

paired.design = (~ fb + ft)

所以现在我的 phenoData 看起来像这样:

    index  Treatment BioRep
WT1    WT    control     B1
WT2    WT    control     B2
WT3    WT    control     B3
WT4    WT    control     B4
LT1    LT littermate     B1
LT2    LT littermate     B2
LT3    LT littermate     B3
LT4    LT littermate     B4
TG1    TG transgenic     B1
TG2    TG transgenic     B2
TG3    TG transgenic     B3
TG4    TG transgenic     B4```

B1 是生物复制,然后我有野生型、同窝仔和转基因组

比较我的样品我正在尝试这个

colnames(paired.design)=c("Intercept","B4vsB1","B3vsB1","B2vsB1","B4vsB2","B3vsB2","B4vsB3","littermatevscontrol","transgenicvscontrol")

但后来我得到了这个错误: Error in `colnames<-`(`*tmp*`, value = c("Intercept", "WTvsLT", "WTvsTG", : attempt to set 'colnames' on an object with less than two dimensions

我做错了什么,这是比较我的数据的正确方法?

data.fit = lmFit(data.rma,paired.design)
data.fit$coefficients
4

1 回答 1

0

因为我没有你的数据,所以我做了 ft 和 fb:

fb <- structure(c(1L, 2L, 3L, 4L, 1L, 2L, 3L, 4L, 1L, 2L, 3L, 4L), .Label = c("B1", 
"B2", "B3", "B4"), class = "factor")

ft <- structure(c(1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L), .Label = c("control", "littermate", "transgenic"), class = "factor")

这一行有一个错字:

paired.design = (~ fb + ft)
dim(paired.design)
NULL

它应该是:

paired.design = model.matrix(~ fb + ft)
head(paired.design)
  (Intercept) fbB2 fbB3 fbB4 ftlittermate fttransgenic
1           1    0    0    0            0            0
2           1    1    0    0            0            0

如果您查看您的 model.matrix,它有 6 列,这不是您要分配的。因此,当您指定 时~ fb + ft,选择因子的第一个水平作为参考,您会发现其他水平对此的影响。例如,对于 fb,B1 是参考,您估计 B2、B3、B4 的影响。

要进行所需的成对比较,对于所有内容与参考,您只需指定列本身,例如,对于 B4 与 B1,它将只是 fbB4,因为 B4 是相对于 B1 估计的。对于其他人,您可以像以前一样指定“-”:

contrast.matrix = makeContrasts(fbB4,fbB3,fbB2,fbB4-fbB2,
fbB3-fbB2,fbB4-fbB3,ftlittermate,fttransgenic,levels=paired.design)

现在我们可以根据需要分配 col 名称:

colnames(contrast.matrix) <-c("B4vsB1","B3vsB1","B2vsB1","B4vsB2","B3vsB2","B4vsB3",
"littermatevscontrol","transgenicvscontrol")

我为 normData 模拟了一些数据以适应模型和对比:

set.seed(100)
normData = matrix(rpois(100*12,100),100,12)
data.fit = lmFit(normData,paired.design)
data.fit.con = contrasts.fit(data.fit,contrast.matrix)

注意,有一个警告:

Warning message:
In contrasts.fit(data.fit, contrast.matrix) :
  row names of contrasts don't match col names of coefficients

因为“(Intercept)”改名为Intercept。

您可以查看 data.fit.con 下的系数以查看所做比较的倍数变化

  Contrasts
           B4vsB1    B3vsB1      B2vsB1     B4vsB2     B3vsB2     B4vsB3
  [1,]  14.333333 11.000000   5.0000000   9.333333   6.000000   3.333333
  [2,]  -3.333333  2.000000   7.0000000 -10.333333  -5.000000  -5.333333
  [3,]   3.666667 -2.666667   6.3333333  -2.666667  -9.000000   6.333333
  [4,] -22.666667 -1.666667 -10.3333333 -12.333333   8.666667 -21.000000
  [5,]  -8.333333 -3.666667   8.3333333 -16.666667 -12.000000  -4.666667
  [6,]  15.333333 -5.666667  -0.3333333  15.666667  -5.333333  21.000000
      Contrasts
       littermatevscontrol transgenicvscontrol
  [1,]                1.25               -7.50
  [2,]                3.50               10.50
  [3,]               -3.75                2.75
  [4,]                7.75               14.00
  [5,]               16.75               -2.50
  [6,]                2.00                9.50

您可以与原始拟合系数进行交叉检查,以查看是否产生了正确的对比度:

head(data.fit$coefficients)
     (Intercept)        fbB2      fbB3       fbB4 ftlittermate fttransgenic
[1,]    95.41667   5.0000000 11.000000  14.333333         1.25        -7.50
[2,]    94.33333   7.0000000  2.000000  -3.333333         3.50        10.50
[3,]    97.66667   6.3333333 -2.666667   3.666667        -3.75         2.75
[4,]    93.41667 -10.3333333 -1.666667 -22.666667         7.75        14.00
[5,]    98.91667   8.3333333 -3.666667  -8.333333        16.75        -2.50
[6,]    97.16667  -0.3333333 -5.666667  15.333333         2.00         9.50
于 2019-11-28T22:41:00.687 回答