convertToRangedSummarizedExperiment: Convert experimental results retrieved by InterMineR queries...

Description Usage Arguments Details Value Author(s) References See Also Examples

View source: R/convertToRangedSummarizedExperiment.R

Description

convertToRangedSummarizedExperiment constitutes a wrapper function for converting genomic information and experimental data from InterMine to an object of the RangedSummarizedExperiment class.

Usage

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
convertToRangedSummarizedExperiment(
  im,
  dataset,
  SampleColumn,
  GeneColumn,
  ValueColumn,
  OrganismValue,
  colsForSampleMetadata,
  exonsForRowRanges = FALSE
)

Arguments

im

a Service object containing the base URL and API token, created by initInterMine.

dataset

a data.frame retrieved with InterMineR queries that contains experimental data from high-throughput assays.

SampleColumn

a character string or an integer indicating which column of the dataset contains the samples.

GeneColumn

a character string or an integer indicating which column of the dataset contains the genes.

ValueColumn

a character string or an integer indicating which column of the dataset contains the experimental data.

OrganismValue

a character string with the name of the organism from which the genomic information are retrieved.

colsForSampleMetadata

an integer vector indicating the columns of the dataset which will be assigned as sample metadata in the colData argument of the RangedSummarizedExperiment.

exonsForRowRanges

a logical value indicating whether the rowRanges argument of the RangedSummarizedExperiment should be assigned with a GRangeslist containing only annotations about the genes (default) or annotations about all exons of each gene.

Details

The InterMineR package provides a flexible interface to InterMine web services, which allow for rapid retrieval of various genomic information.

convertToRangedSummarizedExperiment function facilitates the conversion of genomic information and experimental data from high-throughput assays, which are retrieved by using InterMineR queries, to an object of the RangedSummarizedExperiment class.

It is noteworthy that the reshape function is used to convert the experimental data from the InterMineR format (long format) to the matrix (wide format) assigned to the assays argument of the SummarizedExperiment.

Value

an object of the RangedSummarizedExperiment class containing genomic information and experimental data of high-throughput assays, which are retrieved with the InterMineR queries system.

Author(s)

InterMine Team

References

SummarizedExperiment for Coordinating Experimental Assays, Samples, and Regions of Interest

See Also

RangedSummarizedExperiment, link{convertToGRanges}, link{initInterMine}

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
# 10 Drosophila melanogaster genes of interest
Drosophila.genes = c("BEAF-32", "Antp", "bcd", "caup", "tup", "E2f2", "dsx", "so", "toy", "Lim1")

# retrieve microarray time course experimental data for Drosophila.genes
# get FlyMine instance
im.fly = initInterMine(listMines()["FlyMine"])

# get FlyMine microarray time course template query
queryForData = getTemplateQuery(im.fly, "Gene_TimeCourseExpression")

test.data = list(NULL)
ind.null = c()

for(i in seq(length(Drosophila.genes))){
  # set value in gene constraint
  queryForData$where[[3]]$value = as.character(Drosophila.genes[i])
  
  # run query and save the results of genes that exist in the Microarray time course dataset
  r = runQuery(im.fly, queryForData)
  
  ind.null = c(ind.null, is.null(r))
  
  test.data[[i]] = r
  
}

# remove genes for which no experimental data were retrieved
test.data = test.data[which(!ind.null)]

# rbind data together
test.data = do.call(rbind, test.data)

# using integer index for columns in arguments
test1 = convertToRangedSummarizedExperiment(
  im = im.fly,
  dataset = test.data,
  SampleColumn = 2,
  GeneColumn = 1,
  ValueColumn = 3,
  OrganismValue = "Drosophila melanogaster",
  colsForSampleMetadata = 4:7,
  exonsForRowRanges = TRUE
)

test1

# using directly column names in arguments
test2 = convertToRangedSummarizedExperiment(
  im = im.fly,
  dataset = test.data,
  SampleColumn = "Gene.microArrayResults.assays.sample2",
  GeneColumn = "Gene.symbol",
  ValueColumn = "Gene.microArrayResults.value",
  OrganismValue = "Drosophila melanogaster",
  colsForSampleMetadata = 4:7,
  exonsForRowRanges = TRUE
)

test2

intermine/rintermine documentation built on Jan. 10, 2022, 4:24 p.m.