## Experiment with illumina microArrays 96 samples. There are 3 groups (X,Y,Z) each group has 96 samples ## Loading all the required libraris library(lumi) library(limma) library(mgcv) library(lumiHumanAll.db) library(lumiHumanIDMapping) library(GOstats) library(annotate) library(GO.db) library(illuminaHumanv3BeadID.db) ## Without converting to NuId data_96_raw <- lumiR("Sample_Probe_Profile.txt",detectionTh =0.05) summary(data_96_raw, 'QC') data_96_E_object <- lumiExpresso(data_96_raw) ## Identifying differentially expressed genes # Removing un-annotated genes data_96_Matrix <- exprs(data_96_E) data_96_probeList <- rownames(data_96_Matrix) x <- illuminaHumanv3BeadIDSYMBOL annotated_ids <- mappedkeys(x) idx<- data_96_probeList %in% annotated_ids data_96_Matrix<-data_96_Matrix[idx,] ## Creating a design matrix data_96_sampleType <- c("Y","Y","Y","Y","Y","Z","Z","Z","Y","X","Z","Z ","X","Y","X","Y","Y","Y","Y","Z","Y","Z","Y","Z","X","X","Z","Y","X", "Z","Y","Z","Z","X","X","Z","Z","Y","Y","Y","Z","Y","Z","Y","Y","X","Y ","X","Z","X","Y","Y","Z","X","Z","Z","Z","X","X","Y","Z","Z","X","Y", "Y","X","X","Z","Z","Z","X","X","Z","Y","X","X","X","Z","X","X","Y","Y ","X","Z","X","X","X","Z","X","Z","X","Y","Y","Z","Y","X") data_96_design <- model.matrix(~0+data_96_sampleType) colnames(data_96_design) <- c('X','Y','Z') Design Matrix X Y Z 1 0 1 0 2 0 1 0 3 0 1 0 4 0 1 0 5 0 1 0 6 0 0 1 7 0 0 1 8 0 0 1 9 0 1 0 10 1 0 0 11 0 0 1 12 0 0 1 13 1 0 0 . . . . 96 1 0 0 data_96_fit1 <- lmFit(data_96_Matrix,data_96_design) data_96.contrast <- makeContrasts (Y-X,Z-Y,Z-X,levels=data_96_design) Contrasts Matrix Levels Y - X Z - Y Z - X C -1 0 -1 I 1 -1 0 S 0 1 1 data_96_fit2 <- contrasts.fit(data_96_fit1,data_96.contrast) data_96_fit2 <- eBayes(data_96_fit2) topTable(data_96_fit2,coef=1, adjust="BH") result <- decideTests(p.value=0.05) vennDiagram(result) ********************************************************************** **** topTable(data_96_fit2,coef=1, adjust="BH") ID logFC AveExpr t P.Value adj.P.Val B 5517 1300397 0.11999169 6.341387 4.819393 5.453441e-06 0.1555921 3.8059087 6069 3870273 -0.31283278 6.881315 -4.555384 1.552881e-05 0.2215262 2.8294648 23026 5310079 -0.13036276 6.815443 -4.280789 4.449540e-05 0.2890393 1.8509379 17388 5390427 0.25070344 6.487644 4.185116 6.363124e-05 0.2890393 1.5194451 17387 5340240 0.21502145 6.470917 4.152405 7.183003e-05 0.2890393 1.4072650 13729 580364 0.13918530 6.616036 4.151444 7.208548e-05 0.2890393 1.4039803 10224 3460242 -0.17331384 7.322021 -4.148558 7.285831e-05 0.2890393 1.3941133 19894 4610626 0.05532714 6.224477 4.078810 9.414832e-05 0.2890393 1.1570851 4131 6130370 0.05371882 6.177268 4.066616 9.843716e-05 0.2890393 1.1159301 6762 4830255 -0.11130830 6.478843 -4.025402 1.143630e-04 0.2890393 0.9774658 topTable(data_96_fit2,coef=2, adjust="BH") ID logFC AveExpr t P.Value adj.P.Val B 1926 70767 -0.13534270 6.625839 -4.250721 0.0000498156 0.846945 1.74279904 197 730671 -0.12151636 7.248196 -3.968852 0.0001402750 0.846945 0.78510309 27388 6040500 -0.06651645 6.145904 -3.918356 0.0001680832 0.846945 0.61838206 6961 1110037 -0.05087288 6.238527 -3.889455 0.0001862930 0.846945 0.52364436 2994 7610202 -0.09374225 7.130537 -3.753356 0.0003004434 0.846945 0.08432137 20403 7210564 -0.05199980 6.220342 -3.742198 0.0003123001 0.846945 0.04880959 17651 4200692 -0.11728191 7.258627 -3.644029 0.0004376044 0.846945 -0.26025698 18244 2000458 -0.05455462 6.210703 -3.591672 0.0005226312 0.846945 -0.42258864 17655 110458 -0.05546275 6.227990 -3.574337 0.0005540729 0.846945 -0.47594401 9425 60239 -0.04479364 6.171825 -3.557422 0.0005864732 0.846945 -0.52782091 ********************************************************************** ******** Array Weigh of Raw Data :: 1 2 3 4 5 6 7 8 3.2861102 0.6091478 2.6927503 2.2778201 2.4391841 0.7802437 1.1378995 0.5714981 9 10 11 12 13 14 15 16 0.8091165 0.5260511 0.7565424 0.8354660 0.2929977 0.4950789 0.4512856 0.1283279 17 18 19 20 21 22 23 24 0.6713553 0.6275333 0.5441217 0.2711778 0.7278274 0.3323277 0.5786381 0.7375983 25 26 27 28 29 30 31 32 0.5793944 0.4018671 0.6658730 1.6805341 0.7210032 1.2802383 0.8902326 0.6715769 33 34 35 36 37 38 39 40 0.2959048 0.8854763 0.7587111 1.5067705 0.1766820 0.1639351 0.1739882 0.1833592 41 42 43 44 45 46 47 48 0.2443871 0.2188281 0.2945594 0.2607380 0.3419603 0.5932775 0.3396003 0.7229979 49 50 51 52 53 54 55 56 0.8274153 1.7944800 0.8306203 1.0738534 1.5052515 2.9798751 2.0836544 1.5010941 57 58 59 60 61 62 63 64 2.3827591 2.7575819 1.5967670 3.3397084 2.0576009 2.9660262 3.8417611 1.3469622 65 66 67 68 69 70 71 72 1.7334743 2.7855044 4.5417142 3.2551286 3.9634064 3.8054377 3.8020479 4.2732205 73 74 75 76 77 78 79 80 0.2659556 0.3440704 0.6404601 0.8697107 1.0736648 0.8835502 0.9912724 0.4334270 81 82 83 84 85 86 87 88 0.8801461 1.0309162 2.5352693 1.8366521 2.5501944 4.5388081 4.2603196 1.7665255 89 90 91 92 93 94 95 96 2.5511084 2.3867529 4.0237262 2.1647555 3.3101943 2.5863558 1.8618118 0.2397050 ********************************************************************** ********* Array Weight of normalized Data :: 1 2 3 4 5 6 7 8 0.8908085 1.1978748 1.1454450 1.0215452 1.3605699 0.8141813 1.4172765 1.2958270 9 10 11 12 13 14 15 16 1.0522608 1.2562065 1.3531495 1.1489293 1.0596229 1.1217643 0.9156484 0.6705476 17 18 19 20 21 22 23 24 0.9643412 1.2120423 1.0672521 0.9983735 0.8096782 0.9379575 1.0493924 0.7614746 25 26 27 28 29 30 31 32 1.1573255 1.2735353 1.3795986 0.9499952 1.3602425 1.2549726 1.1772558 1.4158351 33 34 35 36 37 38 39 40 1.3659316 1.0658774 1.3185647 0.9017821 0.7915704 0.6326567 1.0325512 0.6812818 41 42 43 44 45 46 47 48 1.0142990 1.1921695 1.1476346 0.8255798 1.2012711 1.1893762 0.9367947 1.1594407 49 50 51 52 53 54 55 56 1.1072266 1.0561352 0.8488650 1.1756689 1.0554286 0.9326934 1.2836555 0.9665945 57 58 59 60 61 62 63 64 1.3293081 1.3027377 1.3459009 1.2834152 0.9657114 1.0041629 0.8823282 0.6843356 65 66 67 68 69 70 71 72 0.8883121 0.9374396 0.9799943 0.9733364 1.2045700 1.0693506 0.7730626 0.9408430 73 74 75 76 77 78 79 80 0.6685933 1.0492289 0.9952222 0.9690452 1.0150076 1.0148188 0.6687192 0.5641711 81 82 83 84 85 86 87 88 0.7872238 0.8793660 0.8553373 0.9662603 0.6143880 1.0340769 0.9842437 0.6626941 89 90 91 92 93 94 95 96 1.0270788 0.8590296 1.0807271 0.8162435 1.0398548 0.8993595 1.2094240 0.5715857

Please look at some of my previous comments on the interpretation of FDR. If e.g. fewer than 1% of your genes are significant at p<.01, then the FDR is 100%. --Naomi At 09:01 AM 9/8/2009, Md.Mamunur Rashid wrote: >-------- Original Message -------- >Subject: result of linear model >Date: Tue, 8 Sep 2009 13:59:41 +0100 >From: Md.Mamunur Rashid <mamunur.rashid at="" kcl.ac.uk=""> >To: smyth at wehi.edu.au <smyth at="" wehi.edu.au="">, >bioconductor at stat.math.ethz.ch <bioconductor at="" stat.math.ethz.ch=""> > > > >Dear Gordon, > >I am working with a set of illumina microarray data (96 samples) >from three groups >(i.e. group-1(X),group-2(Y),group-3(Z)). 32 samples from each group. > >I have read the data using lumiR method and processed the data using >lumi Methods(lumiExpresso). > >Now I need to identify the differentially expressed genes by >comparing each of these groups >with each other. I am using linear model in limma package and >topTable method to identify >top N differentially expressed genes. Below is the code that I have >used to design the test >in linear model. If you look at the result of the top-table all the >adjusted p-values are very >high. as a result none of the genes passes through the cut-off >p-value 0.05. I have tried all >the tree combination and in all cases adjusted p-values are. > >Now I have tested the same code on another set of data that has 20 >samples that has been analyzed >before (i.e I have the results before hand). In that case also >adjusted p-values are very high but genes in >the top-table are correct( i.e matches with the result of the >previous analysis). > > >so in this situation, I will be really grateful if you can give me >some suggestion on > >1. Is there anything wrong in my code which is making the adjusted >p-value so high.??? >2. Or it might be a problem in the data pre-process phase.??? > >*** I have attached the code , the result and the array-weight in a text file. > Please have a look. > >thanks in advance. If anybody else want to have a suggestion, you >are most welcome. > >regards, >Md.Mamunur Rashid > > > >######## code ############# > > >## norm_object is the normalized object > >d_Matrix<- exprs(norm_object) >probeList<- rownames(d_Matrix) >library(illuminaHumanv3BeadID.db) > >## filtering out the un-annotated genes > >x<- illuminaHumanv3BeadIDSYMBOL >annotated_ids<- mappedkeys(x) >d_Matrix<- exprs(norm_object) >probeList<- rownames(d_Matrix) >idx<- probeList %in% annotated_ids >d_matrix<-d_matrix[idx,] > > >## 32 samples from each group without any pair >sampleType<- >c("X","X","Y","Y",..........96 >samples..........,"Z","X","Y","Y","X","X","Z","I","X") >design<-model.matrix(~0+sampleType) >colnames(design_norm_test)<- c('X','Y','Z') > >fit1<- lmFit(d_Matrix,design) >constrast.matrix<- makeContrasts (Y-X , Z-Y , Z-X, levels=design) >fit1_2<- contrasts.fit(fit1,contrast.matrix) >fit1_2<- eBayes(fit1_2) >topTable(fit1_2,coef=1, adjust="BH") > > > >Result of the toptable method: > >6284 ILMN_1111111 0.11999169 6.341387 4.828711 5.237786e-06 0.2325975 >12919 ILMN_2222222 -0.05966259 6.187268 -4.678886 9.532099e-06 0.2325975 >6928 ILMN_3333333 -0.31283278 6.881315 -4.561366 1.513503e-05 0.2462115 >42428 ILMN_4444444 -0.13036276 6.815443 -4.288051 4.321272e-05 0.3964163 >36153 ILMN_5555555 0.25070344 6.487644 4.190735 6.220719e-05 0.3964163 >36152 ILMN_6666666 0.21502145 6.470917 4.158153 7.019901e-05 0.3964163 >28506 ILMN_7777777 0.13918530 6.616036 4.158140 7.020219e-05 0.3964163 >11763 ILMN_8888888 -0.17331384 7.322021 -4.154668 7.110990e-05 0.3964163 >38906 ILMN_9999999 0.05532714 6.224477 4.093623 8.903425e-05 0.3964163 >4728 ILMN_0000000 0.05371882 6.177268 4.081921 9.293339e-05 0.3964163 > > *B* > >12919 3.236579 >6928 2.801832 >42428 1.818263 >36153 1.477781 >36152 1.364969 >28506 1.364927 >11763 1.352940 >38906 1.143329 >4728 1.103392 > > > > > >_______________________________________________ >Bioconductor mailing list >Bioconductor at stat.math.ethz.ch >https://stat.ethz.ch/mailman/listinfo/bioconductor >Search the archives: >http://news.gmane.org/gmane.science.biology.informatics.conductor Naomi S. Altman 814-865-3791 (voice) Associate Professor Dept. of Statistics 814-863-7114 (fax) Penn State University 814-865-1348 (Statistics) University Park, PA 16802-2111