microarray analysis in R

Question

[bioinformatics]

Hi,

I'm trying to perform microarray analysis in R using the following commands:

BiocManager::install("affy")
BiocManager::install("oligo")
BiocManager::install("Biobase")
BiocManager::install("GEOquery")
 install.packages("splitstackshape")
BiocManager::instal

library("arrayQualityMetrics")
library(oligo)
library(Biobase)
library(affy)
library("split stackshape")
library("tidyr")
library("dplyr")
library("arrayQualityMetrics")

celFiles <- list.celfiles(CelFiles)
affyRaw <- read.celfiles(celFiles)

However, I get an error message after the last command

Error message: All the CEL files must be of the same type

Error in read.celfiles(celfiles) : 
  checkChipTypes(filenames, verbose, "affymetrix", TRUE) is not TRUE

Does anyone know how I might correct this?

Also the working directory has been set to the folder containing the data files.

Thank you

Answer 1

[Kevin Blighe]

It indicates that your CEL files are from multiple Affymetrix platforms. You will have to stratify these [CEL files] by platform and then process them separately.

Note that you should only need one of either oligo or affy. oligo should be okay for all Affymetrix cDNA experiments.

Kevin

[bioinformatics]

Sorry would you be able to explain what it meant to stratify these files by platform and process them separately? The two CEL files I want to compare were both generated using Affymetrix Microarray Suite v5.

https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSM766640 https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSM766537

I managed to get it to work, thanks for your help.

[Kevin Blighe]

Okay, they are all from the same study, which utilised the Affymetrix U133A. I am guessing that your variable, CelFiles, contained non-CEL files, and these triggered the error.

[bioinformatics]

I have formulated a table that has p values for comparison. Do you know how might I convert probes to gene ID?

[Kevin Blighe]

You can create an annotation lookup table for this array for your probe IDs, held in probes, via:

require(hgu133a.db)
annotLookup <- select(hgu133a.db, keys = probes,
  columns = c('PROBEID', 'ENSEMBL', 'SYMBOL'))

[bioinformatics]

Ok thanks for your insight! When I use the command:

annotLookup <- select(hgu133a.db, keys = probes,
  columns = c('PROBEID', 'ENSEMBL', 'SYMBOL'))

I get the following error message: Error in .testForValidKeys(x, keys, keytype, fks) : 'keys' must be a character vector

[Kevin Blighe]

Okay. Please re-check what I wrote. The variable, probes, must contain the probe IDs from your experiment - it should be a character vector.

[bioinformatics]

Ok thanks, how do I easily create a character vector with 16384 probe names?

[Kevin Blighe]

You should already have it as a character vector. Please display all of your code so far so that we can replicate what you are doing. We are not making any progress otherwise.

[bioinformatics]

I have the probe names listed as a column in the excel file with p values for differentially expressed genes. This is the code I used:

library("arrayQualityMetrics")
library(GEOquery)
library(oligo)
library(Biobase)
library(affy)
library("splitsackshape")
library("tidyr")
library("dplyr")


celFiles <- list.celfiles()
affyRaw <- oligo::rma(affyraw)
eset <- oligo::rma(affyRaw)
library(limma)
pData(eset)
Groups <- c("DDLPS", "DDLPS", "WDLPS", "WDLPS")
design <- model.matrix(~factor(Groups))
colnames(design) <- c("DDLPS", "DDLPSvsWDLPS")
fit <- lmFit(eset, design)
fit <- eBayes(fit)
option (digits =2)
res <- topTable (fit, number = Inf, adjust.method = "none", coef = 1)
write.table(res, "diff_exp.txt", sep= "\t)
require(hgu133a.db)
annotLookup <- select(hgu133a.db, keys = probes,
  columns = c('PROBEID', 'ENSEMBL', 'SYMBOL'))

Error in .testForValidKeys(x, keys, keytype, fks) : 'keys' must be a character vector

[Kevin Blighe]

A few points:

First point

why are you running these two lines - they are doing exactly the same thing:

affyRaw <- oligo::rma(affyraw)
eset <- oligo::rma(affyRaw)

It follows that the affyRaw variable's name is misleading because it will not contain raw data based on how you are creating you. It will contain RMA-normalised data.

Second point

The error about probes now makes sense - you never create this variable anywhere. It should be fine to use:

probes <- rownames(eset)

Third point

These lines are risky:

Groups <- c("DDLPS", "DDLPS", "WDLPS", "WDLPS")
design <- model.matrix(~factor(Groups))
colnames(design) <- c("DDLPS", "DDLPSvsWDLPS")

You should ensure that you set the order of your categorical variables. You should do something like:

Groups <- c("DDLPS", "DDLPS", "WDLPS", "WDLPS")
Groups <- factor(Group, levels = c('DDLPS','WDLPS'))
design <- model.matrix( ~ Groups)
colnames(design) <- c("DDLPS", "DDLPSvsWDLPS")

[bioinformatics]

Sorry I think I made a few typing errors in the post that I had tried to correct. Ok thankyou for your help, that command worked. I got a list of gene symbols. I will re-do the analysis to set the order of the categorical variables.

[bioinformatics]

Do you know if these commands are correct for a successful microarray differential analysis?

[Kevin Blighe]

My apologies - I am out of time and am now traveling internationally again.

[bioinformatics]

Ok no worries, thanks anyway.