首先安装xps和oligo。
安装oligo比较简单。直接
source("http://www.bioconductor.org/biocLite.R") biocLite("oligo") |
就可以了。
而安装xps可能麻烦一点，首先要安装ROOT，其主页为http://root.cern.ch/drupal/
而后是安装xps，为了确保安装效果，请使用：
source("http://www.bioconductor.org/biocLite.R") biocLite("xps",type="source") |
安装成功之后，就可以使用library(xps)调用了。
好了，我们先使用oligo来分析。
> library(oligo) Loading required package: oligoClasses Loading required package: Biobase Welcome to Bioconductor Vignettes contain introductory material. To view, type 'browseVignettes()'. To cite Bioconductor, see 'citation("Biobase")' and for packages 'citation("pkgname")'. Attaching package: 'Biobase' ..... Attaching package: 'oligo' > geneCELs<-list.celfiles(full.names=T) > geneCELs [1] "./PAO-1_(RaGene-1_0-st-v1).CEL" "./PAO-2_(RaGene-1_0-st-v1).CEL" "./PAO-3_(RaGene-1_0-st-v1).CEL" [4] "./PAR-1_(RaGene-1_0-st-v1).CEL" "./PAR-2_(RaGene-1_0-st-v1).CEL" "./PAR-3_(RaGene-1_0-st-v1).CEL" > affyGeneFS <- read.celfiles(geneCELs) Loading required package: pd.ragene.1.0.st.v1 Loading required package: RSQLite Loading required package: DBI Platform design info loaded. Reading in : ./PAO-1_(RaGene-1_0-st-v1).CEL Reading in : ./PAO-2_(RaGene-1_0-st-v1).CEL Reading in : ./PAO-3_(RaGene-1_0-st-v1).CEL Reading in : ./PAR-1_(RaGene-1_0-st-v1).CEL Reading in : ./PAR-2_(RaGene-1_0-st-v1).CEL Reading in : ./PAR-3_(RaGene-1_0-st-v1).CEL > affyGeneFS GeneFeatureSet (storageMode: lockedEnvironment) assayData: 1102500 features, 6 samples element names: exprs protocolData rowNames: PAO-1_(RaGene-1_0-st-v1).CEL PAO-2_(RaGene-1_0-st-v1).CEL ... PAR-3_(RaGene-1_0-st-v1).CEL (6 total) varLabels: exprs dates varMetadata: labelDescription channel phenoData rowNames: PAO-1_(RaGene-1_0-st-v1).CEL PAO-2_(RaGene-1_0-st-v1).CEL ... PAR-3_(RaGene-1_0-st-v1).CEL (6 total) varLabels: index varMetadata: labelDescription channel featureData: none experimentData: use 'experimentData(object)' Annotation: pd.ragene.1.0.st.v1 > # RMA - probeset level > genePS <- rma(affyGeneFS, target = "probeset") Background correcting Normalizing Calculating Expression > # RMA - transcript level > geneCore <- rma(affyGeneFS, target = "core") Background correcting Normalizing Calculating Expression > # For Exon arrays, there are three possible > # options for transcript level summarization: core, full and extended. For Gene arrays, only > # summaries for core probesets is available. > # Retrieving NetAffx Biological Annotation > genePS ExpressionSet (storageMode: lockedEnvironment) assayData: 213067 features, 6 samples element names: exprs protocolData rowNames: PAO-1_(RaGene-1_0-st-v1).CEL PAO-2_(RaGene-1_0-st-v1).CEL ... PAR-3_(RaGene-1_0-st-v1).CEL (6 total) varLabels: exprs dates varMetadata: labelDescription channel phenoData rowNames: PAO-1_(RaGene-1_0-st-v1).CEL PAO-2_(RaGene-1_0-st-v1).CEL ... PAR-3_(RaGene-1_0-st-v1).CEL (6 total) varLabels: index varMetadata: labelDescription channel featureData: none experimentData: use 'experimentData(object)' Annotation: pd.ragene.1.0.st.v1 > featureData(genePS) <- getNetAffx(genePS, "probeset") > genePS ExpressionSet (storageMode: lockedEnvironment) assayData: 213067 features, 6 samples element names: exprs protocolData rowNames: PAO-1_(RaGene-1_0-st-v1).CEL PAO-2_(RaGene-1_0-st-v1).CEL ... PAR-3_(RaGene-1_0-st-v1).CEL (6 total) varLabels: exprs dates varMetadata: labelDescription channel phenoData rowNames: PAO-1_(RaGene-1_0-st-v1).CEL PAO-2_(RaGene-1_0-st-v1).CEL ... PAR-3_(RaGene-1_0-st-v1).CEL (6 total) varLabels: index varMetadata: labelDescription channel featureData featureNames: 10700001 10700002 ... 10940690 (213067 total) fvarLabels: probesetid seqname ... probesettype (39 total) fvarMetadata: labelDescription experimentData: use 'experimentData(object)' Annotation: pd.ragene.1.0.st.v1 > featureData(geneCore) <- getNetAffx(geneCore, "transcript") > geneCore ExpressionSet (storageMode: lockedEnvironment) assayData: 29214 features, 6 samples element names: exprs protocolData rowNames: PAO-1_(RaGene-1_0-st-v1).CEL PAO-2_(RaGene-1_0-st-v1).CEL ... PAR-3_(RaGene-1_0-st-v1).CEL (6 total) varLabels: exprs dates varMetadata: labelDescription channel phenoData rowNames: PAO-1_(RaGene-1_0-st-v1).CEL PAO-2_(RaGene-1_0-st-v1).CEL ... PAR-3_(RaGene-1_0-st-v1).CEL (6 total) varLabels: index varMetadata: labelDescription channel featureData featureNames: 10700001 10700002 ... 10940690 (29214 total) fvarLabels: transcriptclusterid probesetid ... category (18 total) fvarMetadata: labelDescription experimentData: use 'experimentData(object)' Annotation: pd.ragene.1.0.st.v1 |
至此，我们手头就有了基因水平的表达结果，其后就可以使用limma来获得差异表达的基因了。
我们使用xps来分析。
> library(xps) Welcome to xps version 1.12.1 an R wrapper for XPS - eXpression Profiling System (c) Copyright 2001-2011 by Christian Stratowa Attaching package: 'xps' The following object(s) are masked from 'package:stats': residuals, weights > dir<-"/Users/jianhongou/Documents/Fumihiko/2011-08-23" > setwd(dir) > celfiles<-dir(pattern=".cel",ignore.case=T) > celfiles [1] "PAO-1_(RaGene-1_0-st-v1).CEL" "PAO-2_(RaGene-1_0-st-v1).CEL" "PAO-3_(RaGene-1_0-st-v1).CEL" "PAR-1_(RaGene-1_0-st-v1).CEL" [5] "PAR-2_(RaGene-1_0-st-v1).CEL" "PAR-3_(RaGene-1_0-st-v1).CEL" > libdir <- "/opt/local/apache2/htdocs/eam/RaGene-1_0-st-v1" > scheme.ragene10stv1r4 <- import.exon.scheme("Scheme_RaGene10stv1r4", + layoutfile=paste(libdir,"RaGene-1_0-st-v1.r4.clf",sep="/"), + schemefile=paste(libdir,"RaGene-1_0-st-v1.r4.pgf",sep="/"), + probeset=paste(libdir,"RaGene-1_0-st-v1.na31.rn4.probeset.csv",sep="/"), + transcript=paste(libdir,"RaGene-1_0-st-v1.na31.rn4.transcript.csv",sep="/")) Creating new file </Users/jianhongou/Scheme_RaGene10stv1r4.root>... Importing </opt/local/apache2/htdocs/eam/RaGene-1_0-st-v1/RaGene-1_0-st-v1.r4.clf> as <RaGene-1_0-st-v1.cxy>... <1102500> records imported...Finished New dataset <RaGene-1_0-st-v1> is added to Content... Importing </opt/local/apache2/htdocs/eam/RaGene-1_0-st-v1/RaGene-1_0-st-v1.na31.rn4.probeset.csv> as <RaGene-1_0-st-v1.anp>... Note: The following header columns are missing or in wrong order: <mouse_fl> <mouse_mrna> <rat_fl> <rat_mrna> Number of probesets is <213067>. Number of annotated transcripts is <29214>. <29214> records read...Finished14>. <28041> records imported...Finished > celnames<-gsub("\\(RaGene-1_0-st-v1\\).CEL","",celfiles) > celnames [1] "PAO-1_" "PAO-2_" "PAO-3_" "PAR-1_" "PAR-2_" "PAR-3_" > data <- import.data(scheme.ragene10stv1r4,filename="tmpdt_Data",filedir=dir,celdir=dir,celfiles=celfiles,celnames=celnames) Opening file </Users/jianhongou/Documents/Fumihiko/2011-08-23/Scheme_RaGene10stv1r4.root> in <READ> mode... Creating new temporary file </Users/jianhongou/Documents/Fumihiko/2011-08-23/tmpdt_Data_cel_20110824_173559.root>... Importing </Users/jianhongou/Documents/Fumihiko/2011-08-23/PAO-1_(RaGene-1_0-st-v1).CEL> as <PAO_1_.cel>... hybridization statistics: 1 cells with minimal intensity 23 3 cells with maximal intensity 65534 New dataset <DataSet> is added to Content... Importing </Users/jianhongou/Documents/Fumihiko/2011-08-23/PAO-2_(RaGene-1_0-st-v1).CEL> as <PAO_2_.cel>... hybridization statistics: 1 cells with minimal intensity 24 1 cells with maximal intensity 63204 Importing </Users/jianhongou/Documents/Fumihiko/2011-08-23/PAO-3_(RaGene-1_0-st-v1).CEL> as <PAO_3_.cel>... hybridization statistics: 1 cells with minimal intensity 22 1 cells with maximal intensity 36841 Importing </Users/jianhongou/Documents/Fumihiko/2011-08-23/PAR-1_(RaGene-1_0-st-v1).CEL> as <PAR_1_.cel>... hybridization statistics: 2 cells with minimal intensity 23 1 cells with maximal intensity 65533 Importing </Users/jianhongou/Documents/Fumihiko/2011-08-23/PAR-2_(RaGene-1_0-st-v1).CEL> as <PAR_2_.cel>... hybridization statistics: 2 cells with minimal intensity 23 2 cells with maximal intensity 65534 Importing </Users/jianhongou/Documents/Fumihiko/2011-08-23/PAR-3_(RaGene-1_0-st-v1).CEL> as <PAR_3_.cel>... Warning message: In import.data(scheme.ragene10stv1r4, filename = "tmpdt_Data", filedir = dir, : characters [](){}.:# etc in 'celnames' will be replaced with '_' hybridization statistics: 5 cells with minimal intensity 26 17 cells with maximal intensity 65534 > data <- attachInten(data) > temp <- intensity(data) > head(temp) X Y PAO_1__RaGene_1_0_st_v1_.cel_MEAN PAO_2__RaGene_1_0_st_v1_.cel_MEAN PAO_3__RaGene_1_0_st_v1_.cel_MEAN 1 0 0 8485 6825 8915 2 1 0 290 197 233 3 2 0 9184 7000 8774 4 3 0 279 191 184 5 4 0 52 58 66 6 5 0 174 127 137 PAR_1__RaGene_1_0_st_v1_.cel_MEAN PAR_2__RaGene_1_0_st_v1_.cel_MEAN PAR_3__RaGene_1_0_st_v1_.cel_MEAN 1 7683 9576 7932 2 201 259 212 3 7111 9466 7366 4 183 265 152 5 56 86 57 6 173 162 220 > data<-attachMask(data) > setwd(dir) > png("hist.png",width=1024,height=1024) > hist(data) > dev.off() quartz 2 > png("boxplot.png",width=1024,height=1024) > boxplot(data) > dev.off() quartz 2 > png("pmplot.png",width=1024,height=1024) > pmplot(data) > dev.off() quartz 2 > data<-removeInten(data) > data<-removeMask(data) > data.rma<-rma(data,"tmpdt_rma") |
其余略。。。。
https://github.com/benilton/oligo/wiki/Getting-the-grips-with-the-oligo-Package
博主，您好，看了您的文章，想请教您，“而安装xps可能麻烦一点，首先要安装ROOT，其主页为http://root.cern.ch/drupal/”我下载了”root_v5.34.14.source.tar”,这个解压了就相当于安装了,是吗？还有这个要安装在哪里？是和R一个文件夹里吗？还有“糗世界”在哪里注册啊？
对，直接安装就好了。不需要在一个文件夹中，但是需要让系统能找到你安装的目录，具体在windows下我无法测试。糗世界无法注册。
至此，我们手头就有了基因水平的表达结果，其后就可以使用limma来获得差异表达的基因了。
博主，有没有博文讲下一步limma的分析啊
你可以参考相关系列的博文。
博主您好，为什么执行代码genePS <- rma(affyGeneFS, target = "probeset")时候出现Error in .local(object, …) : unused argument (target = "probeset")
请检查你的sessionInfo(), oligo是否正常载入。
欧博您好，
我的出现错误了> affyGeneFS <- read.celfiles(geneCELs)
All the CEL files must be of the same type.
Error: checkChipTypes(filenames, verbose, "affymetrix", TRUE) ist nicht TRUE
怎么办啊
检查你的CEL文件是不是都是相同的affymetrix公司的平台的数据。你可以直接使用文本文档打开CEL文件，查看最初几行，看看它们的型号是否一致。如果是来自于不同型号的芯片，分开分析。
欧博，上个问题解决了，但是出现了新的错误>
> affyGeneFS
GeneFeatureSet (storageMode: lockedEnvironment)
assayData: 1102500 features, 39 samples
element names: exprs
protocolData
rowNames: GSM1085665_HE32H001.CEL GSM1085666_HE32H003.CEL …
GSM1085703_HE32H040.CEL (39 total)
varLabels: exprs dates
varMetadata: labelDescription channel
phenoData
rowNames: GSM1085665_HE32H001.CEL GSM1085666_HE32H003.CEL …
GSM1085703_HE32H040.CEL (39 total)
varLabels: index
varMetadata: labelDescription channel
featureData: none
experimentData: use ‘experimentData(object)’
Annotation: pd.hugene.1.0.st.v1
> # RMA – probeset level
> genePS eset<-rma(affyGeneFS)
新的错误》
> # RMA – probeset level
> genePS <- rma(affyGeneFS, target = "probeset")
Error in (function (classes, fdef, mtable) :
unable to find an inherited method for function ‘probeNames’ for signature ‘"GeneFeatureSet"’
你的pd.hugene.1.0.st.v1有没有加载上？把sessionInfo()传上来看一下。
欧博士，我用的是Mac系统（10.10.1),已经启用root用户，已下载root_v5.34.17.source.tar.gz，接下来应该如何安装root？
root的安装比较复杂。我建议您参照root的安装说明一步一步来。有些时候，可以试着使用port, fink等工具。
欧博士，
geneCELs<-list.celfiles(full.names=T)
data.raw <- read.celfiles(geneCELs)
data.rma sessionInfo()
R version 3.1.1 (2014-07-10)
Platform: x86_64-apple-darwin13.1.0 (64-bit)
locale:
[1] zh_CN.UTF-8/zh_CN.UTF-8/zh_CN.UTF-8/C/zh_CN.UTF-8/zh_CN.UTF-8
attached base packages:
[1] parallel stats graphics grDevices utils methods base
other attached packages:
[1] hugene10sttranscriptcluster.db_8.1.0 hugene10stprobeset.db_8.1.0
[3] hugene10stv1probe_2.14.0 org.Hs.eg.db_2.14.0
[5] AnnotationDbi_1.26.1 GenomeInfoDb_1.0.2
[7] affy_1.42.3 BiocInstaller_1.14.3
[9] hugene10stv1cdf_2.14.0 pd.hugene.1.0.st.v1_3.8.0
[11] RSQLite_1.0.0 DBI_0.3.1
[13] RColorBrewer_1.1-2 oligo_1.28.3
[15] Biostrings_2.32.1 XVector_0.4.0
[17] IRanges_1.22.10 oligoClasses_1.26.0
[19] GEOquery_2.30.1 Biobase_2.24.0
[21] BiocGenerics_0.10.0
loaded via a namespace (and not attached):
[1] affxparser_1.36.0 affyio_1.32.0 bit_1.1-12
[4] bitops_1.0-6 codetools_0.2-10 ff_2.2-13
[7] foreach_1.4.2 GenomicRanges_1.16.4 iterators_1.0.7
[10] preprocessCore_1.26.1 RCurl_1.95-4.5 splines_3.1.1
[13] stats4_3.1.1 tools_3.1.1 XML_3.98-1.1
[16] zlibbioc_1.10.0
欧博士，
geneCELs<-list.celfiles(full.names=T)
data.raw <- read.celfiles(geneCELs)
data.rma <- rma(data.raw, target = "core")
错误于(function (classes, fdef, mtable) :
unable to find an inherited method for function ‘probeNames’ for signature ‘"GeneFeatureSet"’
出现了和楼上一样的错误，请解惑？？
sessionInfo()
R version 3.1.1 (2014-07-10)
Platform: x86_64-apple-darwin13.1.0 (64-bit)
locale:
[1] zh_CN.UTF-8/zh_CN.UTF-8/zh_CN.UTF-8/C/zh_CN.UTF-8/zh_CN.UTF-8
attached base packages:
[1] parallel stats graphics grDevices utils methods base
other attached packages:
[1] hugene10sttranscriptcluster.db_8.1.0 hugene10stprobeset.db_8.1.0
[3] hugene10stv1probe_2.14.0 org.Hs.eg.db_2.14.0
[5] AnnotationDbi_1.26.1 GenomeInfoDb_1.0.2
[7] affy_1.42.3 BiocInstaller_1.14.3
[9] hugene10stv1cdf_2.14.0 pd.hugene.1.0.st.v1_3.8.0
[11] RSQLite_1.0.0 DBI_0.3.1
[13] RColorBrewer_1.1-2 oligo_1.28.3
[15] Biostrings_2.32.1 XVector_0.4.0
[17] IRanges_1.22.10 oligoClasses_1.26.0
[19] GEOquery_2.30.1 Biobase_2.24.0
[21] BiocGenerics_0.10.0
loaded via a namespace (and not attached):
[1] affxparser_1.36.0 affyio_1.32.0 bit_1.1-12
[4] bitops_1.0-6 codetools_0.2-10 ff_2.2-13
[7] foreach_1.4.2 GenomicRanges_1.16.4 iterators_1.0.7
[10] preprocessCore_1.26.1 RCurl_1.95-4.5 splines_3.1.1
[13] stats4_3.1.1 tools_3.1.1 XML_3.98-1.1
[16] zlibbioc_1.10.0
你load一下pd.ragene.1.0.st.v1.
如何load pd.ragene.1.0.st.v1.
library(pd.ragene.1.0.st.v1)。如果没有安装，就安装一下：library(BiocInstaller); biocLite(“pd.ragene.1.0.st.v1”); library(pd.ragene.1.0.st.v1)
我的错误和上面差不多，您说可以load，但是应该是我方式错了，然后没有成功，请问可以详细一点么？
> class(AEset)
[1] “ExpressionFeatureSet”
attr(,”package”)
[1] “oligoClasses”
> qc(AEset)
Error in (function (classes, fdef, mtable) :
unable to find an inherited method for function ‘qc’ for signature ‘”ExpressionFeatureSet”’
> sessionInfo()
R version 3.3.0 (2016-05-03)
Platform: i386-w64-mingw32/i386 (32-bit)
Running under: Windows 7 (build 7601) Service Pack 1
locale:
[1] LC_COLLATE=Chinese (Simplified)_People’s Republic of China.936
[2] LC_CTYPE=Chinese (Simplified)_People’s Republic of China.936
[3] LC_MONETARY=Chinese (Simplified)_People’s Republic of China.936
[4] LC_NUMERIC=C
[5] LC_TIME=Chinese (Simplified)_People’s Republic of China.936
attached base packages:
[1] stats4 parallel stats graphics
[5] grDevices utils datasets methods
[9] base
other attached packages:
[1] oligo_1.36.1
[2] Biostrings_2.40.2
[3] XVector_0.12.0
[4] IRanges_2.6.0
[5] S4Vectors_0.10.1
[6] arrayQualityMetrics_3.28.2
[7] BiocInstaller_1.22.2
[8] hgu95av2cdf_2.18.0
[9] simpleaffy_2.48.0
[10] gcrma_2.44.0
[11] genefilter_1.54.2
[12] CLL_1.12.0
[13] affy_1.50.0
[14] Biobase_2.32.0
[15] BiocGenerics_0.18.0
[16] oligoClasses_1.34.0
loaded via a namespace (and not attached):
[1] Rcpp_0.12.5
[2] Formula_1.2-1
[3] affxparser_1.44.0
[4] knitr_1.13
[5] magrittr_1.5
[6] cluster_2.0.4
[7] GenomicRanges_1.24.1
[8] splines_3.3.0
[9] zlibbioc_1.18.0
[10] bit_1.1-12
[11] xtable_1.8-2
[12] colorspace_1.2-6
[13] lattice_0.20-33
[14] foreach_1.4.3
[15] latticeExtra_0.6-28
[16] iterators_1.0.8
[17] preprocessCore_1.34.0
[18] setRNG_2013.9-1
[19] gridSVG_1.5-0
[20] Matrix_1.2-6
[21] acepack_1.3-3.3
[22] limma_3.28.6
[23] BeadDataPackR_1.24.2
[24] scales_0.4.0
[25] Hmisc_3.17-4
[26] annotate_1.50.0
[27] AnnotationDbi_1.34.3
[28] munsell_0.4.3
[29] affyPLM_1.48.0
[30] SVGAnnotation_0.93-1
[31] GenomeInfoDb_1.8.1
[32] hwriter_1.3.2
[33] plyr_1.8.4
[34] stringr_1.0.0
[35] tools_3.3.0
[36] base64_2.0
[37] grid_3.3.0
[38] nnet_7.3-12
[39] SummarizedExperiment_1.2.2
[40] data.table_1.9.6
[41] gtable_0.2.0
[42] ff_2.2-13
[43] beadarray_2.22.2
[44] DBI_0.4-1
[45] openssl_0.9.4
[46] yaml_2.1.13
[47] survival_2.39-4
[48] affyio_1.42.0
[49] RJSONIO_1.3-0
[50] gridExtra_2.2.1
[51] ggplot2_2.1.0
[52] reshape2_1.4.1
[53] RColorBrewer_1.1-2
[54] codetools_0.2-14
[55] rsconnect_0.4.3
[56] rpart_4.1-10
[57] RSQLite_1.0.0
[58] stringi_1.1.1
[59] XML_3.98-1.4
[60] vsn_3.40.0
[61] foreign_0.8-66
[62] chron_2.3-47
[63] illuminaio_0.14.0
[64] Cairo_1.5-9
现在有针对st array的qc了吗？
博主好，这个后边接的limma分析，上边评论里面贴的那个链接失效了。
所以能不能提供个可供参考的，好知道后边接着limma怎么搞。
参见limma的说明书，或者http://blog.qiubio.com:8080/archives/2083/4,这里只有非常粗浅的解释。
老师你好，按照您的指示，下载了ROOT文件（全是一堆文件，没有安装的执行文件），并解压了，做到biocLite(“xps”,”source=type”),一直发出警告：* installing *source* package ‘xps’ …
Warning: 运行命令’sh ./configure.win’的状态是127
ERROR: configuration failed for package ‘xps’
* removing ‘D:/R-3.4.2/library/xps’
这是为什么啊？？？哪里有问题啊
看情况，你还是使用oligo更合适一些。安装ROOT对于你可能有些过于困难了。