Mercurial > repos > lnguyen > id_converter
comparison id_converter_UniProt.R @ 0:02c549457875 draft default tip
planemo upload
author | lnguyen |
---|---|
date | Fri, 15 Sep 2017 06:02:03 -0400 |
parents | |
children |
comparison
equal
deleted
inserted
replaced
-1:000000000000 | 0:02c549457875 |
---|---|
1 # Read file and return file content as data.frame | |
2 readfile = function(filename, header) { | |
3 if (header == "true") { | |
4 # Read only first line of the file as header: | |
5 headers <- read.table(filename, nrows = 1, header = FALSE, sep = "\t", stringsAsFactors = FALSE, fill = TRUE) | |
6 #Read the data of the files (skipping the first row) | |
7 file <- read.table(filename, skip = 1, header = FALSE, sep = "\t", stringsAsFactors = FALSE, fill = TRUE) | |
8 #And assign the header to the data | |
9 names(file) <- headers | |
10 } | |
11 else { | |
12 file <- read.table(filename, header = FALSE, sep = "\t", stringsAsFactors = FALSE, fill = TRUE) | |
13 } | |
14 return(file) | |
15 } | |
16 | |
17 # Mapping IDs using file built from Uniprot file source (ftp://ftp.uniprot.org/pub/databases/uniprot/current_release/knowledgebase/taxonomic_divisions/uniprot_sprot_human.dat.gz) | |
18 # Available databases: | |
19 # UNIPROT_AC: Uniprot accession number (e.g. P31946) | |
20 # UNIPROT_ID: Uniprot identifiers (e.g 1433B_HUMAN) | |
21 # GeneID_EntrezGene: Entrez gene ID (serie of digit) (e.g. 7529) | |
22 # RefSeq: RefSeq (NCBI) protein (e.g. NP_003395.1; NP_647539.1; XP_016883528.1) | |
23 # GI_number: GI (NCBI GI number) ID (serie of digits) assigned to each sequence record processed by NCBI (e.g; 21328448; 377656701; 67464627; 78101741) | |
24 # PDB: Protein DataBank Identifiers (e.g. 2BR9:A; 3UAL:A; 3UBW:A) | |
25 # GO_ID: GOterms (Gene Ontology) Identifiers (e.g. GO:0070062; GO:0005925; GO:0042470; GO:0016020; GO:0005739; GO:0005634) | |
26 # PIR: Protein Information Resource ID (e.g. S34755) | |
27 # OMIM: OMIM (Online Mendelian Inheritance in Man database) ID (serie of digits) (e.g: 601289) | |
28 # UniGene: Unigene Identifier (e.g. Hs.643544) | |
29 # Ensembl_ENSG: Ensembl gene identifiers (e.g. ENSG00000166913) | |
30 # Ensembl_ENST: Ensembl transcript identifiers (e.g. ENST00000353703; ENST00000372839) | |
31 # Ensembl_ENSP: Ensembl protein identifiers (e.g. ENSP00000300161; ENSP00000361930) | |
32 | |
33 mapping = function() { | |
34 # Extract arguments | |
35 args = commandArgs(trailingOnly = TRUE) | |
36 #print(args) | |
37 if (length(args) != 7) { | |
38 stop("Not enough/Too many arguments", call. = FALSE) | |
39 } | |
40 else { | |
41 input_id_type = args[1] | |
42 list_id = args[2] | |
43 list_id_input_type = args[3] | |
44 options = strsplit(args[4], ",")[[1]] | |
45 output = args[5] | |
46 uniprot_map_file = args[6] | |
47 np_uniprot_file = args[7] | |
48 | |
49 # Extract ID maps | |
50 uniprot_map = read.table(uniprot_map_file, header = TRUE, sep = "\t", stringsAsFactors = FALSE, fill = TRUE) | |
51 np_uniprot = read.table(np_uniprot_file, header = TRUE, sep = "\t", stringsAsFactors = FALSE, fill = TRUE) | |
52 | |
53 # Extract input IDs | |
54 if (list_id_input_type == "list") { | |
55 list_id = strsplit(args[2], " ")[[1]] | |
56 } | |
57 else if (list_id_input_type == "file") { | |
58 filename = as.character(strsplit(list_id, ",")[[1]][1]) | |
59 column_number = as.numeric(gsub("c", "" ,strsplit(list_id, ",")[[1]][2])) | |
60 header = strsplit(list_id, ",")[[1]][3] | |
61 file_all = readfile(filename, header) | |
62 list_id = c() | |
63 list_id = sapply(strsplit(file_all[,column_number], ";"), "[", 1) | |
64 } | |
65 names = c() | |
66 | |
67 # Map IDs | |
68 res = matrix(nrow=length(list_id), ncol=0) | |
69 | |
70 for (opt in options) { | |
71 names = c(names, opt) | |
72 # Map to neXtProt ID | |
73 if (opt == "neXtProt_ID") { | |
74 if (input_id_type == "UNIPROT_AC") { | |
75 mapped = sapply(strsplit(np_uniprot[match(list_id, np_uniprot$Uniprot_AC),]$neXtProt_ID, ";"), "[", 1) | |
76 } | |
77 else if (input_id_type == "neXtProt_ID") { | |
78 mapped = matrix(list_id) | |
79 } | |
80 else { | |
81 uniprot = sapply(strsplit(uniprot_map[match(list_id, uniprot_map[input_id_type][,]),]$UNIPROT_AC, ";"), "[", 1) | |
82 mapped = sapply(strsplit(np_uniprot[match(uniprot, np_uniprot$Uniprot_AC),]$neXtProt_ID, ";"), "[", 1) | |
83 } | |
84 } | |
85 # Map to other ID types | |
86 else { | |
87 if (input_id_type == "neXtProt_ID") { | |
88 uniprot = sapply(strsplit(np_uniprot[match(list_id, np_uniprot$neXtProt_ID),]$Uniprot_AC, ";"), "[", 1) | |
89 mapped = sapply(strsplit(uniprot_map[match(uniprot, uniprot_map$UNIPROT_AC),][opt][,], ";"), "[", 1) | |
90 } | |
91 else { | |
92 mapped = sapply(strsplit(uniprot_map[match(list_id, uniprot_map[input_id_type][,]),][opt][,], ";"), "[", 1) | |
93 } | |
94 } | |
95 res = cbind(res, matrix(mapped)) | |
96 } | |
97 | |
98 # Write output | |
99 if (list_id_input_type == "list") { | |
100 res = cbind(as.matrix(list_id), res) | |
101 names = c(input_id_type, names) | |
102 colnames(res) = names | |
103 write.table(res, output, row.names = FALSE, sep = "\t", quote = FALSE) | |
104 } | |
105 else if (list_id_input_type == "file") { | |
106 names(res) = options | |
107 names = c(names(file_all), names) | |
108 output_content = cbind(file_all, res) | |
109 colnames(output_content) = names | |
110 write.table(output_content, output, row.names = FALSE, sep = "\t", quote = FALSE) | |
111 } | |
112 } | |
113 } | |
114 | |
115 mapping() |