Dear Tim,
Currently xps does not directly support vsn, however it is possible to
combine xps with vsn, as the following code for the Affymetrix HuGene
dataset demonstrates:
1, create ROOT scheme file and import HuGene CEL-files:
### new R session: load library xps
> library(xps)
# directory containing Affymetrix library files
> libdir <- "/Volumes/GigaDrive/Affy/libraryfiles"
# directory containing Affymetrix annotation files
> anndir <- "/Volumes/GigaDrive/Affy/Annotation"
# directory to store ROOT scheme files
> scmdir <- "/Volumes/GigaDrive/CRAN/Workspaces/Schemes"
# create ROOT scheme file for whole genome array HuGene
> scheme.genome <-
import.exon.scheme("Scheme_HuGene10stv1r4_na29_hg18",filedir=scmdir,
paste(libdir,"HuGene-1_0-st-v1.r4.analysis-lib-files/HuGene-
1_0-st-v1.r4.clf",sep="/"),
paste(libdir,"HuGene-1_0-st-v1.r4.analysis-lib-files/HuGene-
1_0-st-v1.r4.pgf",sep="/"),
paste(anndir,"Version09Jul/HuGene-
1_0-st-v1.na29.hg18.probeset.csv",sep="/"),
paste(anndir,"Version09Jul/HuGene-
1_0-st-v1.na29.hg18.transcript.csv",sep="/"))
# directory containing Tissues CEL files
> celdir <- "/Volumes/GigaDrive/ChipData/Exon/HuGene"
# directory to store ROOT raw data files
> datdir <- "/Volumes/GigaDrive/CRAN/Workspaces/ROOTData"
# import tissues from Affymetrix Exon Array Dataset for HuGene-
1_0-st-v1
> celfiles <- c("TisMap_Breast_01_v1_WTGene1.CEL",
"TisMap_Breast_02_v1_WTGene1.CEL", "TisMap_Breast_03_v1_WTGene1.CEL",
"TisMap_Prostate_01_v1_WTGene1.CEL",
"TisMap_Prostate_02_v1_WTGene1.CEL",
"TisMap_Prostate_03_v1_WTGene1.CEL")
# rename CEL files
> celnames <- c("Breast01", "Breast02", "Breast03", "Prostate01",
"Prostate02", "Prostate03")
# import CEL files
> data.genome <- import.data(scheme.genome, "HuTissuesGenome6",
filedir=datdir, celdir=celdir, celfiles=celfiles, celnames=celnames)
Note: It is suggested to rename the CEL-files and/or save copies of
the
original CEL-files, since importing the results from "vsn" requires
conversion to CEL-files which may overwrite the original CEL-files!
2, convert raw data to expression levels using vsn and xps:
### new R session:
# need to load xps after vsn (to be able to use intensity<-)
> library(vsn)
> library(xps)
### first, load ROOT scheme file and ROOT data file
> scmdir <- "/Volumes/GigaDrive/CRAN/Workspaces/Schemes"
> scheme.genome <-
root.scheme(paste(scmdir,"Scheme_HuGene10stv1r4_na29_hg18.root",sep="/
"))
> datdir <- "/Volumes/GigaDrive/CRAN/Workspaces/ROOTData"
> data.genome <- root.data(scheme.genome,
paste(datdir,"HuTissuesGenome6_cel.root",sep="/"))
# attach intensity to data
> data.tmp <- attachInten(data.genome)
> str(data.tmp)
# get intensities as data.frame (includes X,Y)
> tmp <- intensity(data.tmp)
> head(tmp)
# apply vsn2 to intensities only (w/o X,Y)
> data.vsn2 <- vsn2(as.matrix(tmp[,3:ncol(tmp)]))
# convert result from log2
> value <- as.data.frame(2^data.vsn2 at hx)
# attach (X,Y) coordinates
> value <- cbind(tmp[,1:2], value)
> head(value)
# replace data with value and save in ROOT file:
# note: this step creates modified CEL-files!
> intensity(data.tmp, "tmp_DataVSN", TRUE) <- value
> str(data.tmp)
# apply rma
> datdir <- getwd()
> data.norm <- rma(data.tmp, "tmp_data_norm", filedir=datdir,
tmpdir="", background="none", normalize=FALSE, exonlevel="core+affx")
> str(data.norm)
> boxplot(data.norm)
# get data.frame
> expr.norm <- validData(data.norm)
> head(expr.norm)
I hope this example shows you how to use xps with vsn. In the long
term
I would like to implement vsn in xps, however this requires converting
the relevant code to C++ and it is not easy to find out which parts of
package vsn I need to convert.
Best regards
Christian
_._._._._._._._._._._._._._._._._._
C.h.r.i.s.t.i.a.n S.t.r.a.t.o.w.a
V.i.e.n.n.a A.u.s.t.r.i.a
e.m.a.i.l: cstrato at aon.at
_._._._._._._._._._._._._._._._._._
Tim Rayner wrote:
> Hi,
>
> I'm curious to know whether there are any plans to support the use
of
> the VSN algorithm (i.e. the vsn package) with the classes defined by
> the oligo package. The reason I ask is that it seems as though the
> latter package is now touted as one of the supported methods of
> handling Affymetrix GeneST expression arrays (e.g. HuGene, MoGene),
> but as far as I can see the only expression
> normalisation/summarisation method currently supported for objects
of
> the oligo GeneFeatureSet class seems to be RMA. I'm not even sure if
> GCRMA is supported for GeneFeatureSet objects? I guess the same
> questions also apply to the xps package, although I suspect that VSN
> support would turn out to be more work to implement there given the
> rather different underlying object structures used by xps. Do please
> correct me if I'm wrong about any of this, but until now I've been
> forced to shoehorn our data into the old AffyBatch workflow using
> custom CDFs and I'd like to know if there's a better (i.e.
supported)
> workflow I could be using.
>
> Best regards,
>
> Tim Rayner
>
>
>
> # Current setup for reference:
>
>> sessionInfo()
>>
> R version 2.10.0 Patched (2009-11-03 r50305)
> x86_64-apple-darwin9.8.0
>
> locale:
> [1] en_GB.UTF-8/en_GB.UTF-8/C/C/en_GB.UTF-8/en_GB.UTF-8
>
> attached base packages:
> [1] stats graphics grDevices utils datasets methods base
>
> other attached packages:
> [1] oligo_1.10.0 preprocessCore_1.8.0 oligoClasses_1.8.0
> [4] gcrma_2.18.0 vsn_3.14.0 affy_1.24.2
> [7] Biobase_2.6.1
>
> loaded via a namespace (and not attached):
> [1] affxparser_1.18.0 affyio_1.14.0 Biostrings_2.14.8 DBI_0.2-4
> [5] grid_2.10.0 IRanges_1.4.9 lattice_0.17-26
limma_3.2.1
> [9] splines_2.10.0
>
>
>
>