frma normalization function not working on ExpressionFeatureSet object

Question

[asalimih]

Hello,
I have an ExpressionFeatureSet object which is created by ae2bioc function from ArrayExpress package. Now i want to perform a Frozen Robust Multiarray Analysis using frma function from frma package but i get this error:

object must be of class AffyBatch, ExonFeatureSet, or GeneFeatureSet.

as i wanted to change some of srdf files of ArrayExpress manually I downloaded all raw and annotation data and used getAE function following by ae2bioc function (instead of main ArrayExpress function). and as a result it provided me an ExpressionFeatureSet object instead of AffyBatch. So my question is:

• can i convert ExpressionFeatureSet to any of those classes or is there any other way to solve my problem?

any help would be appreciated.

Here is my Code:

AEset = getAE("E-MEXP-950", path = "Data/E-MEXP-950" ,type='full',local = T,sourcedir = "Data/E-MEXP-950")
rawset= ae2bioc(mageFiles = AEset)
A = frma(rawset$`A-AFFY-33`)
B = frma(rawset$`A-AFFY-34`)

[ATpoint]

Please show some code that allows reproduction. [ Please read before posting a question ] -- How To Ask A Good Question

[asalimih]

I edited my post and added my code.

Answer 1

[Kevin Blighe]

I can reproduce the error by running the following code:

require(ArrayExpress)
AEset <- getAE("E-MEXP-950", type='full')
rawset <- ae2bioc(mageFiles = AEset)

require(frma)
A <- rawset$`A-AFFY-33`
frma(A)

Error in frma(A) : 
  object must be of class AffyBatch, ExonFeatureSet, or GeneFeatureSet.

There is a way around this issue, though.

You should have the CEL files on your local drive. Even with my command, above, the raw files are downloaded automatically and their paths are stored here:

as.character(A@protocolData@data@.Data[[1]])
 [1] "C:/Users/kevin/Documents/PN5538_U133A_1_1_1.CEL"        "C:/Users/kevin/Documents/PN72_U133A_1_1_1.CEL"         
 [3] "C:/Users/kevin/Documents/kiel1_U133A_1_1_1.CEL"         "C:/Users/kevin/Documents/PDN15Normal_U133A_1_1_1.CEL"  
 [5] "C:/Users/kevin/Documents/PDN16Normal_U133A_1_1_1.CEL"   "C:/Users/kevin/Documents/PDN17Normal_U133A_1_1_1.CEL"  
 [7] "C:/Users/kevin/Documents/PDN8_3_U133A_1_1_1.CEL"        "C:/Users/kevin/Documents/PKN13Normal_U133A_1_1_1.CEL"  
 [9] "C:/Users/kevin/Documents/PKN15Normal_U133A_1_1_1.CEL"   "C:/Users/kevin/Documents/PKN16Normal_U133A_1_1_1.CEL"  
[11] "C:/Users/kevin/Documents/PKN9Normal_U133A_1_1_1.CEL"    "C:/Users/kevin/Documents/pt10_U113A_1_1_1.CEL"         
[13] "C:/Users/kevin/Documents/33B_U133A_1_1_1.CEL"           "C:/Users/kevin/Documents/35B_U133A_1_1_1.CEL"          
[15] "C:/Users/kevin/Documents/39B_U133A_1_1_1.CEL"           "C:/Users/kevin/Documents/PT55_U133A_1_1_1.CEL"         
[17] "C:/Users/kevin/Documents/56A_U133A_1_1_1.CEL"           "C:/Users/kevin/Documents/PDT121-2Tumor_U133A_1_1_1.CEL"
[19] "C:/Users/kevin/Documents/PDT14Tumor_U133A_1_1_1.CEL"    "C:/Users/kevin/Documents/PD51_3_U133A_1_1_1.CEL"       
[21] "C:/Users/kevin/Documents/PDT81-2_U133A_1_1_1.CEL"       "C:/Users/kevin/Documents/PDT91-2_U133A_1_1_1.CEL"      
[23] "C:/Users/kevin/Documents/PKT4_1Tumor_U133A_1_1_1.CEL"   "C:/Users/kevin/Documents/PKT5Tumor_U133A_1_1_1.CEL"    
[25] "C:/Users/kevin/Documents/PKT9Tumor_U133A_1_1_1.CEL"

So, you just need to read those guys / gals into an AffyBatch object via the affy package:

files <- as.character(A@protocolData@data@.Data[[1]])
myAffyObject <- affy::ReadAffy(filenames = files)

Then, you can use frma():

myAffyObjectNormalised <- frma::frma(myAffyObject)
myAffyObjectNormalised
ExpressionSet (storageMode: lockedEnvironment)
assayData: 22283 features, 25 samples 
  element names: exprs, se.exprs 
protocolData: none
phenoData
  sampleNames: PN5538_U133A_1_1_1.CEL PN72_U133A_1_1_1.CEL ... PKT9Tumor_U133A_1_1_1.CEL (25 total)
  varLabels: sample
  varMetadata: labelDescription
featureData: none
experimentData: use 'experimentData(object)'
Annotation: hgu133a

It should automatically detect the annotation that is required. When you run frma(), you may have to install a new package, hgu133afrmavecs, and re-run the command (I did).

Kevin

[asalimih]

Thank you very much. yes it works but there is another problem. if you try your code on A-AFFY-34 you'll face an error saying hgu133bfrmavecs package is required. and actually we have no package named hgu133bfrmavecs .

[Kevin Blighe]

That is unfortunate. I may take a look later, however, you could also try to email the package maintainer: https://bioconductor.org/packages/release/bioc/html/frma.html

[Kevin Blighe]

I took another look but I am not sure how to circumvent that second issue. If you look at the hgu133afrmavecs package data, though, it seems to contain experimentally-derived values (I had planned to create a custom one for HGU133b, but now not sure):

require(hgu133afrmavecs)
data(hgu133afrmavecs)
str(hgu133afrmavecs)
List of 6
 $ normVec        : Named num [1:247965] 15.3 15.5 15.6 15.7 15.7 ...
  ..- attr(*, "names")= chr [1:247965] "129340" "213420" "396671" "82246" ...
 $ probeVec       : Named num [1:247965] -0.239 0.736 1.332 2.768 3.005 ...
  ..- attr(*, "names")= chr [1:247965] "129340" "213420" "396671" "82246" ...
 $ probeVarWithin : Named num [1:247965] 0.161 0.124 0.148 0.127 0.128 ...
  ..- attr(*, "names")= chr [1:247965] "129340" "213420" "396671" "82246" ...
 $ probeVarBetween: Named num [1:247965] 1.199 0.439 0.754 0.78 0.948 ...
  ..- attr(*, "names")= chr [1:247965] "129340" "213420" "396671" "82246" ...
 $ probesetSD     : Named num [1:22283] 1.008 0.619 0.568 0.455 0.334 ...
  ..- attr(*, "names")= chr [1:22283] "1007_s_at" "1053_at" "117_at" "121_at" ...
 $ medianSE       : Named num [1:22283] 0.196 0.171 0.184 0.235 0.219 ...
  ..- attr(*, "names")= chr [1:22283] "1007_s_at" "1053_at" "117_at" "121_at" ...

Unless you really need to use frma(), just use affy::rma(), which is AOK for both of these microarrays.

[asalimih]

Thank you for your time and help. Here is the package maintainer's reply:

Unfortunately, I never created the vectors for the hgu133b platform. You could create your own using the frmaTools package, but that's fairly involved.

for now I'll just use rma. but I'll try frmaTools in future too.

[Kevin Blighe]

Great - thanks for the feedback. I took a quick look and it seems possible to create your own vectors. Take a look at this reproducible example:

require(ArrayExpress)
AEset <- getAE("E-MEXP-950", type='full')
rawset <- ae2bioc(mageFiles = AEset)

B <- rawset$`A-AFFY-34`
files <- as.character(B@protocolData@data@.Data[[1]])
myAffyObject <- affy::ReadAffy(filenames = files)

require(frmaTools)
hgu133bfrmavecs <- makeVectorsAffyBatch(
  files,
  batch.id = rep('batch1', length(files)))

require(frma)
myAffyObjectNormalised <- frma(
  myAffyObject,
  input.vecs = hgu133bfrmavecs)

The only problem is that I run into an error eventually when running makeVectorsAffyBatch(), relating to the medianSE vector. Let me know if you also run into an error.

[asalimih]

yes i could reproduce the error :

Beginning Median SE Calculation ... 
Error in if (any(w < 0)) { : missing value where TRUE/FALSE needed

[Kevin Blighe]

If still a problem, perhaps post on Bioconductor Support site and ask for James, who is an expert on arrays