(!("VERBOSE")) VERBOSE = F
.pardefault <- (no.readonly = T)
.parpin <-()$pin
# .pardefault <- par()
stackedGeneBar<-(tcga_som,
="Mutation types for the top 20 most mutated genes",
colorMatchesFile="./reference_data/MAFcolorMatches.txt"){
ufilt = tcga_som
gs = (=tcga_som,"Hugo_Symbol"], display=F)
top20 = (gs[order(gs,2],decreasing=T),],20)
gsnames = top20$types
toprows = ufilt[ufilt$Hugo_Symbol%in%gsnames,]
slimrows = toprows,("Hugo_Symbol","Variant_Classification")]
slimrows = (x=slimrows, y=top20, by.x="Hugo_Symbol", by.y="types")
slimrows = slimrows[order(slimrows$counts),]
slimrows$Hugo_Symbol<-(slimrows$Hugo_Symbol, =top20$types)
colcols = getColorSequence(cnames=(slimrows$Variant_Classification),
fname=colorMatchesFile)
= ()[colcols]
()<-(colcols)
p1 = ggplot(=slimrows, aes(x=Hugo_Symbol,
fill=(Variant_Classification)))+
geom_bar(color="black")+
coord_flip()+
scale_fill_manual(values=)+
theme_bw()+
theme(legend.title=element_blank())+
ggtitle()
(p1)
}
simpleGGHist<-(dataSet, xlab, ylab, mainTitle, showPlot=T){
(ggplot2)
p = qplot(x=(dataSet), geom="histogram")+
theme_bw()+
xlab(xlab)+
ylab(ylab)+
ggtitle(mainTitle)
(showPlot){
(p)
} {
(p)
}
}
promptNumeric <- function () {
(T){
=()
(=="")
= (expr=())
(!())
("Sorry, your input could not be understood, please try again")
}
()
}
longTextBarPlot<-(, lab, main=""){
bpdata = (, horiz=T, main=main)
xmin = ("usr")[1]
xrange = ("usr")[2] - ("usr")[1]
xcoord = xmin + (xrange/40)
nlev = ()
(x=xcoord, y=bpdata, =lab, pos=4)
}
###############################################################################
#
# helperFunctions.R: This file contains the all helper functions not directly related to any other source file.
# author: Frank Kramer <dev@frankkramer.de>
#
# This is released under GPL-2.
#
# Documentation was created using roxygen
#
###############################################################################
#' Replace factors/levels in a data.frame and use plain strings instead
#'
#' This function takes a data.frame as argument and returns it with strings instead of factors.
#'
#' @param df any data.frame with factor levels in at least one column
#' @return The data.frame is returned using strings instead of factors.
#' @author Frank Kramer
<- function () {
(!("data.frame" %in% ())) { ("Error: unfactorize: data.frame as argument expected") }
for (i 1:()) {
if ((,i]) == "factor") ,i] <- (,i])
}
}
save.plot<-(pname){
(AutoSavePlot(pname))
}
AutoSavePlot<-(pname){
((pname)){
pname = "tmp_image_file"
}
(!(=".png$", x=pname)) pname = (pname,".png",sep="")
pfname = plotFileName(pname)
(,pfname)
()
(pfname)
}
plotFileName<-(pname, addTime=F){
#takes a name, makes sure there's a directory there, adds time stamp to make it unique, adds a postfix to make sure
autoincrament = F
root="./output/imageTemp/"
((="[/]",x=pname)){#the a file path was passed, find the root
root = (((x=pname, ="[/]")[1]][1:(((x=pname, ="[/]")[1]])-1)],collapse="/"),
"/",
sep="")
pname = (x=pname, ="[/]")[1]][length((x=pname, ="[/]")[1]])]#the last part of the path is the actual file name
}
(path=root, recursive=T, showWarnings=F)
ctime = ""
(addTime) ctime = (=" ", replacement=".", x=((()),"."))
fullPath = (root, "graphic", (="[-:]", replacement=".", x=ctime), pname, sep="")
#check if it exists, if it does, append a number
((fullPath)&autoincrament){
pname = ("1", pname, sep=".")
fullPath = (root, "graphic",(="[-:]", replacement=".", x=(())),".", pname, sep="")
}
(fullPath)
}
#pgvmFromStacked makes a patient gene value matrix from data in stacked format
#rows = genes, columns = patients, values = value for gene, ie, for a gene which has a missense mutation,
# the value would be "missense"
#if more than one value are given for a particular patient-gene combination
#takes: stackedData: data frame with patient gene data in stacked format
# valueCol: the name of the column with the gene values
# patientIdCol: the name of the column with patient ids
# geneCol: the name of the column with the gene identifiers
# blankSymbol: what is placed in the pvm when the stacked file does not give a value
# for a given patient-gene combination
#returns: patient gene value matrix
pgvmFromStacked<-(stackedData, valueCol, patientIdCol, geneCol, blankSymbol=){
#1: find list of unique pids
pids = (stackedData,patientIdCol])
#2: find lis of unique genes
ugenes = (stackedData,geneCol])
#case to handle: if patient does not have a value for a gene,
#it will be set to NULL
pvm = ( = blankSymbol, =(pids), =(ugenes), =(ugenes, pids))
#3: put values
(p pids){
prows = stackedData[stackedData,patientIdCol] == p,]
(g ugenes){
pgrows = prows[prows,geneCol]==g,]
((pgrows) == 1){
pvm[g,p] = pgrows[1,valueCol]
} ((pgrows) > 1){
#case to handle: if there are more than one value for a gene..
("\nNote: multiple values found for this patient-gene combo:", p, g, "\n")
pvm[g,p] = (pgrows,valueCol], collapse=" ")
}
}
}
(pvm)
}
PGMFromVector<-(genevector){
((genevector)) genevector = (genevector)
((genevector)) genevector = (genevector)
pgmout = (=T,
=(genevector),
=1,
=(genevector, "Genes_in_set"))
(pgmout)
}
#selectionList()
#prompts the user with a numbered list of selection options
#takes a vector of options to choose from (does not have to be unique set)
#returns logic vector indicating which items in the input vector match the/those
# item(s) chosen by user
selectionList<-(valcol, verbose=T){
usel = (valcol)
useltmp = (="_",replacement=" ",x=usel)
uselmat = (=1:(usel), =1, =(useltmp, "selection number"))
(uselmat)
uin = ("Please enter the number(s) corresponding to your selection(s)\n")
uin = ((x=uin, =" ")[1]])
lvout = valcol%in%usel[uin]
(lvout)
}
filePrompt<-(defaultfile){
#prompts user to select file for data input
#provides default file option
#returns file name
fsel = (("\nTo select a data file, Enter s\n",
"To load the default data from \n",defaultfile,",\njust press enter \n",sep=""))
(fsel=="s"){
pfile = ()
}{
pfile = defaultfile
}
("\nLoading data from:\n",pfile,"\n")
(pfile)
}
#'@title openPGM
#'@description Opens a patient gene matrix file.
#'@param fname The file path to the patient gene matrix file.
#'@return The patient gene matrix (logic matrix with gene identifiers given as row names and patient identifiers gien as column names)
#'@export
<-(fname = ){
#opens a patient gene matrix file
#takes: 1) file name
# or 2) no arg (in this case, user will be prompted to select a file)
#returns: patient gene matrix
((fname)){
fname = ()
}
pgmdat = (=fname,check.names=T,
header=T,
=1,
sep="\t",
stringsAsFactors=F)
pgmdat = (pgmdat)
(pgmdat)
}
#extract_pid_w_matchnormal
#like extract_pid but leaves the sample type code at the end (match/normal)
#takes: n: a string that might contain a TCGA id, ex: TCGA-IQ-7632-01 in unc.edu__IlluminaHiSeq_RNASeqV2__TCGA-IQ-7632-01A-11R-2081-07__expression_rsem_isoforms_normalized.txt
#returns: the TCGA id (with -'s replaced by .'s) if an id is found; NULL if there is not a TCGA id in the input string
<-(n=)
{
((n)){n="unc.edu__IlluminaHiSeq_RNASeqV2__TCGA-IQ-7632-01A-11R-2081-07__expression_rsem_isoforms_normalized.txt"}
#first replace all special characters with spaces
#en = gsub(pattern="[[:punct:]]", replacement=" ", x=n)
sen = (n, ="[[:punct:]]")[1]]
#find where "TCGA" is
i = (="TCGA", x=sen)
((i)<1){
()
}
#take out the TCGA part, plus the next two parts
= sen[c(i, i+1, i+2, i+3)]
[4] = (="[a-zA-Z]*", "", x=[4])
out = (, sep=".", collapse=".")
(out)
}
toPGMWithCoverage<-(sds, coverageSet){
pgm = toPGM(sds)
additionalRowNames = (coverageSet, (pgm))
newChunk = (=F,
=(additionalRowNames),
=(pgm),
=(additionalRowNames, (pgm)))
}
##toPGM
# converts a stacked format to a patient gene matrix
#
#takes: a data frame with columns "Symbol" and "Ids"
#returns: patient gene matrix: rownames = gene names,
# colnames = patientIDs
# cell values: logical, T if gene is active in that patient
toPGM<-(sds){
#find unique patient ids
pids = (sds$Ids)
pids = (pids)
sds$Ids = (sds$Ids)
#find unique genes affected
syms = (sds$Symbol)
#make the matrix
out = (=F, =(syms), =(pids), =(syms,pids))
#for each patient, find the set of genes
(p pids){
#get the current patient's active genes
cursyms = sds$Symbol[sds$Ids == p]
#assign the cells for that patient, for those genes
out[cursyms, p] = T
}
(out)
}
# toBipartateGraph
# takes a data frame with two columns
# makes a bipartate graph from the two columns
toBipartateGraph<-(dfin){
(dfin)<-("Ids", "Symbol")
#find unique patient ids
pids = (dfin$Ids)
pids = (pids)
dfin$Ids = (dfin$Ids)
#find unique genes affected
syms = (dfin$Symbol)
#make the matrix
out = (=F, =(syms), =(pids), =(syms,pids))
#for each patient, find the set of genes
(p pids){
#get the current patient's active genes
cursyms = dfin$Symbol[dfin$Ids == p]
#assign the cells for that patient, for those genes
out[cursyms, p] = T
}
(out)
}
paths<-(smd){
(smd@paths)
}
#rounds columns in table to "figs" number of significant digits
#and replaces underscores in col names with spaces
cleanTables<-(tab,figs=3, verbose=F){
("\nCleaning table..")
(!((tab)|(tab))){
("WAIT, that's not a table..\n")
(tab)
}
((tab)){
(i 1:(tab)){
ccol = tab,i]
cname = (tab)[i]
(verbose){
(cname)
(tab[1,i])
((ccol))
((ccol))
}
((tab,i])%in%("integer","double")){
tab,i] = (tab,i],digits=figs)
}
}
}{
("table has no rows")
}
((tab)){
#replace underscores:
(tab)<-(="[_]", replacement=" ", x=(tab))
}{
("table has no columns")
}
("cleaned\n")
(tab)
}
#qinfo: gives you quick info about a data structure
#takes: ds: vector, dataframe or matrix
qinfo<-(ds,mx=25){
((ds))
((ds))
((ds))
((ds)){
((ds))
((ds))
} ((ds)|(ds)){
((ds))
((ds)<5){
((ds))
((ds))
}{
(ds[1:((ds),5),1:4])
((ds,1:4]))
((ds)[1:((ds),mx)])
((ds)[1:((ds),mx)])
}
} (!(ds)){
(ds)
} ((ds)){
((ds)[1:(mx,(ds))])
(i 1:(ds)){
((ds)[i])
qinfo(ds[i]])
}
}
}
#extract_pid
#takes: n: a string that might contain a TCGA id, ex: TCGA-IQ-7632-01 in unc.edu__IlluminaHiSeq_RNASeqV2__TCGA-IQ-7632-01A-11R-2081-07__expression_rsem_isoforms_normalized.txt
#returns: the TCGA id (with -'s replaced by .'s) if an id is found; NULL if there is not a TCGA id in the input string
<-(n=)
{
((n)){
n="unc.edu__IlluminaHiSeq_RNASeqV2__TCGA-IQ-7632-01A-11R-2081-07__expression_rsem_isoforms_normalized.txt"
("\nNote: extracting example pid from\n", n, "\n")
}
#first replace all special characters with spaces
#en = gsub(pattern="[[:punct:]]", replacement=" ", x=n)
sen = (n, ="[[:punct:]]")[1]]
#find where "TCGA" is
i = (="TCGA", x=sen)
((i)<1){
("** Note: TCGA id not found, returning barcodes as pids.")
(n)
}
#take out the TCGA part, plus the next two parts
= sen[c(i, i+1, i+2)]
out = (, sep=".", collapse=".")
(out)
}
test.summarize_by<-(){
tset1 = (("thirty", times=30), ("twenty", times=20), ("not ten", times=11))
tset1res = (=tset1, display=T, barPlotTitle="test summarize by", left_margin_factor=1.1)
checkEquals(checkNames=F,
current=tset1res,
target=(stringsAsFactors=F,
types=("not ten", "thirty","twenty"), counts=(11,30,20)))
}
###summarize_by
#summarize table by categories in a particular column
#used for making a quick summary of different types of somatic mutations in TCGA somatic mutation tables
#takes: the table
# the column to be summarized
# whether or not you want to display the summary to screen and bar graph
#returns: a summary table: category column, count column
#'@title summarize_by
#'@description summarizes the data levels in a column of categorical data and optionally provides a barplot output
#'@param col Vector, the categorical values to be summarized
#'@param display A flag indicating if the data should be displayed in a bar plot.
#'@param barPlotTitle The title that should be given to the bar plot
#'@param left_margin_factor A numerical factor expanding the left margin of the bar plot. Used to adjust in situations where there is insufficient room to display the text of each bar's label
#'@return A two column data.frame with columns "types", the different levels and "counts", the different numer of times each level occurs
#'@import ggplot2
<-(, display=F, barPlotTitle="Counts across types", left_margin_factor=1)
{
ctab = ()
out = ((ctab), (ctab), stringsAsFactors=F)
(out)<-("types", "counts")
out = out[order(out$types),]
out$counts = (out$counts)
types = out$types
counts = out$counts
(display){
p = ggplot(out, aes(x=types, y=counts))+
geom_bar(stat='identity')+
coord_flip()+
ggtitle(barPlotTitle)
(p)
# oldmar <- par()$mar
# while(T){
# res = try({
# par(mar=c(5.1, max(4.1,max(left_margin_factor*nchar(types))/2.5) ,4.1 ,2.1))
# # try(displayGraph(w), silent=T)
# barplot(counts, horiz=T, las=2, main = barPlotTitle, xlab="Number found in data set", names.arg=types)
# par(oldmar)
# }, silent=T)
# # if(!grepl(pattern="Error", x=res)) break
# if(is.null(res)) break
# par(oldmar)
# readline(prompt="There seems to have been an error with plotting the bar graph.\nPlease increase the size of the plot window, the press enter")
# }
# par(oldmar)
# print(out)
}
(out)
}
#join_tables
#takes: orig: a table to which columns from new_t are to be joined (if this is null t2 will be returned with adjusted names)
# new_t: joins columns from new_t to orig, filling 0s in where new_t doesn't have a value
# name_prefix: prefix to be appended to columns from new_t
# fill: value to be placed where new_t contains no value (typically 0 or NA)
#returns: t1 columns joined to t2 columns, new rows filled with the value in "fill"
joinTables<-(t1,t2, name_prefix=, fill=0){
("inside joinTables()")
((name_prefix)){#invent a name prefix
name_prefix = (((t1), (t2)))
}
#change t2 column names
t2n = (t2)
newt2n = (name_prefix, t2n, sep=".")
(newt2n)
(t2)<-newt2n
((t1)){
(t2)
}
allrownames=((t1), (t2))
#make new chunk
newchunkrows = (allrownames, (t2))
tmpnew = ((=(t2), =(newchunkrows), =fill), =newchunkrows)
(tmpnew)<-(t2)
newchunk = (t2, tmpnew)
rowstoold = (allrownames, (t1))
tmpold = ((=(t1), =(rowstoold), =fill), =rowstoold)
(tmpold)<-(t1)
updatedold = (t1, tmpold)
#join the old and the new
out = (updatedold, newchunk[rownames(updatedold),])
("\nTables Joined.. \n")
(out)
}
#getEmpP
#get empirical p-value
getEmpP<-(path,pgm,reps=1000,paths){
#some number of genes
#draw random genes for pathway
#for each gene: see if it's in the pgm
# if it is, set it to the score from the pgm
# if it's not, set it to zero
#score = sum for all genes in pathway
#add score to set of scores for ECDF
flatpgm = pgm%*%(T, (pgm))
ngenes = (path) #the number of genes in the path
vals = (0,times=reps)
(i 1:reps){
res = (x=0,times=ngenes) # results from this iteration
#select random genes for pathway
curPath = ((paths),size=ngenes,=F)
curInSomatic = (pgm)%in%curPath
#now get the scores for the genes that are actually in the pathway
curScores = flatpgm[curInSomatic] #extract the scores for all the genes in the current randomized pathway
curScore = (curScores)
vals[i] = curScore
}
genesInPath = (flatpgm)%in%path
trueScore = (flatpgm[genesInPath])
curFunction <- (vals)
pout = 1-curFunction(trueScore)
(pout)
}#getEmpP()
#goal: provided this set of data for each piece of the project:
# Panel used: Described in Tyner et al Oct. 2012*
# Total number of unique genes targeted by drug or mutation : 290
# Unique genes from targeted genes found in the provided set of pathways : 180
# genes targeted, but not found in the current set of pathways : 110
# Pathways used : KEGG PATHWAYS / 2011-Mar-14**
# Number of unique genes found in all pathways : 5894
# Total number of unique pathways examined : 229
# Number of unique pathways targeted by aberrations or drugs : 98
# List of genes targeted, but not found in current set of pathways
# AAK1 ACVRL1 ADCK3 ADCK4 ALK ANKK1 NUAK1 AURKB AURKC AXL BMP2K BLK BMX BRSK1 BRSK2 CAMK1 CAMK1D CAMK1G CDK11B CDK11A CDK19 CDK3 CDK8 CDK9 CIT CLK1 CLK2 CLK3 CLK4 DCLK1 DCLK2 DCLK3 DDR1 DDR2 DMPK CDC42BPG STK17A STK17B DYRK1B MAPK6 MAPK4 MAPK15 FRK GAK LATS1 LATS2 STK10 LTK MAP4K5 MARK1 MARK2 MARK3 MARK4 MELK MERTK MAP3K9 MAP3K10 CDC42BPA CDC42BPB MST1 STK24 MST4 MUSK MYO3A MYO3B STK38L NEK1 NEK2 NEK5 NEK6 NEK7 NEK9 CDK16 CDK17 CDK18 CDK14 PIM3 PKN1 PKN2 PLK3 PLK4 PRKD1 PRKD2 PRKD3 PTK6 RIOK3 ROS1 MYLK4 NUAK2 SIK1 SIK2 SRMS SRPK1 SRPK2 STK16 STK33 TESK1 TIE1 TLK1 TLK2 TNIK TNK1 TNK2 TNNI3K TSSK1B TYRO3 STK32B STK32C STK25 ZAK
#
startHTMLPlugIns<-(){
(!("hwriterPlus")){
("Trying to install hwriterPlus so that HTML output can be generated")
("hwriterPlus")
(("hwriterPlus")){
("hwriterPlus installed and loaded")
} {
("could not install hwriterPlus")
}
}
(!("hwriter")){
("Trying to install hwriterPlus so that HTML output can be generated")
("hwriter")
(("hwriter")){
("hwriter installed and loaded")
} {
("could not install hwriter")
}
}
(!("xtable")){
("Trying to install hwriterPlus so that HTML output can be generated")
("hwriterPlus")
(("xtable")){
("xtable installed and loaded")
} {
("could not install xtable")
}
}
}
WriteSettings<-(sumset, p){
dwu = sumset"settings"]]
((dwu)) ()
((dwu)&!(dwu)) dwu = saveSettings(set=dwu)
("Inside WriteSettings()... ")
#WriteDataWorkUp()
#builds the data work up table and outputs it to the HTML page
#takes: sumset: the summary set
# p: the handle for the html page being created
#returns: nothing
hwrite(("<A name=\"section",".002","\">","Settings used for this data set", "</A>", sep=""),p,heading=2)
hwrite("<A href=\"#section0.1\">Back to top</A><BR>", p, heading=4)
("\nWriting settings to HTML..\n")
hwrite(dwu, p, table.style= 'border-collapse:collapse', table.cellpadding='5px', col.names=F)
}
WriteDataWorkUp<-(sumset, p, fpath){
dwu = sumset"Data_work_up_notes"]]
((dwu)|!(dwu)) ()
(!(dwu)&(dwu)){
dwu = (set=dwu, fpath=fpath)
}
("Inside WriteDataWorkUp()... ")
#WriteDataWorkUp()
#builds the data work up table and outputs it to the HTML page
#takes: sumset: the summary set
# p: the handle for the html page being created
#returns: nothing
hwrite(("<A name=\"section",".001","\">","Data work-up notes", "</A>", sep=""),p,heading=2)
hwrite("<A href=\"#section0.1\">Back to top</A><BR>", p, heading=4)
("\nOutputting data work up notes.\n")
((dwu))
((dwu))
((dwu)|(dwu)){
dwu_tab = dwu
}{
(dwu)
dwu_names = (dwu)
dwu_vals =()
(i 1:(dwu)){
((dwu[i]])){
((dwu[i]][1]])==1&!(dwu[i]])){
# dwu_names = c(dwu_names,names(dwu)[i])
dwu_vals = (dwu_vals, dwu[i]][1]])
}
}
}
("bout to bind the data fram")
dwu_tab = (dwu_names, dwu_vals,stringsAsFactors=F)
("bound")
}
("\nwriting work up notes to HTML..\n")
hwrite(dwu_tab, p, table.style= 'border-collapse:collapse', table.cellpadding='5px', col.names=F)
}
#htmlSummary outputs summary information from a summaryTable output list to an HTML file
#the input sumset must have, at the very least:
# $summarystats: list of summary stats (path coverage, etc. .)
# $pathsummary: the generalized output table
# $genesummary: counts of times each gene is targeted across cohort
htmlSummary<-(sumset,
fname="./test.html",
plimit=.05,
sortColumn=5,
pagetitle=,
path_detail=){
("\nCreating HTML summary. . . \n")
("\nFile name of HTML summary:\n")
(fname)
#startHTMLPlugIns()
targetType=sumset$summarystats[2,2]
headings = ()
#check if there's a coverage summary
baseDir = ((fname, "/")[1]][1:(((fname, "/")[1]])-1)], collapse="/")
("coverage_summary"%in%(sumset)){
("Found coverage summary")
htmlSummary(path_detail=path_detail, sumset=sumset$coverage_summary, fname=(baseDir, "coverage_summary.html", sep="/"), path_detail)
(!"summarystats"%in%(sumset)){#if a coverage summary is all there is, then exit
()
}
("Back to writing",fname,"\n" )
}
(!(fname)){
("\nCreating file folder..\n")
(path=baseDir, recursive=T,showWarnings=F)
}
("names in sumset:\n")
((sumset))
((sumset$patientGeneMatrix))
(((sumset$summarystats))) () #if it's blank, get outta there
("just about to look at the image slots...")
multiPatients = F
(!(sumset$patientsums)){
multiPatients = ((sumset$patientsums))>1
}
("\nOpening HTML page..\n")
(fname)
p = openPage(fname)
(!(pagetitle)) hwrite((="_",replacement=" ",x=pagetitle), p, heading=1)
headings = (("Summary of",sumset$summarystats[1,2],"used in enrichment testing"),
("Distribution of",targetType,"genes across patients"),
("Top 20 most commonly",targetType,"genes, across cohort"),
"Analysis of pathways",
"Individual patient enrichments for the patients with the top 5 highest mutation rates ",
headings)
("headings:")
(headings)
(multiPatients)
(("imageSlots"%in%(sumset)&!(sumset$imageSlots$imageSlots))|multiPatients){
#establish image directory: any images will go here
imageDir = (baseDir, "/imageSlots/", sep="")
("\nCreating image directory..\n")
(imageDir,recursive=T,showWarnings=F)
headings = (headings, "Path images")
}
#for all the headings, turn the heading into an anchor
links = headings
hwrite("<A name=\"section0.1\">Table of contents</A>", p, heading=2)
("Data_work_up_notes"%in%(sumset)) hwrite(("<A href=\"#section",".001","\">","Data work-up notes", "</A>", sep=""),p,heading=3)
("settings"%in%(sumset)) hwrite(("<A href=\"#section",".002","\">","Settings", "</A>", sep=""),p,heading=3)
(i 1:(headings)){
links[i] = ("<A href=\"#section",i,"\">",headings[i], "</A>", sep="")
headings[i] = ("<A name=\"section",i,"\">",headings[i], "</A>", sep="")
hwrite(links[i],p,heading=3)
}
WriteDataWorkUp(sumset, p, baseDir)
WriteSettings(sumset, p=p)
("\nTarget type:",targetType,"\n")
hwrite(headings[1],p,heading=2)
hwrite("<A href=\"#section0.1\">Back to top</A><BR>", p, heading=4)
("just wrote heading 1")
#output the summary stats
sum_stats=sumset$summarystats
((sumset))
(sum_stats)
hwrite(sum_stats, p, table.style= 'border-collapse:collapse', table.cellpadding='5px')
("Just wrote sum_stats")
#output most targeted genes
(multiPatients){#if more than one patient's worth of gene data is being summarized here,
("\nOutputting distribution of",targetType,"genes across all patients\n")
#output information about the counts for each gene across the cohort
hwrite(headings[2],p,heading=2)
hwrite("<A href=\"#section0.1\">Back to top</A><BR>", p, heading=4)
#put a histogram into the summary
imageName = ("dist_",targetType,"_genes_per_patient.jpg",sep="")
imgloc = (baseDir,"/imageSlots/",imageName,sep="")
(imgloc)
(sumset$patientsums,1],breaks=20,main=("Distribution of",targetType,"gene counts for all patients:"),xlab=("Number of",targetType,"genes in patient"))
()
hwriteImage(("./imageSlots/", imageName, sep=""), p, br=)
genelist = sumset$genesummary
geneindex = (genelist,1],decreasing=T)
# ordgenelist = as.matrix(genelist[geneindex[ 1:min( 20,length(geneindex) ) ],],ncol=1)
ordgenelist = genelist[geneindex 1:( 20,(geneindex) ) ],,=F]
(!(sumset$preFilteringGeneSummary)){
("!is.null(sumset$preFilteringGeneSummary)")
genelistpre = sumset$preFilteringGeneSummary
geneindexpre = (genelistpre,2],decreasing=T)
ordgenelistpre = (genelistpre[geneindexpre[1:20],],=1)
ordgenelist = ((ordgenelist),ordgenelist,1])
ordgenelistpre = GenesInPaths(ordgenelistpre, path_detail)
ordgenelist = GenesInPaths(ordgenelist, path_detail)
ordgenelist = (ordgenelistpre, ordgenelist)
((ordgenelist))
(ordgenelist)<-("Top mutated genes, before filtering",
"Number of patients with mutated gene, pre filtering",
"In the current set of pathways?",
"Top mutated genes, after filtering",
"Number of patients with mutated gene, after filtering",
"In the current set of pathways?")
}{
ordgenelist = GenesInPaths(ordgenelist, path_detail)
}
("Bout to print heading 3")
hwrite(headings[3],p,heading=2)
hwrite("<A href=\"#section0.1\">Back to top</A><BR>", p, heading=4)
hwrite(ordgenelist, p, table.style= 'border-collapse:collapse', table.cellpadding='1px', =(sumset$preFilteringGeneSummary),col.names=T)
}
hwrite("<br><br>",p)
#pathway analysis/coverage/enrichment
hwrite(headings[4],p,heading=2)
hwrite("<A href=\"#section0.1\">Back to top</A><BR>", p, heading=4)
#output the significant paths
ptable = sumset$pathsummary
tmpTest = F
((ptable)){
tmpTest = T
} (!(ptable)|!(ptable)){
tmpTest=T
}
(tmpTest){
(sumset$summarystats[1,2])
# readline("empty table")
hwrite(("There are no",targetType,"genes in this patient.",p))
}{
ptable = cleanTables(ptable)
(ptable) = (="_",replacement=" ",x=(ptable))
((ptable)<6){
("Coverage summary.")
sortColumn = 4
}
hwrite(("Pathways sorted by",(ptable)[sortColumn],"."),p,heading=2)
torder = (ptable,sortColumn],decreasing=T)
hwrite(ptable[torder,], p, table.style= 'border-collapse:collapse', table.cellpadding='1px', =F)
}
#determine, for the patients with the most mutation, in the paths that are most significantly
#mutated, which genes are mutated the most
#find genes and pathways mutated in patients with highest overall mutational load
((sumset$path_summary_each_patient)){
hwrite(headings[5],p,heading=2)
hwrite("<A href=\"#section0.1\">Back to top</A><BR>", p, heading=4)
("About to write info about top 5 most affected patients.. . ")
((sumset$patientsums)>1){
#first: select to 5 patients
ordmurates = sumset$patientsums[order(sumset$patientsums,decreasing=T),]
((ordmurates))
top5pat = (ordmurates[1:5])#patients with the top 5 highest mutation rates
#second: select pathways most enriched paths in those patients
(ci 1:(top5pat)){
cp = top5pat[ci]
("\nPatient: ",cp,"\n")
#Error in if (nrow(sigpaths)) { : argument is of length zero
cpat = sumset$path_summary_each_patient[cp]]
(!(cpat$pathsummary)){
hwrite(("<br><br>In patient",cp,"(",ci,"of 5 patients with the largest number of",targetType,"genes):"),p,heading=2)
#third: find genes in those pathways which are significantly aberrational
sigpaths = cpat$pathsummary[cpat$pathsummary$hyperg_p_w_FDR<0.05,,=F]
((sigpaths)){
hwrite((targetType,"genes in pathways significantly enriched (hypergeometric, with BD FDR correction) for",targetType,"genes"),p,heading=2)
}{
hwrite(("No paths were found to be significantly enriched for",targetType,"genes as to the hypergeometric test with BH FDR"),p,heading=2)
sigpaths = cpat$pathsummary[cpat$pathsummary$hyperg_p_w_FDR<0.05,]
((sigpaths)){
hwrite((targetType,"genes in pathways significantly enriched (hypergeometric test) for",targetType,"genes"),p,heading=2)
}{
hwrite(("No paths were found to be significantly enriched for",targetType,"genes as to BH FDR"),p,heading=2)
hwrite((targetType,"genes in top 5 pathways most enriched (proportion enrichment) for ",targetType,"genes"),p,heading=2)
sigpaths = cpat$pathsummary[order(cpat$pathsummary,5],decreasing=T),][1:5,]
}
}
sigpaths2 = (sigpaths,cpat$active_genes_ea_path[rownames(sigpaths),])
(sigpaths2) = (="_",replacement=" ",x=(sigpaths2))
(sigpaths2)[ncol(sigpaths2)] = (targetType,"genes in this pathway")
hwrite(sigpaths2, p, table.style= 'border-collapse:collapse', table.cellpadding='1px', =F)
hwrite(((sigpaths2),"pathways total."),p,heading=2)
hwrite("<br>",p)
#find commonality of genes in top paths
tmp = (cpat$active_genes_ea_path[rownames(sigpaths),],=" ")
(tmp) = (sigpaths)
tmpmat = (tmp)
gcounts = ((T,(tmpmat))%*%tmpmat,=1)
(gcounts) = (tmpmat)
hwrite(("Number of the above",(tmpmat),"pathways each gene is found",targetType,"in."),p,heading=2)
(gcounts) = ("Number of the above pathways the gene is found in")
hwrite(gcounts, p, table.style= 'border-collapse:collapse', table.cellpadding='1px', =T)
hwrite(((tmpmat),"genes total."),p,heading=2)
}
}
}
}
("About to write image tags... ")
((sumset))
#this is defined above: imageDir = paste(baseDir, "/imageSlots/")
#search the summary for imageSlots
("imageSlots"%in%(sumset)){
((sumset$imageSlots)){
("looking at image slots")
(fname)
hwrite(x=headings[7], =p, heading=2)
hwrite("<A href=\"#section0.1\">Back to top</A><BR>", p, heading=4)
(sumset$imageSlots)
images = (sumset$imageSlots$imageSlots, =1)
(images)
((images)){
(i 1:(images)){
iname = (images)[i]
(iname)
ifile = images[i,1]
hwrite(iname,p,heading=2)
tmpImageLoc = ("./imageSlots/",ifile, sep="")
hwriteImage(tmpImageLoc, p, br=)
}
}
}{
hwrite(x="Sorry, no images were found", =p, heading=2)
}
}
closePage(p)
("\nSaved summary to HTML file: ",fname,"\n")
()
}# htmlSummary()
# htmlSummary(sumset=somatic_summary,
# fname="./quick_check_aml_somatic_nci.html")
#htmlSummary(sumset=somatic_summary, fname="./output/testSavesomaticsummary.html")
#toHTML outputs an HTML file
#the input table_list cat have:
# tables
#
# image file names
#
toHTML<-(table_list,
fname="./test.html",
plimit=.05,
maxrows=100,
reorderTables=,
pagetitle=,
limit_col=,
path_detail=){
("____________________________Saving HTML page to :",fname,"\n")
#startHTMLPlugIns()
fnameroot = ((fname, ="/")[1]][1:(((fname, ="/")[1]])-1)],sep="/",collapse="/")
(fnameroot)
(fname)
p = openPage(fname)
(!(pagetitle)) hwrite(pagetitle, p, heading=1)
sl1 = "<A href=\"#section"
sl2 = "\">"
sl3 = "</A><BR>"
hwrite("<A name=\"section0.1\">Table of contents</A>", p, heading=2)
#first, print the list of elements
hwrite("These are the sections of this HTML page:", p, heading=3)
(i 1:(table_list)){
curlink = (sl1,i,sl2, (table_list)[i], sl3, sep="")
hwrite(curlink, p, heading=3)
}
trg1 = "<A name=\"section"
trg2 = "\">"
trg3 = "</A>"
((table_list)){
("settings"%in%(table_list)){
#if it's from the overlap analysis
((table_list$settings)=="list"){
(( (table_list$settings)%in%("basic_overlap", "defaultSummaryTable", "combinedAberrations") ) ){
nsnames = ("setting_for_",(table_list$settings))
(table_list$settings)<-nsnames
table_list[nsnames] = table_list$settings
table_list$settings=
}
}
}
(i 1:(table_list)){
(!(table_list)) (table_list)
curname = (table_list)[i]
(VERBOSE) (curname)
((curname)) (table_list)
(VERBOSE) ("names(table_list)")
(VERBOSE) ((table_list))
(VERBOSE) (i)
(VERBOSE) ("stops here boo hoo")
postfix = (curname,="\\." )[1]] ((curname,="\\.")[1]]) ]
( (="setting", x=curname) & (table_list[i]])!="data.frame" ){
table_list[i]] = listToDf(lst=table_list[i]])
}
( postfix %in% ("png","jpeg","bmp")){#if it's an image
curnameanchor = (trg1, i, trg2, curname, trg3, sep="")
(VERBOSE) (curnameanchor)
hwrite(curnameanchor, p, heading=2)
hwrite("<A href=\"#section0.1\">Back to top</A><BR>", p, heading=4)
("\nImage found")
hwriteImage(curname, p, br=)
} (curname=="imageSlots"){
#add anchor
curnameAnchor = (trg1, i, trg2, "Network diagrams of pathways which are significantly enriched in aberration and contain drug targets", trg3, sep="")
(VERBOSE) (curnameAnchor)
hwrite(curnameAnchor, p, heading=2)
hwrite("<A href=\"#section0.1\">Back to top</A><BR>", p, heading=4)
("Image list found...\n")
("(file:",fname,")\n")
("\nCurrent item being output:",curname,"\n")
#if the curname is imageSlots,
(VERBOSE) ((table_list))
#get the matrix out and print the images
imageSet = table_list"imageSlots"]]"imageSlots"]]
#check if there are images in there
((imageSet)){
(VERBOSE) ("no images found in image folder")
hwrite("Sorry, images not available", p, heading=2)
}{
(VERBOSE) ("images found..")
(VERBOSE) (imageSet)
imageSetNames = (imageSet)
((imageSet)) imageSetNames = (imageSet)
(im 1:(imageSet)){
imageName = imageSetNames[im]
(VERBOSE) (imageName)
imageFile = ("./imageSlots/",imageSet[im], sep="")
(VERBOSE) (imageFile)
hwrite(imageName, p, heading=3)
hwriteImage(imageFile, p, br=)
}
}
} (("_each_patient",x=curname, ignore.case=T)){#it's the ... for each patient
("OLA eap")
eaPatSet = table_list[i]]
(patsetname (eaPatSet)){
("**********************************************************")
ePatRoot = (fnameroot,curname,patsetname,("overlap_analysis/"),sep="/")
SaveToHTML_inner(results=eaPatSet[patsetname]],
study_name=patsetname,
fileRoot=ePatRoot,
selNames=(eaPatSet[patsetname]]))
}
} (("_summary",x=curname,ignore.case=T)){
("Found a summary:::::::::::::::::::::::::::+:::::::::::::::::::::::::")
sumfname =(fnameroot,curname,(curname,"_.html",sep=""),sep="/")
ctab = table_list[curname]]
htmlSummary(path_detail=path_detail,
sumset=ctab,
fname=sumfname,
pagetitle=curname)
} ((table_list[i]])&!(table_list[i]])){
("is.list(table_list[[i]])&!is.data.frame(table_list[[i]])")
#if it's a list, then call toHTML again
listfname = (fnameroot,curname,sep="/")
toHTML(path_detail=path_detail,
table_list=table_list[i]],
limit_col=limit_col,
pagetitle=(pagetitle, curname,sep=":"),
reorderTables=reorderTables,
fname=listfname,
plimit=plimit,
maxrows=maxrows)
}{ # it's not an image, summary table, folder of images or 'analysis each patient' set
curname = (trg1, i, trg2, curname, trg3, sep="")
hwrite(curname, p, heading=2)
hwrite("<A href=\"#section0.1\">Back to top</A><BR>", p, heading=4)
initial_table = table_list[i]]#if it's turned into a vector clean table will not work; needs to be data frame
((initial_table)&!(initial_table)){
initial_table = (initial_table,=1)
}
((initial_table))
((initial_table)){
curtab = cleanTables(tab=initial_table)#makes sure there are only 3 sig figs instead of 9
((initial_table))
(!(limit_col)&(limit_col%in%(curtab))){#if there are columns to limit the output by
("There are columns by which to limit the number of table rows printed.\n")
( limit_col){#for each type of limit that is to be set, scan the table for the limit's column,
(%in%(curtab)){#if the column is found, limit and sort the table based on that column
curindex = (curtab,])<plimit
tmptab=curtab[curindex,]
curorder = 1:(tmptab)
(!(reorderTables)){
("reorderTables")
curorder = (tmptab,],decreasing=reorderTables)
}
(tmptab)<-(="_",replacement=" ",x=(tmptab))
hwrite((curname,"limited to ",,"<",(plimit),"(", ((tmptab)), "records exist, showing top", ((maxrows,(tmptab))) ," )"), p, heading=2)
hwrite(x=tmptab[curorder[1:(maxrows,(tmptab))],], =p,
=!((tmptab)==tmptab,1]),
table.style= 'border-collapse:collapse',
table.cellpadding='1px')
("header and table written")
}
}
}{ #if there are no columns to limit the output by
curtab = curtab[1:(maxrows,(curtab)),,=F]#limit the number of rows
(!(reorderTables)){#order the rows if wanted
curtab = curtab[order(curtab,5],decreasing=T),,=F]
}
= ("(", ((initial_table)), "records, showing top", ((curtab)) ,")")
hwrite(, p, heading=4)
(curtab)<-(="_",replacement=" ",x=(curtab))
hasRowNames = !((curtab)==curtab,1])
hwrite(curtab,=hasRowNames,
p,
table.style= 'border-collapse:collapse',
table.cellpadding='1px')
#readline("press enter to continue")
("~")
}
}{
hwrite(("Zero",curname,"were found."), p, heading=4)
}#if/else there are rows in curtab
}#if/else is image
}#for each element in the table_list
}
("\nHTML summary written to",fname,"\n")
closePage(p)
}
#'@title Save results to an HTML page.
#'@description Writes data from the \code{results} slot of a \code{Study} object to an HTML document.
#'@param study A \code{Study} object with results sets from one or more study arm.
#'@import xtable
#'@import hwriter
#'@export
<-(study){
# if(!require("hwriterPlus")){
# print("Trying to install hwriterPlus so that HTML output can be generated")
# install.packages("hwriterPlus")
# if(require("hwriterPlus")){
# print("hwriterPlus installed and loaded")
# } else {
# stop("could not install hwriterPlus")
# }
# }
# if(!require("xtable")){
# print("Trying to install hwriterPlus so that HTML output can be generated")
# install.packages("hwriterPlus")
# if(require("xtable")){
# print("xtable installed and loaded")
# } else {
# stop("could not install xtable")
# }
# }
SaveToHTML_inner(study_name=study@studyMetaData@studyName,
results=study@results,
path_detail=study@studyMetaData@paths)
}
SaveToHTML_inner<-(results,
study_name,
selNames=,
fileRoot=,
path_detail=){
("Saving summary to HTML..")
(!(="^test", ignore.case=T, x=study_name)){
study_name = ("study_", study_name)
}
((fileRoot)) fileRoot = ("./output/",study_name,"/results",sep="")
("file root:")
(fileRoot)
isels = (results)
badi = (="file", x=(results), ignore.case=T)
((badi)){
isels = (results)-badi]
}
#pick out all the indexes except the path_file_name
((selNames)){
#let the user decide
selvector.li = selectionList(valcol=isels)
#pullOutTheCorrectNames
selNames = isels[selvector.li]
}
(n selNames){
("selnames:")
(selNames)
(n)
#save em
ssum = results[n]]
((="summary", x=n, ignore.case=T)){
("chose summary")
root2 = (fileRoot,"/",n,"/",sep="")
sumfname = (root2,n,"_ViewableSummary.html",sep="")
# readline("enter to continue")
htmlSummary(path_detail=path_detail, sumset=ssum, fname=sumfname, pagetitle=n)
htmlPerPatient(results=results, summ=n, fileroot=root2, overlap=F, path_detail=path_detail)
} ((="overlap_analysis$", x=n, ignore.case=T)){
("overlap analysis selected")
root2 = (fileRoot,"/",n,"/",sep="")
(root2,recursive=T,showWarnings=F)
outHTMLname = (root2, "overlap_analysis.html",sep="")
(outHTMLname)
(n)
# readline("enter to continue")
toHTML(path_detail=path_detail,
table_list=results$overlap_analysis,
fname=outHTMLname,
maxrows=10000,
pagetitle=n)
}
}
}
htmlPerPatient<-(results, summ, fileroot, overlap=T, path_detail=){
("\nWriting an HTML summary for each patient...........\n")
fileroot = (fileroot,"path_summary_each_patient/",sep="")
psep = results[summ]]$path_summary_each_patient
(patient (psep)){
("\n***************************",patient,"*******************************\n")
#for each patient, establish a summary
perPatRoot = (fileroot,patient,"/",summ,"_ViewableSummary.html",sep="")
(summ)
ss = psep[patient]]
(overlap) ss = psep[patient]][summ]]
htmlSummary(path_detail=path_detail, sumset=ss, fname=perPatRoot)
}
}
test.twoHistOnePlot<-(){
prefilt = (muts = (x=100:200, =T, size=100), stringsAsFactors=F)
postfilt = (muts=(x=50:150, =T, size=100), stringsAsFactors=F)
(dataset1=prefilt,
dataset2=postfilt,
x_label="Number of mutations",
y_label="Number of patients",
legend_titles=("Before filtering", "After filtering"),
main_title="Distributions of mutations in patients before and after filtering")
}
#twoHistOnePlot
#
#'@title Puts two histograms on the same plot.
#'@description Puts two histograms on the same plot, shows each with a different color, and overlap with a third color.
#'@param dataset1 A numeric data set; \code{data.frame or matrix} with one column.
#'@param dataset2 A numeric data set; \code{data.frame or matrix} with one column.
#'@param frequency_not_density Depricated.
#'@param main_title The main title for the plot.
#'@param x_label The label to be put along the x axis.
#'@param y_label The label to be put along the y axis.
#'@param legend_titles Character vector, length 2, giving the names of the data in dataset1 and dataset2.
#'@param breaks Depricated.
#'@import ggplot2
<-(dataset1, dataset2,
frequency_not_density=T,
main_title="main_title not set",
x_label="x_label not set",
y_label="Frequency",
legend_titles=("dataset 1","dataset 2"),
breaks=(30,30)){
("\nReporting distributions of ", legend_titles[1],"and",legend_titles[2],"...")
#first combine the two data sets
(dataset1)<-"dat"
(dataset2)<-"dat"
stackedSet = (dataset1, dataset2)
idCol = ((legend_titles[1],times=(dataset1)), (legend_titles[2], times=(dataset2)))
stackedSet = (stackedSet, idCol)
(stackedSet)<-("value","idCol")
# pmin = min(c(dataset1, dataset2))
# pmax = max(c(dataset1, dataset2))
#
# p1 <- hist(dataset1, freq=frequency_not_density,breaks=breaks[1],)
# p2 <- hist(dataset2, freq=frequency_not_density,breaks=breaks[2])
#
# hist = plot( p2, col=rgb(0,0,1,1/4), xlim=c(pmin,pmax),
# freq=frequency_not_density,
# main=main_title,
# xlab=x_label,
# ylab = y_label) # first histogram
# hist2 = plot( p1, col=rgb(1,0,0,1/4), xlim=c(pmin,pmax),
# add=T,
# freq=frequency_not_density) # second
#
# if(!is.null(legend_titles)){
# legend("topright", inset=.05,
# legend_titles,
# fill=c(rgb(1,0,0,1/4),rgb(0,0,1,1/4)),
# horiz=F)
# }
bwidth = (test=((stackedSet$value))>10, yes=3, no=1)
p1 = ggplot(stackedSet, aes(x=value, fill=idCol)) +
geom_histogram(alpha=0.5, position="identity", binwidth=bwidth)+
ggtitle(main_title)+
theme_bw()+
theme(legend.title=element_blank())+
xlab(x_label)
# if(min(stackedSet[,1])<=0){
# rng = (range(stackedSet$value)[2]-range(stackedSet$value)[1])*.1
# p1 = p1+scale_x_continuous(limits=c( (min(stackedSet[,"value"])-rng), (max(stackedSet[,"value"])+rng ) ))
# }
(p1)
("...done\n")
(p1)
}
scatterhist <- (x, y, xlab="", ylab="", main=""){
zones=((2,0,1,3), =2, byrow=)
(zones, widths=(4/5,1/5), heights=(1/5,4/5))
xhist = (x, =)
yhist = (y, =)
top = ((xhist$counts, yhist$counts))
(mar=(3,3,1,1))
(x,y, main="")
(mar=(0,3,1,1))
(xhist$counts, axes=, =(0, top), space=0)
(mar=(3,0,1,1))
(yhist$counts, axes=, =(0, top), space=0, horiz=)
(oma=(3,3,0,0))
(=main,side=3,=T,padj=2)
(xlab, side=1, =1, =, adj=0,
at=.8 * ((x) - (x))/((x)-(x)))
(ylab, side=2, =1, =, adj=0,
at=(.8 * ((y) - (y))/((y) - (y))))
(mfrow = (1,1))
(mar=(5, 4, 4, 2) + 0.1 )
}
scatterhist2 <- (x1, y1, x2, y2, xlab="", ylab="", main="", legendTxt=){
zones=((2,0,1,3), =2, byrow=)
(zones, widths=(4/5,1/5), heights=(1/5,4/5))
xhist = (dataset1=x1, dataset2=x2,
legend_titles=,
main_title="",
y_label="",
x_label="")#hist(x, plot=FALSE)
yhist = (dataset1=y1, dataset2=y2,
legend_titles=,
main_title="",
y_label="",
x_label="")#hist(y, plot=FALSE)
top = ((xhist$p1$, yhist$p1$, xhist$p2$, yhist$p2$))
(mar=(3,3,1,1))
(x1,y1, main=main, =(0,0,1,1/4), =46)
(x2,y2,=(1,0,0,1/4), =46)
(!(legendTxt)){
("topright", (legendTxt[1], legendTxt[2]), fill=((1,0,0,1/4),(0,0,1,1/4)))
}
(mar=(0,3,1,1))
(xhist$p1$, axes=, =(0, top), space=0, =(0,0,1,1/4))
(xhist$p2$, axes=, =(0, top), space=0, add=T, =(1,0,0,1/4))
(mar=(3,0,1,1))
(yhist$p1$, axes=, =(0, top), space=0, horiz=, =(0,0,1,1/4))
(yhist$p2$, axes=, =(0, top), space=0, horiz=, add=T, =(1,0,0,1/4))
(oma=(3,3,0,0))
(xlab, side=1, =1, =, adj=0,
at=.8 * ((x1) - (x1))/((x1)-(x1)))
(ylab, side=2, =1, =, adj=0,
at=(.8 * ((y1) - (y1))/((y1) - (y1))))
(mfrow = (1,1))
(mar=(5, 4, 4, 2) + 0.1 )
}
addGeneInfo<-(genelist, path_detail){
#adds additional gene info to a list of genes
#indicates genes which did not have hugo records
#initilize variables
symtab = path_detail$symtable
# #assure there are no symbols given as "\\N"
# genelist = genelist[genelist!="\\N"]
#attempt to make sure the symbols match
genelist = (genelist, symref=symtab, verbose=T)
#extract the needed columns
xsymtab = symtab,("Approved.Symbol", "Approved.Name", "Locus.Type", "Chromosome", "Previous.Names")]
(xsymtab) = xsymtab$Approved.Symbol
#figure out which gene symbols we actually have info on:
have = genelist%in%xsymtab$Approved.Symbol
#reduce the gene list to those we have info on and report those we dont have info on
goodset = genelist[have]
badset = genelist!have]
out=()
out$'Genes and their records' = xsymtab[goodset,]
out$'Gene symbols for which records could not be found' = badset
(out)
}
