#pathway_functions.R
####################
#'@title Retreive the default set of cellular pathways.
#'@description Retreives the default set of cellular pathways. These are the gene sets for all Reactome pathways as downloaded on 4/6/14 .
#'@param path_file string giving the relative or absolute file path to an accessible set of cellular pathways. Set of pathways must first have been recorded in the reference_data/paths/pathMetaData.txt file.
#'@return A \code{Path_Detail} object.
#'@export
#'@examples
#'cellularPathwayGeneSets = getDefaultPaths()
<-( path_file=("extdata/Reactome.2014.04.06.12.52.27.txt", package = "packageDir") ){
path_detail_tmp<-(path_file=path_file, verbose=F)
(path_detail_tmp)
}
checkPathsMatch<-(STUDY){
fileInMetaData = STUDY@studyMetaData@paths$
allInputArms = (STUDY@results)[grep(="functional_|aberration_", x=(STUDY@results))]
(carm allInputArms){
( !"summarystats" %in% (STUDY@results[carm]]) ){
( "coverage_summary" %in% (STUDY@results[carm]]) ){
ss = STUDY@results[carm]]$coverage_summary$summarystats
}
}{
ss = STUDY@results[carm]]$summarystats
}
sspathfile = ss[ss,1]=="Path source:",2]
( sspathfile!=fileInMetaData ){
("Warning: pathways currently loaded do not match pathways used to analyze input data types!")
("Current input arm: ", carm)
("File for the pathways that are currently loaded:\n", fileInMetaData)
("File for the pathways that were used in the data analysis:\n", sspathfile)
uin = ("Press any key to continue")
}
}
}
#pathsFromGraphite function
#takes path name
#returns nodes as HUGO symbols
getNodes<-(path){
#returns nodes in a pathway with their identifiers converted to hugo symbols
(nodes(convertIdentifiers(pathway=path, ="symbol")))
}
# importAllGraphite()
#'@title Import all pathway repositories from the graphite package to this package.
#'@description Gives options to import pathway repositories from Reactome, Spike, NCI, KEGG and Biocarta.
#'@param repositories A \code{list} object, each slot containing a pathway repository from graphite that is to be imported
#'@param symtab A HUGO symbol look up table, as provided by the package HGNChelper or from genenames.org.
#'@param verbose A flag indicating if this function should be run in an interactive mode. This allows individual pathway repositories to be selected, insted of automatically selecting all pathways.
#'@import graphite
#'@import tools
<-(repositories = (reactome=reactome, spike=spike, nci=nci, kegg=kegg, biocarta=biocarta),
symtab=,
verbose=F){
((symtab)) symtab = ()
(verbose){
((=(repositories), =1, =(1:(repositories), "Repository Name")))
(T){
= ("Please enter the number(s) of the pathway set(s) you would like to import from graphite\n (separated by spaces): ")
= (((x=, =" ")[1]]), silent=T)
(!(())&!(!%in%1:(repositories)))
("\nSorry, that input was not understood, please try again.\n")
}
repositories = repositories[line]
}
(pws repositories){
preped_paths = importFromGraphite(db=pws)
#correct gene symbols
(preped_paths$paths)<-(symbol_set=(preped_paths$paths),
symref=symtab,
verbose=T)
#setPathMetaData
preped_paths = manualPathMetaData(preped_paths=preped_paths, symtab=symtab)
#establish path set
recordPathSet(ps=preped_paths)
}
}
importFromGraphite<-(db){
outset = ()
= ()
dbname = ()
(i 1:(db)){
pname = (db)[i]
(pname, (nodes(db[pname]])), "\n")
outset[pname]] = getNodes(path=db[pname]])
= (, db[pname]]@)
dbname = (dbname, db[pname]]@database)
}
bip = (pth=outset)
# neededSlots = c("name", "date", "source")
preped_paths = ()
preped_paths"date"]] = (db[1]]@, sep=" ", collapse=" ")
preped_paths"source"]] = (db[1]]@database, "via the bioconductor graphite package", sep=" ", collapse=" ")
preped_paths"name"]] = db[1]]@database
preped_paths"paths"]] = (bip)
preped_paths"original_file_or_source"]] = (db[1]]@database, "from Bioconductor package 'graphite', version", ("graphite") )
preped_paths"original_file_creation_date"]] = ("graphite utilized on",(()))
(preped_paths)
}
#parseReactomeUniprot()
# loads paths from file in this format: <path name><\t><gene><space><gene>....
# everything after and including "firstGeneColum" is assumed to be a gene"
#the first column is assumed to be the path name or other unique identifier for the pathway
#
#takes: firstGeneColum: the first tab stop with gene names in it, usually 2 or 3
# keepLowerCaseGenes: keep genes who are all lower case
parseReactomeUniprot<-(firstVectorColum, fname = ){
pset = (fname, sep="\t",comment.char="", stringsAsFactors=F, ="")
outlist = ()
(i 1:(pset)){
cur = pset[i,]
#pull out the "[19 processes]:" part
paths = (="[[0-9]+ processes]: ",replacement="",x=cur[firstVectorColum])
gene = cur[1,1]
pathsV = (paths,="; ")[1]]
outlist[gene]] = pathsV
}
(outlist)
}
#switchIds
#takes set of ids (vector)
#returns set of ids with as many switched as possible
switchIds<-(idv,symtab=){
#see which can be switched
#idv%in%symtab$
dict = ("testSymbol")
(dict) = ("UniProt:P12956")
(i 1:(idv)){
swap = dict[idv[i]]
(!(swap)){
idv[i] = swap
}
}
(idv)
}
switchPathSymbols<-(po){
#first switch the node names
po@nodes = switchIds(idv=po@nodes)
#second switch the names in the edges
po@edges$src = switchIds(idv=po@edges$src)
po@edges$dest = switchIds(idv=po@edges$dest)
(po)
}
#'@title Convert list of vectors to bipartate graph format.
#'@description Used to restructure list of vectors to bipartate graph format.
#'@param pth A list object, with each slot containing a vector.
#'@param uni The full, unique set of elements in all the vectors in the list passed in the argument, pth. If not provided this will be automatically recalculated.
#'@return \code{matrix} object with logical values representing a bipartate graph; rows are set as the slot names of the input list. Columns are set as all the unique elements in the vectors the list contains.
#'@export
<-(pth, uni=){
#switch from list pathway format to logic table pathway format
#takes: pathways in list format
#returns: pathways in logic table format
((uni)){uni = getUniverse(pth)}
#initilize the bipartate graph table with all Fs (non membership)
pmat = (=(pth), = (uni), = (x=F, times=(pth)*(uni)))
(pmat) <- (pth)
(pmat) <- uni
(i 1:(pth)){
= pth[i]]
#remove any missing data (ie, the NAs)
= [which(!())]
#set appropriate membership of columns (gene names) in current path
pmat[i,] = T
}
(pmat)
}
#like matrix to list but does not produce empty lists
#takes a bipartite graph/logical matrix with row and column names the pathways and elements, respectively
#outputs a list, with each member a vector path members
matrix_to_list<-(m){
out = ()
(i 1:(m)){
= (m)[i]
members = ((m)[m[i,]])
((members)){
out[name]] = members
}
}
(out)
}
getUniverse<-(gsv)
{
####getUniverse
#takes pathway list (list$pathname = [node1,node2,node3])
#returns "universe", the list of unique nodes in all pathways
num = ((gsv))
uni = ()
tot = 0
(i 1:num){
tot = tot + (gsv[i]])
uni = (uni, gsv[i]])
}
#print(paste("sum of numbers of nodes in all pathways:", tot))
#print(paste("total number of unique nodes:", length(uni)))
(uni)
}
test.getPaths<-(){
symtab = ()
path_file=
referenceFileName = ("extdata/pathMetaData.txt", package = "packageDir")
testPths = (symtab=symtab, referenceFileName=referenceFileName)
((testPths$paths))
(testPths$paths[1:2,1:2])
}
#getPaths
#main function to obtain pathways
#requires path meta data be available in file "./reference_data/paths/pathMetaData.txt"
# takes: path_file: the file name of the pathways to be loaded; if provided, interactive pathway selection will not be used
# referenceFileName: the name of the path meta data file
# symtab: the symbol look up table for hugo symbols
# returns path_detail object
#'@title Get set of cellular pathways
#'@description Get \code{Path_Detail} object, housing gene sets for cellular pathways.
#'@param path_file String giving path to file from which pathways should be loaded. Alternatively, a \code{Study} object can be provided here and pathways will automatically be retreived from it.
#'@param referenceFileName Optional. String giving path to file where pathway meta data (date loaded, path repository name, etc..) can be found.
#'@param symtab A gene identifier lookup table. If not provided, one will be automatically obtained (may take longer)
#'@param verbose Logical flag indicating if interactive user prompts and additional information should be displayed.
#'@return \code{Path_Detail} object.
#'@details
#'Pathways are stored in the local working directory in the folder \code{./reference_data/paths/}.
#'Pathways must be imported into this package so that meta-data and import records can be maintained.
#'Path meta data is stored in the file \code{./reference_data/paths/pathMetaData.txt} . Any time a path
#'repository is imported, a record of its import is added to the path meta data file. When loading a path
#'repository for use in a study, this program checks the path meta data file, for available path repositories
#'and for path repository file locations.
#'@export
#'@examples
#'
#'\dontrun{
#'#select and/or load pathways interactively
#'path_detail = getPaths()
#'}
#'
#'#Get pathways from a study object:
#'study=getTestStudyObject()
#'path_detail=getPaths(study)
#'
#'#select pathways from their file name:
#'path_detail = getPaths(path_file="./reference_data/paths/Reactome.2014.04.06.12.52.27.txt",
#' verbose=FALSE)
<-(path_file=,
referenceFileName = "./reference_data/paths/pathMetaData.txt",
symtab=,
verbose=T){
((path_file)=="Study") (path_file@studyMetaData@paths)
#initilize
((symtab)) symtab = (verbose=verbose)
(!(path_file)) path_file=
(!(path_file)){
(!(path_file)){
("\n\nThe path file\n",path_file,"\ncould not be found.\n\n")
# blnk = readline("Press enter to continue and select a different pathway gene set file")
path_file =
}
}
#provide user selections and path file selection
(T){#allows repeated loading of pathways
pRecord = choosePaths(ref=referenceFileName, path_file=path_file)
(!(pRecord))
#if excecution gets here, no pathway was selected and the option to import new pathway sets will be provided
importPathways(symtab)
}
#load the actual paths
paths = loadPaths(pRecord=pRecord, symtab=symtab, firstGeneColumn=3)
(paths)
}#getPaths
getPathMetaData<-(ref="./reference_data/paths/pathMetaData.txt"){
((ref)){
fullTab = (=ref, header=T, sep="\t", stringsAsFactors=F, comment.char="")
}{
fullTab = (=("extdata/pathMetaData.txt",package = "packageDir"),
header=T, sep="\t", stringsAsFactors=F, comment.char="")
(!((path=ref))) (path=(path=ref), recursive=T, showWarnings=F)
(from=("extdata/pathMetaData.txt",package = "packageDir"), to=ref)
(from=("extdata/Reactome.2014.04.06.12.52.27.txt",package = "packageDir"),
to=((ref), "/Reactome.2014.04.06.12.52.27.txt"))
# write.table(x=fullTab, file=ref, quote=F, sep="\t", col.names=T, row.names=F)
(fullTab)
}
}
#choosePaths()
#allows interactive selection of pathways
#takes: ref: path meta data file (default:"./reference_data/paths/pathMetaData.txt")
# path_file: the file name of a pathway set (must be in GSEA format)
#returns: list: 1 row from the path meta data file given by the ref argument
choosePaths<-(ref="./reference_data/paths/pathMetaData.txt",
path_file=){
(!(path_file)) path_file =
outRecord =
fullTab = getPathMetaData(ref=ref)
((path_file)){
refTab=fullTab,(1:4)]
(refTab)<-(="[._]", replacement=" ",x=(refTab))
(refTab)<-("Repository source", "procurement date", "Number of paths", "Number of genes")
refTab = (refTab, 1:(refTab))
(refTab)[ncol(refTab)]<-"Selection Number"
("\nCurrently available pathway repositories:\n")
(refTab)
(T){
= (("Please enter the selection number for the pathway repository you would like to use.\n",
"Or, if you would like to import a different pathway repository, (for example: a custom repository)\n",
"please enter i: ", sep=""))
(=="i"){
()
} (!(())){
(()%in%1:(refTab)){
outRecord = fullTab[as.numeric(),]
}
}
("Sorry, that input was not recognized, please try again.")
}
}{
#check if the path file can be found..
outreci = (fullTab$==path_file)[1]
((outreci)) outreci = ((fullTab$Original.file.name)==(path_file))[1]
((outreci)) outreci = ((fullTab$)==(path_file))[1]
((outreci)){
(("Pathway file could not be found.... \nPath file:",path_file,"Files in pathway records:\n",(fullTab$,sep="\n",collapse="\n")))
(choosePaths())
}
outRecord = fullTab[outreci,]
}
(outRecord)
}#choosePaths
getPathObject<-(symtab, pRecord, pData, repset=){
#getPathObject
#creates path list object with all possible slots filled
# symtab: the hugo lookup table
# pData: the bipartate graph of pathways
# pRecord: the row of data from the path meta data file
# repset: depricated: the set of repositories
# pd=list()
pd = $()
pData = (pData)
(pData) = (="\\.",replacement="-",x=(pData))#replace periods with dashes
(pRecord$Gene_Id_Type=="HUGO"){
(pData) = (symbol_set=(pData),verbose=F,symref=symtab)
}{
("\nNon-HUGO symbols used in paths thus skipping gene symbol correction\n")
}
pd = setPathMetaData(pd=pd,
symbol_type=pRecord$Gene_Id_Type,
p.file=pRecord$,
symtab=symtab,
p.info = (pRecord$Set_name, pRecord$date_procured),
p.date = pRecord$date_procured,
p.source=pRecord$Set_name,
p.paths=pData)
pd$paths = pData
(pd)
}
setPathMetaData<-(symtab,
p.paths,
p.file,
p.date,
p.source,
p.name=,
p.info=,
pd=(),
symbol_type="HUGO"){
#setPathMetaData
#establishes path meta data
#takes: symtab: the symbol lookup table
# p.paths: the set of pathways in bipartate graph format
# p.file: the file containing the pathways
# p.date: the date the pathways were produced
# p.source: the source of the pathways
# p.name: optional, the name of the set of pathways
# (if not provided, a formatted version of the pathname will be used)
# p.info: identifying data about the pathways
# pd: a nascent path object; providing this will preserve any custom list slots that have already been filled.
#returns: path_details object with meta data filled in
p.info = ((p.info),
yes=p.file,
no=p.info)
p.name = ((p.name),
yes=(=" [0-9]*.[0-9]*.[0-9]* [0-9]*.[0-9]*.[0-9]*.txt", replacement="", x=(path=p.file)),
no=p.name)
p.paths = (p.paths)
pd"name"]] = p.name
pd"file"]] = p.file
pd"info"]] = p.info
pd"date"]] = p.date
pd"source"]] = p.source
pd"gene_overlap_counts"]] = (T,(p.paths))%*%p.paths
pd"full_path_length"]] = p.paths%*%(T,(p.paths))
pd"symtable"]] = symtab
pd"paths"]] = p.paths
pd"symbol_type"]] = symbol_type
pd"graphite"]] = ()
(pd)
}#setPathMetaData
#loadPaths
#used by getPaths()
#loads a set of paths from set of files at ./reference_data/paths/working_paths/
#will try to load path in GSEA format from the file denoted by file_name
#takes: paths pRecord: a list with one slot named "file" and other slots to be used by getPathObject
#returns: pathset
loadPaths<-(pRecord, symtab=, firstGeneColumn=3){
((symtab)) symtab = ()
tmptab = (pth=loadPathsAsSets(fname=pRecord$, firstGeneColum=firstGeneColumn, keepLowerCaseGenes=F))
("getting path object")
pathset = getPathObject(symtab=symtab, pRecord=pRecord, pData=tmptab)
("got path object..")
(pathset)
}
#imports set of pathways in bipartate format
#takes: optional: fname: the file name of the bipartate path set
# if fname is not provided, user will be prompted to select a set of pathways
#returns: bipartate graph of the path set
bipartateImport<-(fname=){
preped_paths = ()
("\nPlease select the tab-delimited file containing the pathway set you would like to use.\n")
("Pathways in file should be in bipartate graph format, with the columns the\n")
("colums the path members, the rows the path names and the values the values all\n")
("T or F indicating path membership.\n")
(T){
((fname)) fname = ()
("\nFile selected: ", fname, "\n")
preped_paths"paths"]] = (expr=(=fname, header=T, =1, sep="\t"), silent=T)
((preped_paths"paths"]]))
(!(fname)) ("Alert! the file provided could not be read as a bipartate graph.")
(!(fname)) ()
blnk=("Sorry, the file chosen could not be read as a bipartate graph.\nPlease press enter and select another file.")
}
("\nFile created on", ((fname)$mtime),"\n")
preped_paths"paths"]] = (=fname, header=T, =1, sep="\t")
preped_paths"paths"]] = (preped_paths"paths"]])
preped_paths"original_file_or_source"]] = fname
preped_paths"original_file_creation_date"]] = ((fname)$mtime)
(preped_paths)
}
#imports set of pathways in bipartate format
#takes: optional: fname: the file name of the GSEA format path set
# if fname is not provided, user will be prompted to select a set of pathways
#returns: bipartate graph of the path set
GSEAImport<-(fname=){
preped_paths = ()
("\nPlease select the tab-delimited file containing the pathway set you would like to use.\n")
("Pathways in file should be in GSEA format, with each row containing a path name followed by\n")
("path members, all sepparated by tabs.\n")
(T){
((fname)) fname = ()
("\nFile selected: ", fname, "\n")
preped_paths"paths"]] = (expr=(pth=loadPathsAsSets(fname=fname)),silent=T)
((preped_paths"paths"]])|(preped_paths"paths"]]))
(!(fname)) ("Alert! the file provided could not be read as a GSEA-format path set.")
(!(fname)) ()
blnk=("\nSorry, the file chosen could not be read as a GSEA-format path set.\n(Ex: <path_name><path_source><gene1><gene2>...)\nPlease press enter and select another file.")
}
("\nFile created on", ((fname)$mtime),"\n")
preped_paths"paths"]] = (preped_paths"paths"]])
preped_paths"original_file_or_source"]] = fname
preped_paths"original_file_creation_date"]] = ((fname)$mtime)
(preped_paths)
}#GSEAImport
#importPathways
#allows import of pathways either interactively or by choice and file name
#takes: symtab : table of symbol lookups
# choice: currently b for bipartate graph or g for GSEA format
# fname: the file name of the pathway repository
# must match the format indicated by the choice argument
#returns: nothing
importPathways<-(symtab=, choice=, fname=){
((symtab)) symtab = ()
#initialize
preped_paths = #will hold the forming path object
#allow selection
(T){
(!(choice)) (choice%in%("b", "g"))
choice = ("To provide a set of pathways in bipartate graph format, enter b\nTo provide a set of pathways in GSEA/'gene set' format, enter g: ")
}
#&open path file
# conversion to bipartate format
(choice == "b"){
preped_paths = bipartateImport(fname=fname)
} (choice=="g"){
preped_paths = GSEAImport(fname=fname)
}
usesHUGO = ("Does the selected pathway set employ HGNC/HUGO symbols? (y/n)")
(usesHUGO=="y"){
symbol_type="HUGO"
#correct gene symbols
(preped_paths$paths)<-(symbol_set=(preped_paths$paths),
symref=symtab,
verbose=T)
}{
symbol_type=("Please enter the type of gene identifiers employed (ex: UniProt) ")
symtab=
}
#setPathMetaData
preped_paths = manualPathMetaData(preped_paths=preped_paths, symtab=symtab, symbol_type=symbol_type)
#establish path set
recordPathSet(ps=preped_paths)
}#importPathways
#manualPathMetaData
#allows user to interactively establish path meta data
#used by importPathways() function
#takes: preped_paths: nascent path_detail object with file name and paths slots filled in.
# symtab: the symbol lookup table
#returns: preped_paths: the path_detail object with meta data filled in
manualPathMetaData<-(preped_paths, symtab=, symbol_type="HUGO"){
# library("tools")
((symtab)) symtab = ()
pathsFolder = "./reference_data/paths/"
# allSlots = c("paths","name", "file", "info", "date", "source", "gene_overlap_counts", "full_path_length", "symtable", "original_file_or_source", "original_file_creation_date")
neededSlots = ("name", "date", "source")
missingSlots = neededSlots!neededSlots%in%(preped_paths)]
addedSlots = ()
(T&(!neededSlots%in%(preped_paths))){
defaultName = (path=preped_paths$original_file_or_source)
defaultName = (=(".",(defaultName),"$",sep=""),replacement="",x=defaultName)
defaultName = (="[._]", replacement=" ", x=defaultName)
("source"%in%missingSlots) addedSlots"source"]] = (("Please enter the source of the pathway set.\nTo use the default,\"",defaultName,"\", just press enter: "))
("name"%in%missingSlots) addedSlots"name"]] = (("Please enter the name of the set of pathways just added: \nTo use the default,\"",defaultName,"\", just press enter: "))
("date"%in%missingSlots) addedSlots"date"]] = ("Please enter the date the new pathway set was generated or obtained\n(from the providers of the pathway repository): ")
(s (addedSlots)) (addedSlots[s]]=="") addedSlots[s]] = defaultName
("\nPlease verify these entries are all correct: \n")
(r (addedSlots)) ("For \"", r, "\" you entered: ", addedSlots[r]], "\n",sep="")
(("If the above entries are all correct, press enter. \nTo change one or more, enter \"c\"")=="")
}
preped_paths = (preped_paths, addedSlots)
fileName = (pathsFolder, preped_paths$," ", (="[-:]", replacement=".", x=(())),".txt", sep="")
("\nusing file name:", fileName,"to save the current set of pathways.\n")
preped_paths = setPathMetaData(symtab=symtab, symbol_type=symbol_type,
pd=preped_paths,
p.info=(preped_paths"name"]], preped_paths"date"]]),
p.name=preped_paths"name"]],
p.source=preped_paths"source"]],
p.paths=preped_paths$paths,
p.file=fileName,
p.date=preped_paths"date"]])
(preped_paths)
}#manualPathMetaData
#establishPathSet
# loads set of pathways into program memory
# a)gets meta data to put in path summary table
# b)corrects gene names
# c)makes GSEA format list of paths
# d)saves GSEA format paths file
# e)adds to path records file
test.recordPathSet<-(){
pathset=path_detail$paths
}
#
# #obtain path meta data
# ps = setPathMetaData(symtab=hugo,
# p.date=p.date,
# p.source=p.source,
# p.file=fileName,
# p.paths=pathset)
recordPathSet<-(ps){
pathset = ps$paths
#make gsea data structure
gsea = toGSEAformat(bip=pathset, psetName=ps$)
#save gsea data structure
(x=gsea, =ps$, =F, sep="\t", =F,col.names=F)
#add to pth records file
prec = addToPathRecords(pfull=ps)
}#recordPathSet
#addToPathRecords
#used by recordPathSet() function to add path repository meta data to
#takes: pfull: the completed path_detail object
# precFname: the name of the path meta data/path records file
#returns: prec
addToPathRecords<-(pfull, precFname = "./reference_data/paths/pathMetaData.txt"){
# data: path set name
# date
# number of paths
# number of genes
# number of paths with zero nodes
# number of unapproved hugo symbols
# pathFileName
#open file
prec =
((precFname)) prec = (=precFname, header=T, sep="\t", ="", stringsAsFactors=F)
#assemble input data
newRow = (pfull$,
pfull$,
(pfull$paths),
(pfull$paths),
(pfull$full_path_length==0),
(!(pfull$paths)%in%pfull$symtable$Approved.Symbol),
pfull$,
pfull$original_file_or_source,
pfull$original_file_creation_date,
pfull$symbol_type)
prec = (prec, newRow)
(prec)<-("Set_name",
"date_procured",
"Number_of_pathways",
"Number of genes in pathways",
"Number of paths with zero genes",
"Number of path nodes with out HUGO symbols",
"file",
"Original file name",
"Path_Import_Date",
"Gene_Id_Type")
#check if more than one record use the same source path file
prec=checkMultipleRecordsOneRepository(prec=prec)
#write file
(x=prec, =precFname, =F, sep="\t", col.names=T, =F)
(prec)
}#addToPathRecords
#checkMultipleRecordsOneRepository()
#assures multiple records weren't produced from single repository files
#takes: the set of pathway records to be written to disk
#returns: the pathway records with any unwanted records removed
checkMultipleRecordsOneRepository<-(prec){
prec2=(prec)
bnames = (prec2"Original file name"]])
duprows = (((x=bnames)), ((x=bnames, fromLast=T)))
(!(duprows)) (prec2) #return if there's nothing to correct
#see if everything except the file name is the same
dupTest = prec,!(prec)%in%"file"]
((dupTest)>((dupTest))){
remFile = prec$[duplicated(prec$,fromLast=T)]
prec2 = prec2!(prec$,fromLast=T),]
(remFile)
}
reprows = prec2[duprows,]
(reprows)<-1:(reprows)
("\nIt appears the same set of pathways was imported multiple times.\n",
"(repeated \"original file name\"s were found)\n",
"These are the records that appear to be repeated:\n")
(reprows)
(T){
= ("Please enter the line numbers of the path sets you would like to keep\n(sepparated by spaces): ")
= (x=, =" ")[1]]
= ((), silent=T)
(!(())&!(!%in%1:(reprows)))
("Error, input not recognized, please try again.")
}
deleteRows = duprows!1:(duprows)%in%]
#pull out the gsea file to be deleated
todel = prec2[deleteRows,]
#remove the GSEA file that was created
(todel$)
# now limit prec2 to the selections
prec3 = prec2!1:(prec2)%in%deleteRows,]
(prec3)
}#checkMultipleRecordsOneRepository()
test.toGSEAformat<-(){
pset = (STUDY@studyMetaData@paths$info)
gsea = toGSEAformat(bip=STUDY@studyMetaData@paths$paths, psetName=pset)
}
toGSEAformat<-(bip, psetName){
out = (="", = (bip), =3)
(r 1:(bip)){
pname = (bip)[r]
genes = ( (bip)[bip[r,]] , collapse="\t")
out[r,] = (pname, psetName, genes)
}
(out)
}
#loadPathsAsSets()
# loads paths from file in this format: <path name><\t><gene><space><gene>....
# everything after and including "firstGeneColum" is assumed to be a gene"
#the first column is assumed to be the path name or other unique identifier for the pathway
#
#takes: firstGeneColum: the first tab stop with gene names in it, usually 2 or 3
# keepLowerCaseGenes: keep genes who are all lower case
loadPathsAsSets<-(firstGeneColum=,
keepLowerCaseGenes=,
fname="./reference_data/paths/ReactomePathways_dl_June_4_13.gmt.txt",
verbose=T){
("\nLoading pathway sets..\n")
#figure out what happens when this gets passed a bipartate graph..
pset = (fname, sep="\n",comment.char="", stringsAsFactors=F, ="")
outlist = ()
excluded = ()
( (firstGeneColum) ){
rone = (pset[1,],="\t")[1]]
rone = (=rone, =1, =(,"row contents"))
("\nThis is the contents of the first row from the file:\n")
(rone[1:((rone),8),,=F])
(T){
uin = ("Please enter the index of the first gene identifier: ")
uin = (uin)
(!(uin)){
firstGeneColum = uin
}
("Error, please enter a number.")
}
}
((keepLowerCaseGenes)){
keepLowerCaseGenes = "Y"==(("Would you like to keep lowercase genes?\n(note: lower case likely indicates gene is from a virus)\nPlease enter y or n"))
}
( i 1:(pset) ){
("\nWorking on row",i,"of",(pset),"rows.\n")
currow = (pset[i,],="\t")[1]]
("\nGenes in current row:",((currow)-firstGeneColum+1))
pname = currow[1]
("\nPath name:",pname,"\n")
((currow)>2){
genes = currow[firstGeneColum:(currow)]
(!keepLowerCaseGenes){
#check gene names to see if any are all lower case
badInds = !(="[A-Z]+",x=genes)
((badInds)){
(genes)
(verbose) ("Lowercase gene names found. \nFrom:\n",pname,
"\nexcluding",(badInds),"gene(s):\n",
genes[which(badInds)],
"\nindexes:",(badInds),"\n")
(" .. ")
excluded = (excluded, genes[badInds])
(" .. ")
}
(".")
genes = genes!badInds]
(".")
}
("*")
outlist[pname]] = genes
("*")
}{
("\nThe pathway \n\"",pname,"\"\nwas found not to have any genes under the currently employed type of gene symbols.\n", sep="")
}
}
excluded = (excluded)
((excluded)){
("\nExcluded gene symbols:\n")
(excluded)
}
(outlist)
}#loadPathsAsSets
####getTargetMatrix
#get a bipartate graph indicating which genes are targets in which pathways
#paths with no targets are removed
#takes tgenes: vector of targeted genes
# paths: set of pathways in logic matrix/bipartite graph format
#returns: bipartate graph: rows = paths columns = genes
getTargetMatrix<-(tgenes, paths){
("\ntargetMatrix(): ")
nmat = paths #make copy of the paths
tcols = (paths)%in%tgenes#find how many of the path genes are in the set of target genes
("genes from current set that are found the pathways: ", (tcols),"\n")
((tcols)==0){
("There are no target genes in the set of paths")
()
} ((tcols)==1){#if there is only one gene in the target genes that is in the set of paths
("There is only one target gene in the set of paths")
#if there is only one gene, make a one column matrix out of it
colreduced = nmat,tcols]
rowreduced = (colreduced[colreduced],=1)
(rowreduced) = (paths)[tcols]
(rowreduced)
}
# tcol = rep(T, nrow(nmat))%*%nmat>0
nmat = nmat,tcols, = ]
trow = nmat%*%(T, (nmat))>0
nmat = nmat[trow, , = ]
(nmat)
}
#'@title Get the names of the pathways containing submitted genes.
#'@description Given a vector of gene identifiers, this function will return a two column matrix: the gene ids and a column indicating which, if any pathways each gene belongs to.
#'@param genevector A character vector of gene identifiers
#'@param path_detail A \code{Path_Detail} object.
#'@return A two column \code{data.frame}, with columns "gene" and "path_membership"
#'@export
#'@examples
#'path_detail = getDefaultPaths()
#'res = getPathsWithGenes(genevector=c("FAT1","TP53","ENAH","HOMER3"), path_detail=path_detail)
<-(genevector, path_detail){
((genevector)|(genevector)){
genevector = genevector,1,=T]
}
membership = genevector%in%(path_detail$paths)
inpaths = genevector[membership]
pnames = ()
number_of_pathways = ()
(g inpaths){
pathsWithGene = (path_detail$paths)[path_detail$paths,g]]
pnames = (pnames, (pathsWithGene, sep="; ", collapse="; "))
number_of_pathways = (number_of_pathways, (pathsWithGene))
}
(pnames)<-inpaths
genematrix = ( gene=genevector,
number_of_pathways = (0,times=(genevector)),
path_membership = ("Not in pathways", times=(genevector)),
stringsAsFactors=)
(genematrix)<-genematrix,1]
genematrix[names(pnames),"path_membership"] = pnames
genematrix[names(pnames),"number_of_pathways"] = number_of_pathways
(genematrix)<-
(genematrix)
}
#appends column indicating if the gene is in the path repository
GenesInPaths<-(ordgenelist, path_detail){
("inside GenesInPaths()")
(path_detail$info)
pmemb = ("No", times=(ordgenelist))
((ordgenelist,1])){
pmemb[rownames(ordgenelist)%in%(path_detail$paths)]="Yes"
ordgenelist = (ordgenelist, pmemb)
(ordgenelist)<-("Nubmer of patients","In the current set of pathways?")
(ordgenelist)
}
pmemb[ordgenelist,1]%in%(path_detail$paths)]="Yes"
ordgenelist = (ordgenelist, pmemb)
(ordgenelist)<-("Gene","Number of patients","In the current set of pathways?")
(ordgenelist)
}
#'@title Get gene identifiers from a list of pathways.
#'@description Takes a list of path identifiers and returns a list of identifiers for the genes in the submitted pathways.
#'@param pids Character vector of cellular pathway identifiers.
#'@param STUDY A study object, a Path_Detail object or a bipartate graph in a "matrix" class R object.
#'@return Character vector of gene identifiers.
#'@export
<-(pids, STUDY){
((STUDY)=="Study"){
pths = STUDY@studyMetaData@paths$paths
} ((STUDY)=="Path_Detail"){
pths = STUDY$paths
} ((STUDY)=="matrix"){
pths = STUDY
}{
("Unrecognized path source. Must be a matrix, or a Study or Path_Detail object")
}
(!(pids)) pids=pids,1,=]
missingP = !pids%in%(pths)
((missingP)){
("Could not find these pathways:\n", (pids[missingP], collapse="\n ", sep="\n "))
}
pthsub = pths[pids,,=]
geneIdentifiers = (pthsub)[rep(, times=(pthsub))%*%pthsub>0]
geneIdentifiers = (geneIdentifiers)
(geneIdentifiers)
}
#getPathsForGenes()
#gives quick report on set of pathways an input set of genes belongs to
#takes: glist: a vector of genes symbols
#path_detail: a path_detail object
getPathsForGenes<-(glist=, path_detail=path_detail, verbose=F){
((glist)){
glist = (("Enter set of genes sepparated by spaces: "),=" ")[1]]
}
tm = getTargetMatrix(tgenes=glist,paths=path_detail$paths)
#add row and column totals
gsums = (T,(tm))%*%tm
psums = tm%*%(T, (tm))###note, this contains the grand total
#reorder it all
gorder = (gsums,decreasing=T)
porder = (psums,decreasing=T)
tm = tm[porder,gorder,=F]
gsums = gsums[gorder]
psums = psums[porder]
tm = (tm,gsums)
tmppsums = (psums,(psums))
tm = (tm,tmppsums)
#correct marginal's names
(tm)[nrow(tm)] = "Marginal total: paths gene belongs to"
(tm)[ncol(tm)] = "Marginal total: input genes found in path"
(verbose){
("\nGenes from input list not found in current pathway repository:\n")
((x=glist, (path_detail$paths)))
("\nSummary of genes from input list found in current pathway repository:\n")
(tm)
}
(tm)
}#getPathsForGenes()
#getFrequencyMatrix
#takes: genecounts: 2 column table <gene_names> <counts across cohort of patients inwhich that gene is active>
# paths: pathway target matrix in bipartate logic matrix form
#returns: frequencymatrix: pseudo bipartate graph: rows: path names,
# columns: genes,
# values: number of times gene is found mutated across cohort
#getFrequencyMatrix(genecounts=genesum, paths=somtargetmatrix)
getFrequencyMatrix<-(genecounts, paths){
(genecounts)<-genecounts,1] #make it like a dictionary so that gene extraction can be done
inpath = genecounts[colnames(paths),] #extract only the genes which are actually in the pathways
out =
(i 1:(paths)){#for each path
= inpath,2]*paths[i,]#multiply every value in the i'th path by the count for that gene. . . NO!
out = (out, )
}
(out)<-(paths)
(out)
}
#getMutCountMatrix: makes pseudo bipartate graph containing number of times gene is found mutated across cohort
# as the values for each member in each pathway
#takes: genecounts: 2 column table <gene_names> <counts across cohort of patients inwhich that gene is active>
# paths: pathway target matrix in bipartate logic matrix form
#returns: frequencymatrix: pseudo bipartate graph: rows: path names,
# columns: genes,
# values: number of times gene is found mutated across cohort
#getFrequencyMatrix(genecounts=genesum, paths=somtargetmatrix)
getMutCountMatrix<-(genecounts, paths){
(" In getMutCountMatrix()")
(genecounts)<-genecounts,1] #make it like a dictionary so that gene extraction can be done
inpath = genecounts[colnames(paths),] #extract only the genes which are actually in the pathways
out = (=0,=(paths),=(paths),=(rn=(paths),cn=(paths)))
(gene (paths)){#for each path
= inpath[gene,2]*paths,gene]#multiply every value in the i'th path by the count for that gene. . . NO!
out,gene] =
}
(out)<-(paths)
(out)
}##getMutCountMatrix:
test.prepPathListForSave<-(){
matrixToSave = (gset=reactome)
(matrixToSave[1:2])
(x=matrixToSave, ="./reference_data/paths/GraphiteReactomeUniprotFormat.txt",
=F,
=F,
col.names=F)
pset = ("./reference_data/paths/GraphiteReactomeUniprotFormat.txt",
header=F,
sep="\n",
comment.char="",
stringsAsFactors=F,
="")
pset = (=pset,1], =1)
checkEquals(target=matrixToSave, current=pset)
GSEAImport(fname="./reference_data/paths/GraphiteReactomeUniprotFormat.txt")
}
#'@title prepPathListForSave()
#'@description preps paths in gene-set format (list with names=path names and values = vector of gene names) for saving as a tab delimited file
#'@param gset The list object of paths
#'@param path_source A string desribing the source of the pathways, to be pasted as the second element in each line of the output in the manner as "Reactome pathway" is seen here:
#' Abacavir metabolism Reactome pathway UniProt:P49902 UniProt:Q16774
#'@return one column matix, with each row containing the tab-delimited set: path name, path source and path members
#'@export
#'@import graphite
#'@examples
#'library(graphite)
#'plistsAsDataFrame = prepPathListForSave(gset=reactome, path_source="Reactome pathway")
<-(gset, path_source="Reactome pathway"){
(gset)
pasted =(X=gset, FUN=(x){(nodes(x), sep="", collapse="\t")})
dout = ()
(i 1:(pasted)){
dout = (dout,((pasted)[i], path_source, pasted[i], sep="\t"))
}
((=dout, =1))
}#prepPathListForSave
