BiomaRt agilent array IDs
0
0
Entering edit mode
7.6 years ago
wiggs38 • 0

I have some gene expression data from a publicly available array data using the SurePrint G3 Human GE 8x60K Microarray platform. I am trying to annotate the Agilent probe ids with entrezIDs using biomaRt in R. However it appears that several Agilent IDs are missing in biomaRt.

As I am not familiar with the Agilent's technology I am not sure whether this issue may be due to the design of the array (i.e. custom probes etc) or whether there is missing data in the biomaRt package, or of course I have an error in my code.

The link provides the table of Agilent probe IDs along with other gene identifiers used by the expression data set

## This is just loading in the table from the link above.
probe = read.delim("GPL15931-probe_annotation.txt", comment.char = '#')

library(biomaRt)
ensembl = useMart("ensembl", dataset = "hsapiens_gene_ensembl")

## Should return all the agilent probe ids in biomaRt
agilent = getBM(attributes=c( 'efg_agilent_sureprint_g3_ge_8x60k' ),values="*", mart= ensembl)

First thing that strikes me is that for a 60K array there is only 31K probe IDs returned. Am I missing something here with either the technology or the code?

If I look for which probes are matched between the two datasets the difference is 11K probes. All the probes in biomaRt match but there is missing 11K that are in the GEO dataset.

table(probe$ID %in% agilent$efg_agilent_sureprint_g3_ge_8x60k)

Is there any agilent bioconductor packages that might have more complete IDs and gene identifiers? Any other thoughts on how to work around this problem?

biomaRt Agilent R • 3.9k views
ADD COMMENT
0
Entering edit mode

Although the question above still stands I did find a work around in case some is having a similar experience. The work around involved using the genomic coordinates of the probes and look for overlap of all genes (with entrez IDs).

## Have to use the same build as the array was done with. 
ensemblGRCh37 = useMart("ENSEMBL_MART_ENSEMBL", host = "grch37.ensembl.org", dataset = "hsapiens_gene_ensembl")
entrez = getBM(attributes = c("chromosome_name","start_position", "end_position", "strand","entrezgene"),values="*", mart= ensemblGRCh37)

## tidy up the entrez dataframe
chr = c(seq(1,22,1),"X","Y", "MT")
entrez = entrez[entrez$chromosome_name %in% chr, ]
entrez = entrez[!is.na(entrez$entrezgene),]
entrez$chromosome_name<-gsub("MT", "M", entrez$chromosome_name)

After splitting the genomic coordinates from the GEO data table into chr, start and stop I can use genomicRanges package. There was also an issue with reverse strand gene coordinates which required me to reverse them with a simple loop.

## loop to reverse the start and stop coords for reverse stranded genes.
for (i in 1:nrow(p1)){
  if (p1$stop[i] - p1$start[i] < 0){
    temp_start = p1$stop[i]
    temp_stop = p1$start[i]
    p1$start[i] = temp_start
    p1$stop[i] = temp_stop
  }
}

## Create GRanges objects and look for coord overlap
gr1 = with(p1, GRanges(chr, IRanges(start=start, end=stop)))
gr2 = with(entrez, GRanges(chromosome_name, IRanges(start=start_position, end=end_position)))

hits = findOverlaps(gr1, gr2)
probe<-cbind(p1[queryHits(hits),],entrez[subjectHits(hits),])
ADD REPLY

Login before adding your answer.

Traffic: 1671 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6