diff NcbigeneEntry.R @ 2:20d69a062da3 draft

planemo upload for repository https://github.com/workflow4metabolomics/lcmsmatching.git commit d4048accde6bdfd5b3e14f5394902d38991854f8
author prog
date Thu, 02 Mar 2017 08:55:00 -0500
parents
children
line wrap: on
line diff
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/NcbigeneEntry.R	Thu Mar 02 08:55:00 2017 -0500
@@ -0,0 +1,110 @@
+#####################
+# CLASS DECLARATION #
+#####################
+
+NcbigeneEntry <- methods::setRefClass("NcbigeneEntry", contains = "BiodbEntry")
+
+###########
+# FACTORY #
+###########
+
+createNcbigeneEntryFromXml <- function(contents, drop = TRUE) {
+
+	entries <- list()
+
+	# Define xpath expressions
+	xpath.expr <- character()
+	xpath.expr[[BIODB.ACCESSION]] <- "//Gene-track_geneid"
+	xpath.expr[[BIODB.KEGG.ID]] <- "/Dbtag_db[text()='KEGG']/..//Object-id_str"
+	xpath.expr[[BIODB.UNIPROT.ID]] <- "//Gene-commentary_heading[text()='UniProtKB']/..//Dbtag_db[text()='UniProtKB/Swiss-Prot']/..//Object-id_str"
+	xpath.expr[[BIODB.LOCATION]] <- "//Gene-ref_maploc"
+	xpath.expr[[BIODB.PROTEIN.DESCRIPTION]] <- "//Gene-ref_desc"
+	xpath.expr[[BIODB.SYMBOL]] <- "//Gene-ref_locus"
+	xpath.expr[[BIODB.SYNONYMS]] <- "//Gene-ref_syn_E"
+
+	for (content in contents) {
+
+		# Create instance
+		entry <- NcbigeneEntry$new()
+	
+		# Parse HTML
+		xml <-  XML::xmlInternalTreeParse(content, asText = TRUE)
+
+		# An error occured
+		if (length(XML::getNodeSet(xml, "//Error")) == 0 && length(XML::getNodeSet(xml, "//ERROR")) == 0) {
+
+			# Test generic xpath expressions
+			for (field in names(xpath.expr)) {
+				v <- XML::xpathSApply(xml, xpath.expr[[field]], XML::xmlValue)
+				if (length(v) > 0) {
+
+					# Eliminate duplicates
+					v <- v[ ! duplicated(v)]
+
+					# Set field
+					entry$setField(field, v)
+				}
+			}
+		
+			# CCDS ID
+			ccdsid <- .find.ccds.id(xml)
+			if ( ! is.na(ccdsid))
+				entry$setField(BIODB.NCBI.CCDS.ID, ccdsid)
+		}
+
+		entries <- c(entries, entry)
+	}
+
+	# Replace elements with no accession id by NULL
+	entries <- lapply(entries, function(x) if (is.na(x$getField(BIODB.ACCESSION))) NULL else x)
+
+	# If the input was a single element, then output a single object
+	if (drop && length(contents) == 1)
+		entries <- entries[[1]]
+
+	return(entries)
+
+	# Get data
+
+}
+
+################
+# FIND CCDS ID #
+################
+
+.find.ccds.id <- function(xml) {
+
+	# 1) Get all CCDS tags.
+	ccds_elements <- XML::getNodeSet(xml, "//Dbtag_db[text()='CCDS']/..//Object-id_str")
+
+	# 2) If all CCDS are the same, go to point 4.
+	ccds <- NA_character_
+	for (e in ccds_elements) {
+		current_ccds <- XML::xmlValue(e)
+		if (is.na(ccds))
+			ccds <- current_ccds
+		else {
+			if (current_ccds != ccds) {
+				ccds <- NA_character_
+				break
+			}
+		}
+	}
+
+	# 3) There are several CCDS values, we need to find the best one (i.e.: the most current one).
+	if (is.na(ccds)) {
+		# For each CCDS, look for the parent Gene-commentary tag. Then look for the text content of the Gene-commentary_label which is situed under. Ignore CCDS that have no Gene-commentary_label associated. Choose the CCDS that has the smallest Gene-commentary_label in alphabetical order.
+		version <- NA_character_
+		for (e in ccds_elements) {
+			versions <- XML::xpathSApply(e, "ancestor::Gene-commentary/Gene-commentary_label", XML::xmlValue)
+			if (length(versions) < 1) next
+			current_version <- versions[[length(versions)]]
+			if (is.na(version) || current_version < version) {
+				version <- current_version
+				ccds <- XML::xmlValue(e)
+			}
+		}
+	}
+
+	return(ccds)
+}