Hi, I have a question concerning the analysis of some affymetrix chips. I downloaded some of the data from GEO GSE11324 (see below). In doing so I'm stuck after I identified the probesets with significant changes. I have problems in assigning probeset specific gene names as well as getting the genomic coordinates. Furthermore I have no clue how to deal with the fact, that most genes have different probesets with differential transcriptional outcomes. I did this based on the affy and limma manuals like: targets file: Name FileName Target 0h1 GSM286031.CEL control 0h2 GSM286032.CEL control 0h3 GSM286033.CEL control 3h1 GSM286034.CEL three 3h2 GSM286035.CEL three 3h3 GSM286036.CEL three 6h1 GSM286037.CEL six 6h2 GSM286038.CEL six 6h3 GSM286039.CEL six library(affy) library(limma) library(vsn) pd <- read.AnnotatedDataFrame("er_for_affy.txt", header = TRUE, row.names = 2) pData(pd) #### load a <- ReadAffy(filenames = rownames(pData(pd)), phenoData = pd, verbose = TRUE) #### normalize x <- expresso(a, bg.correct = FALSE, normalize.method = "vsn", normalize.param = list(subsample = 1000), pmcorrect.method = "pmonly", summary.method = "medianpolish") ### genes with highest variation library(hgu133plus2.db) rsd <- apply(exprs(x), 1, sd) sel <- order(rsd, decreasing = TRUE)[1:250] u<-mget(row.names(exprs(x)[sel,]),hgu133plus2SYMBOL) heatmap(exprs(x)[sel, ], labRow=u) ### in this case it works to extract the gene symbol ### limma contrasts design <- model.matrix(~ -1+factor(c(1,1,1,2,2,2,3,3,3))) colnames(design) <- c("control", "three", "six") fit <- lmFit(x, design) meanSdPlot(x) contrast.matrix <- makeContrasts(three-control, six-control, levels=design) fit2 <- contrasts.fit(fit, contrast.matrix) fit2 <- eBayes(fit2) #### top list topTable(fit2, coef=1, adjust="BH", number=20, sort.by="M") library(hgu133plus2.db) u<-mget(row.names(fit2),hgu133plus2SYMBOL) How can I produce a topTable result with according gene names, somehow I do not understand the genelist argument? Obviously I still do not really understand how to manipulate R-data objects. I tried: d<-topTable(fit2, coef=1, adjust="BH", number=6600, sort.by ="P",genelist=fit$genes) symbol<-mget(d[,1], hgu133plus2SYMBOL) cbind(d, symbol) but I get Error in data.frame(..., check.names = FALSE) : arguments imply differing number of rows: 6600, 1???? anyhow I feel this might be way more complicated than actually necessary. Next, I would like to produce a standard clustering of the "fold changes" observed within (averaged) contrasts 1 (three - control) and 2 (six - control) and a heatmap presentation of the results. How to extract for example all fold-changes of those genes with a p-value<0.001 in at least one of the two contrasts? The coordinates of the genes I seem to get with v<-mget(row.names(fit2),hgu133plus2CHRLOC) v<-mget(row.names(fit2),hgu133plus2CHRLOCEND) But again I do not know, how to implement it into my fit2 object or topTable results. Furthermore there are many probesets with multiple mappings, should these not be excluded from the analysis? Actually, in the end I'd like to get the corresponding genes' coordinates in a way, that the maximum region size is given, eg in case of genes with multiple TSSs and TESs - does this mean I have to write a little function, that checks the strand of a gene and then returns min and max for gene start and end? As mentioned above, I do not know how to deal with the fact, that many genes are represented with mutliple probesets, often with different fold changes - is there a general recipe to deal with this question? Furthermore there are many probesets with multiple mappings, should these not be excluded from the analysis? I know it's a lot of questions, so is there a general source of information, that may help me to overcome the hurdles? Maxim [[alternative HTML version deleted]]