Mercurial > repos > testtool > find_dmr

comparison findDMR/.Rhistory @ 12:083895dbc289 draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Uploaded
author testtool
date Mon, 12 Jun 2017 11:05:04 -0400
parents 7df2b7d79391
children
comparison
equal deleted inserted replaced
11:24ac6f93cc3e 12:083895dbc289
1 clusterSize=2 1 ??fread
2 class(clusterSize)
3 type(clusterSize)
4 require("minfi", quietly = TRUE)
5 ??getGenomicRatioSetFromGEO
6 getGenomicRatioSetFromGEO
7 getwd()
8 GEODataTable <- getGEO(GSE51547)
9 require("BiocGenerics", quietly = TRUE)
10 require("data.table", quietly = TRUE)
11 require("GEOquery", quietly = TRUE)
12 require("rtracklayer", quietly = TRUE)
13 require("FDb.InfiniumMethylation.hg19", quietly = TRUE)
14 GEODataTable <- getGEO(GSE51547)
15 GEODataTable <- getGEO("GSE51547")
16 IlmnIDTable <- Table(GEODataTable)
17 MetaData <- data.frame(Meta(GEODataTable))
18 IlmnIDTable <- Table(GEODataTable)
19 GEODataTable
20 GSE<-"GSE51547"
21 function (GSE = NULL, path = NULL, array = "IlluminaHumanMethylation450k",
22 annotation = .default.450k.annotation, what = c("Beta", "M"),
23 mergeManifest = FALSE, i = 1)
24 {
25 what <- match.arg(what)
26 if (is.null(GSE) && is.null(path))
27 stop("Either GSE or path must be supplied.")
28 if (!is.null(GSE))
29 gset <- GEOquery::getGEO(GSE)
30 else gset <- GEOquery::getGEO(filename = file.path(path,
31 list.files(path, pattern = ".soft")))
32 if (length(gset) == 0)
33 stop("Empty list retrieved from GEO.")
34 if (length(gset) > 1) {
35 warning("More than one ExpressionSet found:\n", names(gset),
36 "\nUsing entry ", i)
37 gset <- gset[[i]]
38 }
39 else gset <- gset[[1]]
40 platform <- annotation(gset)
41 if (platform != "GPL13534")
42 warning(sprintf("%s is not the platform ID associated with IlluminaHumanMethylation450k. Should be GPL13534.",
43 platform))
44 if (what == "Beta" & (min(exprs(gset)[, 1], na.rm = TRUE) <
45 0 | max(exprs(gset)[, 1], na.rm = TRUE) > 1))
46 warning("Values outside [0,1] detected. 'what' argument should not be Beta.")
47 ann <- .getAnnotationString(c(array = array, annotation = annotation))
48 if (!require(ann, character.only = TRUE))
49 stop(sprintf("cannot load annotation package %s", ann))
50 object <- get(ann)
51 gr <- getLocations(object, mergeManifest = mergeManifest,
52 orderByLocation = TRUE)
53 locusNames <- names(gr)
54 sampleNames(gset) <- gset$title
55 common <- intersect(locusNames, fData(gset)$Name)
56 if (length(common) == 0)
57 stop("No rowname matches. 'rownames' need to match IlluminaHumanMethylation450k probe names.")
58 ind1 <- match(common, fData(gset)$Name)
59 ind2 <- match(common, locusNames)
60 preprocessing <- c(rg.norm = paste0("See GEO ", GSE, " for details"))
61 if (what == "Beta") {
62 out <- GenomicRatioSet(gr = gr[ind2, ], Beta = exprs(gset)[ind1,
63 , drop = FALSE], M = NULL, CN = NULL, pData = pData(gset),
64 annotation = c(array = array, annotation = annotation),
65 preprocessMethod = preprocessing)
66 }
67 else {
68 out <- GenomicRatioSet(gr = gr[ind2, ], Beta = NULL,
69 M = exprs(gset)[ind1, , drop = FALSE], CN = NULL,
70 pData = pData(gset), annotation = c(array = array,
71 annotation = annotation), preprocessMethod = preprocessing)
72 }
73 return(out)
74 }
75 <environment: namespace:minfi>
76 function (GSE = NULL, path = NULL, array = "IlluminaHumanMethylation450k",
77 annotation = .default.450k.annotation, what = c("Beta", "M"),
78 mergeManifest = FALSE, i = 1)
79 {
80 what <- match.arg(what)
81 if (is.null(GSE) && is.null(path))
82 stop("Either GSE or path must be supplied.")
83 if (!is.null(GSE))
84 gset <- GEOquery::getGEO(GSE)
85 else gset <- GEOquery::getGEO(filename = file.path(path,
86 list.files(path, pattern = ".soft")))
87 if (length(gset) == 0)
88 stop("Empty list retrieved from GEO.")
89 if (length(gset) > 1) {
90 warning("More than one ExpressionSet found:\n", names(gset),
91 "\nUsing entry ", i)
92 gset <- gset[[i]]
93 }
94 else gset <- gset[[1]]
95 platform <- annotation(gset)
96 if (platform != "GPL13534")
97 warning(sprintf("%s is not the platform ID associated with IlluminaHumanMethylation450k. Should be GPL13534.",
98 platform))
99 if (what == "Beta" & (min(exprs(gset)[, 1], na.rm = TRUE) <
100 0 | max(exprs(gset)[, 1], na.rm = TRUE) > 1))
101 warning("Values outside [0,1] detected. 'what' argument should not be Beta.")
102 ann <- .getAnnotationString(c(array = array, annotation = annotation))
103 if (!require(ann, character.only = TRUE))
104 stop(sprintf("cannot load annotation package %s", ann))
105 object <- get(ann)
106 gr <- getLocations(object, mergeManifest = mergeManifest,
107 orderByLocation = TRUE)
108 locusNames <- names(gr)
109 sampleNames(gset) <- gset$title
110 common <- intersect(locusNames, fData(gset)$Name)
111 if (length(common) == 0)
112 stop("No rowname matches. 'rownames' need to match IlluminaHumanMethylation450k probe names.")
113 ind1 <- match(common, fData(gset)$Name)
114 ind2 <- match(common, locusNames)
115 preprocessing <- c(rg.norm = paste0("See GEO ", GSE, " for details"))
116 if (what == "Beta") {
117 out <- GenomicRatioSet(gr = gr[ind2, ], Beta = exprs(gset)[ind1,
118 , drop = FALSE], M = NULL, CN = NULL, pData = pData(gset),
119 annotation = c(array = array, annotation = annotation),
120 preprocessMethod = preprocessing)
121 }
122 else {
123 out <- GenomicRatioSet(gr = gr[ind2, ], Beta = NULL,
124 M = exprs(gset)[ind1, , drop = FALSE], CN = NULL,
125 pData = pData(gset), annotation = c(array = array,
126 annotation = annotation), preprocessMethod = preprocessing)
127 }
128 return(out)
129 }
130 GSE<-"GSE51547"
131 function (GSE = NULL, path = NULL, array = "IlluminaHumanMethylation450k",
132 annotation = .default.450k.annotation, what = c("Beta", "M"),
133 mergeManifest = FALSE, i = 1)
134 {
135 what <- match.arg(what)
136 if (is.null(GSE) && is.null(path))
137 stop("Either GSE or path must be supplied.")
138 if (!is.null(GSE))
139 gset <- GEOquery::getGEO(GSE)
140 else gset <- GEOquery::getGEO(filename = file.path(path,
141 list.files(path, pattern = ".soft")))
142 if (length(gset) == 0)
143 stop("Empty list retrieved from GEO.")
144 if (length(gset) > 1) {
145 warning("More than one ExpressionSet found:\n", names(gset),
146 "\nUsing entry ", i)
147 gset <- gset[[i]]
148 }
149 else gset <- gset[[1]]
150 platform <- annotation(gset)
151 if (platform != "GPL13534")
152 warning(sprintf("%s is not the platform ID associated with IlluminaHumanMethylation450k. Should be GPL13534.",
153 platform))
154 if (what == "Beta" & (min(exprs(gset)[, 1], na.rm = TRUE) <
155 0 | max(exprs(gset)[, 1], na.rm = TRUE) > 1))
156 warning("Values outside [0,1] detected. 'what' argument should not be Beta.")
157 ann <- .getAnnotationString(c(array = array, annotation = annotation))
158 if (!require(ann, character.only = TRUE))
159 stop(sprintf("cannot load annotation package %s", ann))
160 object <- get(ann)
161 gr <- getLocations(object, mergeManifest = mergeManifest,
162 orderByLocation = TRUE)
163 locusNames <- names(gr)
164 sampleNames(gset) <- gset$title
165 common <- intersect(locusNames, fData(gset)$Name)
166 if (length(common) == 0)
167 stop("No rowname matches. 'rownames' need to match IlluminaHumanMethylation450k probe names.")
168 ind1 <- match(common, fData(gset)$Name)
169 ind2 <- match(common, locusNames)
170 preprocessing <- c(rg.norm = paste0("See GEO ", GSE, " for details"))
171 if (what == "Beta") {
172 out <- GenomicRatioSet(gr = gr[ind2, ], Beta = exprs(gset)[ind1,
173 , drop = FALSE], M = NULL, CN = NULL, pData = pData(gset),
174 annotation = c(array = array, annotation = annotation),
175 preprocessMethod = preprocessing)
176 }
177 else {
178 out <- GenomicRatioSet(gr = gr[ind2, ], Beta = NULL,
179 M = exprs(gset)[ind1, , drop = FALSE], CN = NULL,
180 pData = pData(gset), annotation = c(array = array,
181 annotation = annotation), preprocessMethod = preprocessing)
182 }
183 return(out)
184 }
185 out
186 GSE<-"GSE51547"
187 gset <- getGEO(GSE)
188 ann <- .getAnnotationString(c(array = array, annotation = annotation))
189 platform <- annotation(gset)
190 hm450.hg19 <- getPlatform()
191 IlmnInfo <-
192 data.table(
193 IlmnID = names(hm450.hg19),
194 CHR = as.data.frame(hm450.hg19@seqnames)$value,
195 BP = as.numeric(hm450.hg19@elementMetadata$probeStart)
196 )
197 IlmnIDTable <- Table(gset)
198 ??getGenomicRatioSetFromGEO()
199 getGenomicRatioSetFromGEO()
200 getGenomicRatioSetFromGEO