Dear Members,

I have this dataset extracted from GEO database (GSE16962). According to the data information provided below, I designed a matrix for lmFit. I am not sure if it is correct. Could somebody take a quick look?

> show(pData(phenoData(gse[[1]]))[1:12,c(1,6,8)])
                                title type                                                          source_name_ch1
GSM424759 pSUPER-scramble replicate 1  RNA HUVEC infected by retroviral vectors bearing a control scramble sequence
GSM424760 pSUPER-scramble replicate 2  RNA HUVEC infected by retroviral vectors bearing a control scramble sequence
GSM424761 pSUPER-scramble replicate 3  RNA HUVEC infected by retroviral vectors bearing a control scramble sequence
GSM424762  pSUPER-mir-210 replicate 1  RNA        HUVEC infected by retroviral vectors bearing pre-miR-210 sequence
GSM424763  pSUPER-mir-210 replicate 2  RNA        HUVEC infected by retroviral vectors bearing pre-miR-210 sequence
GSM424764  pSUPER-mir-210 replicate 3  RNA        HUVEC infected by retroviral vectors bearing pre-miR-210 sequence
GSM424765   Anti-scramble replicate 1  RNA                                    HUVEC transfected with a scramble-LNA
GSM424766   Anti-scramble replicate 2  RNA                                    HUVEC transfected with a scramble-LNA
GSM424767   Anti-scramble replicate 3  RNA                                    HUVEC transfected with a scramble-LNA
GSM424768    Anti-mir-210 replicate 1  RNA                                  HUVEC transfected with anti-miR-210-LNA
GSM424769    Anti-mir-210 replicate 2  RNA                                  HUVEC transfected with anti-miR-210-LNA
GSM424770    Anti-mir-210 replicate 3  RNA                                  HUVEC transfected with anti-miR-210-LNA
>treatments=factor(c(1,1,1,2,2,2,3,3,3,4,4,4), labels=c("CTRL","mir210","CTRL1","CTRL2"))
> design=model.matrix(~treatments)
> colnames(design)=c("CTRL","mir210","CTRL1","CTRL2")
> design
   CTRL mir210 CTRL1 CTRL2
1     1      0     0     0
2     1      0     0     0
3     1      0     0     0
4     1      1     0     0
5     1      1     0     0
6     1      1     0     0
7     1      0     1     0
8     1      0     1     0
9     1      0     1     0
10    1      0     0     1
11    1      0     0     1
12    1      0     0     1
>library(limma)
> fit=lmFit(eset,design)
> fit=eBayes(fit)
Warning message:
Zero sample variances detected, have been offset
> results=classifyTestsF(fit,p.value=0.0001)
> summary(results)
     CTRL      mir210   CTRL1    CTRL2
-1     0           1364        733       961
0      0          45192     48699   48787
1  54675       8119       5243    4927

###mir210 and CTRL2 are the two arrays that I really care about. CTRL and CTRL1 are experiment control. CTRL2 may have some impact to true result.
So I decided to make following cont.contrasts

>cont.contrasts=makeContrasts(mir210=c(0,1,0,0),CTRL2=c(0,0,0,1),levels=design)
> cont.contrasts
        Contrasts
Levels   mir210    CTRL2
  CTRL          0            0
  mir210      1            0
  CTRL1        0           0
  CTRL2        0           1
> fit2=contrasts.fit(fit,cont.contrasts)
> fit2=eBayes(fit2)
> results=classifyTestsF(fit2,p.value=0.0001)
> summary(results)
   mir210    CTRL2
-1    331       559
0   54151   53908                      It looks to me that this is the TRUE differentially expressed list for mir210 treatments.
1     193       208

> table(Cmir210=results[,1],CTRL2=results[,2])
                  CTRL2
Cmir210    -1          0     1
     -1           17     272    42
     0          450  53535   166
     1            92     101     0


Many Many thanks

Jing




	[[alternative HTML version deleted]]

