Mercurial > repos > fubar > rglasso
changeset 4:fb4959ed5b2b draft
Fixes to paths in git for deps
| author | fubar |
|---|---|
| date | Sat, 31 Oct 2015 02:26:24 -0400 |
| parents | 01a5391488cd |
| children | 2b7e28499fd8 |
| files | .shed.yml readme.rst rgToolFactory.py rg_nri.xml rglasso_cox.xml test-data/cox_coxlassotest_glmnet_cvdeviance.pdf test-data/cox_coxlassotest_glmnetdev.pdf test-data/cox_test.xls test-data/coxlassotest.html test-data/coxlassotest_modelres.xls test-data/genTest.R test-data/nri_test1.xls test-data/nri_test1_out.html test-data/nri_test1_out.xls tool_dependencies.xml |
| diffstat | 14 files changed, 3126 insertions(+), 0 deletions(-) [+] |
line wrap: on
line diff
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/.shed.yml Sat Oct 31 02:26:24 2015 -0400 @@ -0,0 +1,9 @@ +categories: +- Statistics +description: glmnet and lars for galaxy +long_description: | + Allows hybrid lasso-ridge regression using glmnet +name: rglasso_1_9_8 +owner: fubar +remote_repository_url: https://github.com/galaxyproject/tools-iuc/tree/master/tools/rglasso +type: unrestricted
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/readme.rst Sat Oct 31 02:26:24 2015 -0400 @@ -0,0 +1,39 @@ +glmnet wrappers +=============== + +This is a self installing Galaxy tool exposing the glmnet_ R package which has excellent documentation at +glmnet_ Minimal details are provided in this wrapper - please RTM to get the best out of it. + +The tool exposes the entire range of penalised maximum likelihood +GLM models ranging from pure lasso (set alpha to 1) to pure ridge-regression (set alpha to 0). + +These models can be k-fold internally cross validated to help select an "optimal" predictive or classification +algorithm. Predictive coefficients for each included independent variable are output for each model. + +Predictors can be forced into models to adjust for known confounders or explanatory factors. + +The glmnet_ implementation of the coordinate descent algorithm is fast and efficient even on relatively large problems +with tens of thousands of predictors and thousands of samples - such as normalised microarray intensities and anthropometry +on a very large sample of obese patients. + +The user supplies a tabular file with rows as samples and columns containing observations, then chooses +as many predictors as required. A separate model will be output for each of potentially multiple dependent +variables. Models are reported as the coefficients for terms in an 'optimal' model. +These optimal predictors are selected by repeatedly setting +aside a random subsample, building a model in the remainder and estimating AUC or deviance +using k (default 10) fold internal cross validation. For each of these steps, a random 1/k +of the samples are set aside and used to estiamte performance of an optimal model estimated +from the remaining samples. Plots are provided showing the range of these (eg 10) internal validation +estimates and mean model AUC (binomial) or residual deviance plots at each penalty increment step. + +A full range of link functions are available including Gaussian, Poisson, Binomial and +Cox proportional hazard time to failure for censored data in this wrapper. + +Note that multinomial and multiresponse gaussian models are NOT yet implemented since I have not yet +had use for them - send code! + +.. _glmnet: http://web.stanford.edu/~hastie/glmnet/glmnet_alpha.html + +Wrapper author: Ross Lazarus +19 october 2014 +
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/rgToolFactory.py Sat Oct 31 02:26:24 2015 -0400 @@ -0,0 +1,649 @@ +# rgToolFactory.py +# see https://bitbucket.org/fubar/galaxytoolfactory/wiki/Home +# +# copyright ross lazarus (ross stop lazarus at gmail stop com) May 2012 +# +# all rights reserved +# Licensed under the LGPL +# suggestions for improvement and bug fixes welcome at https://bitbucket.org/fubar/galaxytoolfactory/wiki/Home +# +# july 2014 +# added buffered read of sterror after run +# +# august 2013 +# found a problem with GS if $TMP or $TEMP missing - now inject /tmp and warn +# +# july 2013 +# added ability to combine images and individual log files into html output +# just make sure there's a log file foo.log and it will be output +# together with all images named like "foo_*.pdf +# otherwise old format for html +# +# January 2013 +# problem pointed out by Carlos Borroto +# added escaping for <>$ - thought I did that ages ago... +# +# August 11 2012 +# changed to use shell=False and cl as a sequence + +# This is a Galaxy tool factory for simple scripts in python, R or whatever ails ye. +# It also serves as the wrapper for the new tool. +# +# you paste and run your script +# Only works for simple scripts that read one input from the history. +# Optionally can write one new history dataset, +# and optionally collect any number of outputs into links on an autogenerated HTML page. + +# DO NOT install on a public or important site - please. + +# installed generated tools are fine if the script is safe. +# They just run normally and their user cannot do anything unusually insecure +# but please, practice safe toolshed. +# Read the fucking code before you install any tool +# especially this one + +# After you get the script working on some test data, you can +# optionally generate a toolshed compatible gzip file +# containing your script safely wrapped as an ordinary Galaxy script in your local toolshed for +# safe and largely automated installation in a production Galaxy. + +# If you opt for an HTML output, you get all the script outputs arranged +# as a single Html history item - all output files are linked, thumbnails for all the pdfs. +# Ugly but really inexpensive. +# +# Patches appreciated please. +# +# +# long route to June 2012 product +# Behold the awesome power of Galaxy and the toolshed with the tool factory to bind them +# derived from an integrated script model +# called rgBaseScriptWrapper.py +# Note to the unwary: +# This tool allows arbitrary scripting on your Galaxy as the Galaxy user +# There is nothing stopping a malicious user doing whatever they choose +# Extremely dangerous!! +# Totally insecure. So, trusted users only +# +# preferred model is a developer using their throw away workstation instance - ie a private site. +# no real risk. The universe_wsgi.ini admin_users string is checked - only admin users are permitted to run this tool. +# + +import sys +import shutil +import subprocess +import os +import time +import tempfile +import optparse +import tarfile +import re +import shutil +import math + +progname = os.path.split(sys.argv[0])[1] +myversion = 'V000.2 June 2012' +verbose = False +debug = False +toolFactoryURL = 'https://bitbucket.org/fubar/galaxytoolfactory' +buffsize = 1048576 + + +def timenow(): + """return current time as a string + """ + return time.strftime('%d/%m/%Y %H:%M:%S', time.localtime(time.time())) + +html_escape_table = { + "&": "&", + ">": ">", + "<": "<", + "$": "\$" + } + +def html_escape(text): + """Produce entities within text.""" + return "".join(html_escape_table.get(c,c) for c in text) + +def cmd_exists(cmd): + return subprocess.call("type " + cmd, shell=True, + stdout=subprocess.PIPE, stderr=subprocess.PIPE) == 0 + + +class ScriptRunner: + """class is a wrapper for an arbitrary script + """ + + def __init__(self,opts=None): + """ + cleanup inputs, setup some outputs + + """ + self.useGM = cmd_exists('gm') + self.useIM = cmd_exists('convert') + self.useGS = cmd_exists('gs') + self.temp_warned = False # we want only one warning if $TMP not set + if opts.output_dir: # simplify for the tool tarball + os.chdir(opts.output_dir) + self.thumbformat = 'png' + self.opts = opts + self.toolname = re.sub('[^a-zA-Z0-9_]+', '', opts.tool_name) # a sanitizer now does this but.. + self.toolid = self.toolname + self.myname = sys.argv[0] # get our name because we write ourselves out as a tool later + self.pyfile = self.myname # crude but efficient - the cruft won't hurt much + self.xmlfile = '%s.xml' % self.toolname + s = open(self.opts.script_path,'r').readlines() + s = [x.rstrip() for x in s] # remove pesky dos line endings if needed + self.script = '\n'.join(s) + fhandle,self.sfile = tempfile.mkstemp(prefix=self.toolname,suffix=".%s" % (opts.interpreter)) + tscript = open(self.sfile,'w') # use self.sfile as script source for Popen + tscript.write(self.script) + tscript.close() + self.indentedScript = '\n'.join([' %s' % x for x in s]) # for restructured text in help + self.escapedScript = '\n'.join([html_escape(x) for x in s]) + self.elog = os.path.join(self.opts.output_dir,"%s_error.log" % self.toolname) + self.tlog = os.path.join(self.opts.output_dir,"%s_runner.log" % self.toolname) + if opts.output_dir: # may not want these complexities + art = '%s.%s' % (self.toolname,opts.interpreter) + artpath = os.path.join(self.opts.output_dir,art) # need full path + artifact = open(artpath,'w') # use self.sfile as script source for Popen + artifact.write(self.script) + artifact.close() + self.html = [] + self.cl = (opts.interpreter,self.sfile) + self.outFormats = 'tabular' # TODO make this an option at tool generation time + self.inputFormats = 'tabular' # TODO make this an option at tool generation time + self.test1Input = '%s_test1_input.xls' % self.toolname + self.test1Output = '%s_test1_output.xls' % self.toolname + self.test1HTML = '%s_test1_output.html' % self.toolname + + def makeXML(self): + """ + Create a Galaxy xml tool wrapper for the new script as a string to write out + fixme - use templating or something less fugly than this example of what we produce + + <tool id="reverse" name="reverse" version="0.01"> + <description>a tabular file</description> + <command interpreter="python"> + reverse.py --script_path "$runMe" --interpreter "python" + --tool_name "reverse" --input_tab "$input1" --output_tab "$tab_file" + </command> + <inputs> + <param name="input1" type="data" format="tabular" label="Select a suitable input file from your history"/><param name="job_name" type="text" label="Supply a name for the outputs to remind you what they contain" value="reverse"/> + + </inputs> + <outputs> + <data format="tabular" name="tab_file" label="${job_name}"/> + + </outputs> + <help> + +**What it Does** + +Reverse the columns in a tabular file + + </help> + <configfiles> + <configfile name="runMe"> + +# reverse order of columns in a tabular file +import sys +inp = sys.argv[1] +outp = sys.argv[2] +i = open(inp,'r') +o = open(outp,'w') +for row in i: + rs = row.rstrip().split('\t') + rs.reverse() + o.write('\t'.join(rs)) + o.write('\n') +i.close() +o.close() + + + </configfile> + </configfiles> + </tool> + + """ + newXML="""<tool id="%(toolid)s" name="%(toolname)s" version="%(tool_version)s"> + %(tooldesc)s + %(command)s + <inputs> + %(inputs)s + </inputs> + <outputs> + %(outputs)s + </outputs> + <configfiles> + <configfile name="runMe"> + %(script)s + </configfile> + </configfiles> + %(tooltests)s + <help> + %(help)s + </help> + </tool>""" # needs a dict with toolname, toolid, interpreter, scriptname, command, inputs as a multi line string ready to write, outputs ditto, help ditto + + newCommand="""<command interpreter="python"> + %(toolname)s.py --script_path "$runMe" --interpreter "%(interpreter)s" + --tool_name "%(toolname)s" %(command_inputs)s %(command_outputs)s + </command>""" # may NOT be an input or htmlout + tooltestsTabOnly = """<tests><test> + <param name="input1" value="%(test1Input)s" ftype="tabular"/> + <param name="job_name" value="test1"/> + <param name="runMe" value="$runMe"/> + <output name="tab_file" file="%(test1Output)s" ftype="tabular"/> + </test></tests>""" + tooltestsHTMLOnly = """<tests><test> + <param name="input1" value="%(test1Input)s" ftype="tabular"/> + <param name="job_name" value="test1"/> + <param name="runMe" value="$runMe"/> + <output name="html_file" file="%(test1HTML)s" ftype="html" lines_diff="5"/> + </test></tests>""" + tooltestsBoth = """<tests><test> + <param name="input1" value="%(test1Input)s" ftype="tabular"/> + <param name="job_name" value="test1"/> + <param name="runMe" value="$runMe"/> + <output name="tab_file" file="%(test1Output)s" ftype="tabular" /> + <output name="html_file" file="%(test1HTML)s" ftype="html" lines_diff="10"/> + </test></tests>""" + xdict = {} + xdict['tool_version'] = self.opts.tool_version + xdict['test1Input'] = self.test1Input + xdict['test1HTML'] = self.test1HTML + xdict['test1Output'] = self.test1Output + if self.opts.make_HTML and self.opts.output_tab <> 'None': + xdict['tooltests'] = tooltestsBoth % xdict + elif self.opts.make_HTML: + xdict['tooltests'] = tooltestsHTMLOnly % xdict + else: + xdict['tooltests'] = tooltestsTabOnly % xdict + xdict['script'] = self.escapedScript + # configfile is least painful way to embed script to avoid external dependencies + # but requires escaping of <, > and $ to avoid Mako parsing + if self.opts.help_text: + xdict['help'] = open(self.opts.help_text,'r').read() + else: + xdict['help'] = 'Please ask the tool author for help as none was supplied at tool generation' + coda = ['**Script**','Pressing execute will run the following code over your input file and generate some outputs in your history::'] + coda.append(self.indentedScript) + coda.append('**Attribution** This Galaxy tool was created by %s at %s\nusing the Galaxy Tool Factory.' % (self.opts.user_email,timenow())) + coda.append('See %s for details of that project' % (toolFactoryURL)) + coda.append('Please cite: Creating re-usable tools from scripts: The Galaxy Tool Factory. Ross Lazarus; Antony Kaspi; Mark Ziemann; The Galaxy Team. ') + coda.append('Bioinformatics 2012; doi: 10.1093/bioinformatics/bts573') + xdict['help'] = '%s\n%s' % (xdict['help'],'\n'.join(coda)) + if self.opts.tool_desc: + xdict['tooldesc'] = '<description>%s</description>' % self.opts.tool_desc + else: + xdict['tooldesc'] = '' + xdict['command_outputs'] = '' + xdict['outputs'] = '' + if self.opts.input_tab <> 'None': + xdict['command_inputs'] = '--input_tab "$input1" ' # the space may matter a lot if we append something + xdict['inputs'] = '<param name="input1" type="data" format="%s" label="Select a suitable input file from your history"/> \n' % self.inputFormats + else: + xdict['command_inputs'] = '' # assume no input - eg a random data generator + xdict['inputs'] = '' + xdict['inputs'] += '<param name="job_name" type="text" label="Supply a name for the outputs to remind you what they contain" value="%s"/> \n' % self.toolname + xdict['toolname'] = self.toolname + xdict['toolid'] = self.toolid + xdict['interpreter'] = self.opts.interpreter + xdict['scriptname'] = self.sfile + if self.opts.make_HTML: + xdict['command_outputs'] += ' --output_dir "$html_file.files_path" --output_html "$html_file" --make_HTML "yes" ' + xdict['outputs'] += ' <data format="html" name="html_file" label="${job_name}.html"/>\n' + if self.opts.output_tab <> 'None': + xdict['command_outputs'] += ' --output_tab "$tab_file"' + xdict['outputs'] += ' <data format="%s" name="tab_file" label="${job_name}"/>\n' % self.outFormats + xdict['command'] = newCommand % xdict + xmls = newXML % xdict + xf = open(self.xmlfile,'w') + xf.write(xmls) + xf.write('\n') + xf.close() + # ready for the tarball + + + def makeTooltar(self): + """ + a tool is a gz tarball with eg + /toolname/tool.xml /toolname/tool.py /toolname/test-data/test1_in.foo ... + """ + retval = self.run() + if retval: + print >> sys.stderr,'## Run failed. Cannot build yet. Please fix and retry' + sys.exit(1) + self.makeXML() + tdir = self.toolname + os.mkdir(tdir) + if self.opts.input_tab <> 'None': # no reproducible test otherwise? TODO: maybe.. + testdir = os.path.join(tdir,'test-data') + os.mkdir(testdir) # make tests directory + shutil.copyfile(self.opts.input_tab,os.path.join(testdir,self.test1Input)) + if self.opts.output_tab <> 'None': + shutil.copyfile(self.opts.output_tab,os.path.join(testdir,self.test1Output)) + if self.opts.make_HTML: + shutil.copyfile(self.opts.output_html,os.path.join(testdir,self.test1HTML)) + if self.opts.output_dir: + shutil.copyfile(self.tlog,os.path.join(testdir,'test1_out.log')) + op = '%s.py' % self.toolname # new name + outpiname = os.path.join(tdir,op) # path for the tool tarball + pyin = os.path.basename(self.pyfile) # our name - we rewrite ourselves (TM) + notes = ['# %s - a self annotated version of %s generated by running %s\n' % (op,pyin,pyin),] + notes.append('# to make a new Galaxy tool called %s\n' % self.toolname) + notes.append('# User %s at %s\n' % (self.opts.user_email,timenow())) + pi = open(self.pyfile,'r').readlines() # our code becomes new tool wrapper (!) - first Galaxy worm + notes += pi + outpi = open(outpiname,'w') + outpi.write(''.join(notes)) + outpi.write('\n') + outpi.close() + stname = os.path.join(tdir,self.sfile) + if not os.path.exists(stname): + shutil.copyfile(self.sfile, stname) + xtname = os.path.join(tdir,self.xmlfile) + if not os.path.exists(xtname): + shutil.copyfile(self.xmlfile,xtname) + tarpath = "%s.gz" % self.toolname + tar = tarfile.open(tarpath, "w:gz") + tar.add(tdir,arcname=self.toolname) + tar.close() + shutil.copyfile(tarpath,self.opts.new_tool) + shutil.rmtree(tdir) + ## TODO: replace with optional direct upload to local toolshed? + return retval + + + def compressPDF(self,inpdf=None,thumbformat='png'): + """need absolute path to pdf + note that GS gets confoozled if no $TMP or $TEMP + so we set it + """ + assert os.path.isfile(inpdf), "## Input %s supplied to %s compressPDF not found" % (inpdf,self.myName) + our_env = os.environ.copy() + if not (our_env.get('TMP',None) or our_env.get('TEMP',None)): + our_env['TMP'] = '/tmp' + if not self.temp_warned: + print >> sys.stdout,'## WARNING - no $TMP or $TEMP!!! Please fix - using /tmp temporarily' + self.temp_warned = True + hlog = os.path.join(self.opts.output_dir,"compress_%s.txt" % os.path.basename(inpdf)) + sto = open(hlog,'w') + outpdf = '%s_compressed' % inpdf + cl = ["gs", "-sDEVICE=pdfwrite", "-dNOPAUSE", "-dUseCIEColor", "-dBATCH","-dPDFSETTINGS=/printer", "-sOutputFile=%s" % outpdf,inpdf] + x = subprocess.Popen(cl,stdout=sto,stderr=sto,cwd=self.opts.output_dir,env=our_env) + retval1 = x.wait() + sto.close() + if retval1 == 0: + os.unlink(inpdf) + shutil.move(outpdf,inpdf) + os.unlink(hlog) + else: + x = open(hlog,'r').readlines() + print >> sys.stdout,x + hlog = os.path.join(self.opts.output_dir,"thumbnail_%s.txt" % os.path.basename(inpdf)) + sto = open(hlog,'w') + outpng = '%s.%s' % (os.path.splitext(inpdf)[0],thumbformat) + if self.useGM: + cl2 = ['gm', 'convert', inpdf, outpng] + else: # assume imagemagick + cl2 = ['convert', inpdf, outpng] + x = subprocess.Popen(cl2,stdout=sto,stderr=sto,cwd=self.opts.output_dir,env=our_env) + retval2 = x.wait() + sto.close() + if retval2 <> 0: + x = open(hlog,'r').readlines() + print >> sys.stdout,x + else: + os.unlink(hlog) + retval = retval1 or retval2 + return retval + + + def getfSize(self,fpath,outpath): + """ + format a nice file size string + """ + size = '' + fp = os.path.join(outpath,fpath) + if os.path.isfile(fp): + size = '0 B' + n = float(os.path.getsize(fp)) + if n > 2**20: + size = '%1.1f MB' % (n/2**20) + elif n > 2**10: + size = '%1.1f KB' % (n/2**10) + elif n > 0: + size = '%d B' % (int(n)) + return size + + def makeHtml(self): + """ Create an HTML file content to list all the artifacts found in the output_dir + """ + + galhtmlprefix = """<!DOCTYPE html PUBLIC "-//W3C//DTD XHTML 1.0 Transitional//EN" "http://www.w3.org/TR/xhtml1/DTD/xhtml1-transitional.dtd"> + <html xmlns="http://www.w3.org/1999/xhtml" xml:lang="en" lang="en"> + <head> <meta http-equiv="Content-Type" content="text/html; charset=utf-8" /> + <meta name="generator" content="Galaxy %s tool output - see http://g2.trac.bx.psu.edu/" /> + <title></title> + <link rel="stylesheet" href="/static/style/base.css" type="text/css" /> + </head> + <body> + <div class="toolFormBody"> + """ + galhtmlattr = """<hr/><div class="infomessage">This tool (%s) was generated by the <a href="https://bitbucket.org/fubar/galaxytoolfactory/overview">Galaxy Tool Factory</a></div><br/>""" + galhtmlpostfix = """</div></body></html>\n""" + + flist = os.listdir(self.opts.output_dir) + flist = [x for x in flist if x <> 'Rplots.pdf'] + flist.sort() + html = [] + html.append(galhtmlprefix % progname) + html.append('<div class="infomessage">Galaxy Tool "%s" run at %s</div><br/>' % (self.toolname,timenow())) + fhtml = [] + if len(flist) > 0: + logfiles = [x for x in flist if x.lower().endswith('.log')] # log file names determine sections + logfiles.sort() + logfiles = [x for x in logfiles if os.path.abspath(x) <> os.path.abspath(self.tlog)] + logfiles.append(os.path.abspath(self.tlog)) # make it the last one + pdflist = [] + npdf = len([x for x in flist if os.path.splitext(x)[-1].lower() == '.pdf']) + for rownum,fname in enumerate(flist): + dname,e = os.path.splitext(fname) + sfsize = self.getfSize(fname,self.opts.output_dir) + if e.lower() == '.pdf' : # compress and make a thumbnail + thumb = '%s.%s' % (dname,self.thumbformat) + pdff = os.path.join(self.opts.output_dir,fname) + retval = self.compressPDF(inpdf=pdff,thumbformat=self.thumbformat) + if retval == 0: + pdflist.append((fname,thumb)) + else: + pdflist.append((fname,fname)) + if (rownum+1) % 2 == 0: + fhtml.append('<tr class="odd_row"><td><a href="%s">%s</a></td><td>%s</td></tr>' % (fname,fname,sfsize)) + else: + fhtml.append('<tr><td><a href="%s">%s</a></td><td>%s</td></tr>' % (fname,fname,sfsize)) + for logfname in logfiles: # expect at least tlog - if more + if os.path.abspath(logfname) == os.path.abspath(self.tlog): # handled later + sectionname = 'All tool run' + if (len(logfiles) > 1): + sectionname = 'Other' + ourpdfs = pdflist + else: + realname = os.path.basename(logfname) + sectionname = os.path.splitext(realname)[0].split('_')[0] # break in case _ added to log + ourpdfs = [x for x in pdflist if os.path.basename(x[0]).split('_')[0] == sectionname] + pdflist = [x for x in pdflist if os.path.basename(x[0]).split('_')[0] <> sectionname] # remove + nacross = 1 + npdf = len(ourpdfs) + + if npdf > 0: + nacross = math.sqrt(npdf) ## int(round(math.log(npdf,2))) + if int(nacross)**2 != npdf: + nacross += 1 + nacross = int(nacross) + width = min(400,int(1200/nacross)) + html.append('<div class="toolFormTitle">%s images and outputs</div>' % sectionname) + html.append('(Click on a thumbnail image to download the corresponding original PDF image)<br/>') + ntogo = nacross # counter for table row padding with empty cells + html.append('<div><table class="simple" cellpadding="2" cellspacing="2">\n<tr>') + for i,paths in enumerate(ourpdfs): + fname,thumb = paths + s= """<td><a href="%s"><img src="%s" title="Click to download a PDF of %s" hspace="5" width="%d" + alt="Image called %s"/></a></td>\n""" % (fname,thumb,fname,width,fname) + if ((i+1) % nacross == 0): + s += '</tr>\n' + ntogo = 0 + if i < (npdf - 1): # more to come + s += '<tr>' + ntogo = nacross + else: + ntogo -= 1 + html.append(s) + if html[-1].strip().endswith('</tr>'): + html.append('</table></div>\n') + else: + if ntogo > 0: # pad + html.append('<td> </td>'*ntogo) + html.append('</tr></table></div>\n') + logt = open(logfname,'r').readlines() + logtext = [x for x in logt if x.strip() > ''] + html.append('<div class="toolFormTitle">%s log output</div>' % sectionname) + if len(logtext) > 1: + html.append('\n<pre>\n') + html += logtext + html.append('\n</pre>\n') + else: + html.append('%s is empty<br/>' % logfname) + if len(fhtml) > 0: + fhtml.insert(0,'<div><table class="colored" cellpadding="3" cellspacing="3"><tr><th>Output File Name (click to view)</th><th>Size</th></tr>\n') + fhtml.append('</table></div><br/>') + html.append('<div class="toolFormTitle">All output files available for downloading</div>\n') + html += fhtml # add all non-pdf files to the end of the display + else: + html.append('<div class="warningmessagelarge">### Error - %s returned no files - please confirm that parameters are sane</div>' % self.opts.interpreter) + html.append(galhtmlpostfix) + htmlf = file(self.opts.output_html,'w') + htmlf.write('\n'.join(html)) + htmlf.write('\n') + htmlf.close() + self.html = html + + def run(self): + """ + scripts must be small enough not to fill the pipe! + """ + my_env = os.environ.copy() + if self.opts.output_dir: + ste = open(self.elog,'w') + sto = open(self.tlog,'w') + sto.write('## Toolfactory running %s as %s script\n' % (self.toolname,self.opts.interpreter)) + sto.flush() + p = subprocess.Popen(self.cl,shell=False,stdout=sto,stderr=ste,cwd=self.opts.output_dir,env=my_env) + retval = p.wait() + sto.close() + ste.close() + # get stderr, allowing for case where it's very large + tmp_stderr = open( self.elog, 'rb' ) + stderr = '' + try: + while True: + stderr += tmp_stderr.read( buffsize ) + if not stderr or len( stderr ) % buffsize != 0: + break + except OverflowError: + pass + tmp_stderr.close() + else: + p = subprocess.Popen(self.cl,shell=False,env=my_env) + retval = p.wait() + if self.opts.make_HTML: + self.makeHtml() + return retval + + + def remove_me_runBash(self): + """ + cannot use - for bash so use self.sfile + """ + if self.opts.output_dir: + s = '## Toolfactory generated command line = %s\n' % ' '.join(self.cl) + ste = open(self.elog,'w') + sto = open(self.tlog,'w') + sto.write(s) + sto.flush() + p = subprocess.Popen(self.cl,shell=False,stdout=sto,stderr=ste,cwd=self.opts.output_dir) + else: + p = subprocess.Popen(self.cl,shell=False) + retval = p.wait() + if self.opts.output_dir: + sto.close() + ste.close() + # get stderr, allowing for case where it's very large + tmp_stderr = open(self.elog, 'rb' ) + stderr = '' + try: + while True: + stderr += tmp_stderr.read( buffsize ) + if not stderr or len( stderr ) % buffsize != 0: + break + except OverflowError: + pass + tmp_stderr.close() + return retval + + +def main(): + u = """ + This is a Galaxy wrapper. It expects to be called by a special purpose tool.xml as (eg): + <command interpreter="python">rgToolFactory.py --script_path "$scriptPath" --tool_name "foo" --interpreter "Rscript" + </command> + The tool writes a script to a scriptPath using a configfile. + Some things in the script are templates. + The idea here is that depending on how this code is called, it uses the specified interpreter + to run a (hopefully correct) script/template. Optionally generates a clone of itself + which will run that same script/template as a toolshed repository tarball for uploading to a toolshed. + There's now an updated version which allows arbitrary parameters. + And so it goes. + """ + op = optparse.OptionParser() + a = op.add_option + a('--script_path',default=None) + a('--tool_name',default=None) + a('--interpreter',default=None) + a('--output_dir',default=None) + a('--output_html',default=None) + a('--input_tab',default="None") + a('--output_tab',default="None") + a('--user_email',default='Unknown') + a('--bad_user',default=None) + a('--make_Tool',default=None) + a('--make_HTML',default=None) + a('--help_text',default=None) + a('--tool_desc',default=None) + a('--new_tool',default=None) + a('--tool_version',default=None) + opts, args = op.parse_args() + assert not opts.bad_user,'UNAUTHORISED: %s is NOT authorized to use this tool until Galaxy admin adds %s to admin_users in universe_wsgi.ini' % (opts.bad_user,opts.bad_user) + assert opts.tool_name,'## Tool Factory expects a tool name - eg --tool_name=DESeq' + assert opts.interpreter,'## Tool Factory wrapper expects an interpreter - eg --interpreter=Rscript' + assert os.path.isfile(opts.script_path),'## Tool Factory wrapper expects a script path - eg --script_path=foo.R' + if opts.output_dir: + try: + os.makedirs(opts.output_dir) + except: + pass + r = ScriptRunner(opts) + if opts.make_Tool: + retcode = r.makeTooltar() + else: + retcode = r.run() + os.unlink(r.sfile) + if retcode: + sys.exit(retcode) # indicate failure to job runner + + +if __name__ == "__main__": + main() + +
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/rg_nri.xml Sat Oct 31 02:26:24 2015 -0400 @@ -0,0 +1,633 @@ +<tool id="rg_nri" name="NRI" version="0.03"> + <description>and other model improvement measures</description> + <requirements> + <requirement type="package" version="3.2.2">R_3_2_2</requirement> + <requirement type="package" version="1.3.18">graphicsmagick</requirement> + <requirement type="package" version="9.10">ghostscript</requirement> + <requirement type="package" version="3.2">glmnet_lars_3_2</requirement> + </requirements> + <command interpreter="python"> + rgToolFactory.py --script_path "$runme" --interpreter "Rscript" --tool_name "rg_NRI" + --output_dir "$html_file.files_path" --output_html "$html_file" --make_HTML "yes" + </command> +<configfiles> +<configfile name="runme"> + +<![CDATA[ + +### http://www.researchgate.net/publication/264672640_R_function_for_Risk_Assessment_Plot__reclassification_metrics_NRI_IDI_cfNRI code +### http://cjasn.asnjournals.org/content/early/2012/05/24/CJN.09590911.full is the reference +### lots of little tweaks and but fixes. Using t as a variable name seems fraught to me. +### Ross Lazarus october 2014 for a Galaxy tool wrapper using the toolfactory infrastucture + +############################################################################# +###Functions to create risk assessment plots and associated summary statistics +############################################################################# +### +### (c) 2012 Dr John W Pickering, john.pickering@otago.ac.nz, and Dr David Cairns +### Last modified August 2014 +### +### Redistribution and use in source and binary forms, with or without +### modification, are permitted provided that the following conditions are met: +### * Redistributions of source code must retain the above copyright +### notice, this list of conditions and the following disclaimer. +### * Redistributions in binary form must reproduce the above copyright +### notice, this list of conditions and the following disclaimer in +### the documentation and/or other materials provided with the distribution + +### FUNCTIONS +### raplot +### Produces a Risk Assessment Plot and outputs the coordinates of the four curves +### Based on: Pickering, J. W. and Endre, Z. H. (2012). New Metrics for Assessing Diagnostic Potential of +### Candidate Biomarkers. Clinical Journal of the American Society of Nephrology, 7, 1355–1364. doi:10.2215/CJN.09590911 +### +### statistics.raplot +### Produces the NRIs, IDIs, IS, IP, AUCs. +### Based on: Pencina, M. J., D'Agostino, R. B. and Steyerberg, E. W. (2011). +### Extensions of net reclassification improvement calculations to measure usefulness of new biomarkers. Statistics in Medicine, 30(1), 11–21. doi:10.1002/sim.4085 +### Pencina, M. J., D'Agostino, R. B. and Vasan, R. S. (2008). Evaluating the added predictive ability of a new marker: From area under the ROC curve to reclassification and beyond. +### Statistics in Medicine, 27(2), 157–172. doi:10.1002/sim.2929 +### DeLong, E., DeLong, D. and Clarke-Pearson, D. (1988). Comparing the areas under 2 or more correlated receiver operating characteristic curves - a nonparametric approach. +### Biometrics, 44(3), 837–845. +### +### summary.raplot +### Produces the NRIs, IDIs, IS, IP, AUCs with confidence intervals using a bootstrap or asymptotic procedure. (I prefer bootstrap which is chosed by cis=c("boot")) + + +### Required arguments for all functions: +### x1 is calculated risk (eg from a glm) for the null model, i.e. predict(,type="response") on a glm object +### x2 is calculated risk (eg from a glm) for the alternative model +### y is the case-control indicator (0 for controls, 1 for cases) +### Optional argument +### t are the boundaries of the risks for each group (ie 0, 1 and the thresholds beteween. eg c(0,0,3,0,7,1)). If missing, defaults to c(0, the incidence, 1) + + +### risk assessment plot + +library('e1071') +library('caret') +library('pROC') +library('Hmisc') +library('pracma') + +raplot = function(x1, x2, y, outplot, title) { + + roc.model1 = roc(y, x1) + roc.model2 = roc(y, x2) + sens.model1 = roc.model1\$sensitivities + spec.model1 = 1 - roc.model1\$specificities + n.model1 = length(sens.model1) + thresh.model1 = roc.model1\$thresholds + thresh.model1 = thresh.model1[c(-1,-n.model1)] + sens.model1 = sens.model1[c(-1,-n.model1)] + spec.model1 = spec.model1[c(-1,-n.model1)] + sens.model2 = roc.model2\$sensitivities + spec.model2 = 1 - roc.model2\$specificities + n.model2 = length(sens.model2) + thresh.model2 = roc.model2\$thresholds + thresh.model2[1]=0 + thresh.model2[length(thresh.model2)]=1 + thresh.model2 = thresh.model2[c(-1,-n.model2)] + sens.model2 = sens.model2[c(-1,-n.model2)] + spec.model2 = spec.model2[c(-1,-n.model2)] + + n.model1 = length(sens.model1) + n.model2 = length(sens.model2) + + ### actual plotting + pdf(outplot) + plot(thresh.model1, sens.model1, xlim = c(0, 1), ylim = c(0, 1), type = "n", + lty = 2, lwd = 2, xlab = "Risk of Event", ylab = "", col = "black", main=title) + grid() + + polygon(x = c(thresh.model1, thresh.model2[n.model2:1]), + y = c(sens.model1, sens.model2[n.model2:1]), border = NA, col = gray(0.8)) + polygon(x = c(thresh.model1, thresh.model2[n.model2:1]), + y = c(spec.model1, spec.model2[n.model2:1]), border = NA, col = gray(0.8)) + + lines(thresh.model1, sens.model1, type = "l", lty = 2, lwd = 2, col = "black") + lines(thresh.model2, sens.model2, type = "l", lty = 1, lwd = 2, col = "black") + + lines(thresh.model1, spec.model1, type = "l", lty = 2, lwd = 2, col = "red") + lines(thresh.model2, spec.model2, type = "l", lty = 1, lwd = 2, col = "red") + + text(x = -0.15, y = 0.4, labels = "Sensitivity, ", col = "black", xpd = TRUE, srt = 90) + text(x = -0.15, y = 0.4 + 0.175, labels = "1-Specificity", col = "red", xpd = TRUE, srt = 90) + legend("topleft", c("Event: New model", "Event: Baseline model", + "No Event: New model", "No Event: Baseline model"), + col = c("black", "black", "red", "red"), + lty = c(1,2, 1, 2), lwd = 2, bty = "n") + dev.off() + return(data.frame("Null.p.sens"=thresh.model1, + "Null.sens"=sens.model1, + "Null.p.1spec"=thresh.model1, + "Null.1spec"=sens.model1, + "Alt.p.sens"=thresh.model2, + "Alt.sens"=sens.model2, + "Alt.p.1spec"=thresh.model2, + "Alt.1spec"=sens.model2)) + +} + + + +### statistics from a raplot (is an adaptation of improveProb() from Hmisc) + +statistics.raplot = function(x1, x2, y, threshvec) +{ + + s = is.na(x1 + x2 + y) ###Remove rows with missing data + if (any(s)) { + smiss = sum(s) + s = !s + x1 = x1[s] + x2 = x2[s] + y = y[s] + print.noquote(paste('Warning: removed',smiss,'cases with missing values')) + } + n = length(y) + y = as.numeric(y) + u = sort(unique(y)) + if (length(u) != 2 || u[1] != 0 || u[2] != 1) { + print.noquote("INPUT ERROR: y must have only two values: 0 and 1") + sink() + quit(save="no",status=2) + } + r = range(x1, x2) + if (r[1] < 0 || r[2] > 1) { + print.noquote("INPUT ERROR: x1 and x2 must be in [0,1]") + sink() + quit(save="no",status=3) + } + incidence=sum(y)/n + if (missing(threshvec)) { + threshvec=c(0, incidence,1) + print(paste('threshvec missing. using',paste(threshvec,collapse=','))) + } + a = (y == 1) + b = (y == 0) + na = sum(a) + nb = sum(b) + d = x2 - x1 + ### NRI + n.thresh=length(threshvec)-1 + risk.class.x1.ev=cut2(x1[a],threshvec) + risk.class.x2.ev=cut2(x2[a],threshvec) + thresh=c() + lt = length(threshvec) + for (i in 1:(lt-1)) { + thresh[i] = paste("[",toString(threshvec[i]),",",toString(threshvec[i+1]),"]") + } + levels(risk.class.x1.ev)=thresh + levels(risk.class.x2.ev)=thresh + cM.ev=confusionMatrix(risk.class.x2.ev,risk.class.x1.ev) + pup.ev=0 + pdown.ev=0 + for (i in 1:(n.thresh-1)) { pup.ev = pup.ev + sum(cM.ev\$table[(i+1):n.thresh,i])} + for (i in 2:n.thresh) { pdown.ev = pdown.ev + sum(cM.ev\$table[1:(i-1),i])} + pup.ev=pup.ev/na + pdown.ev=pdown.ev/na + risk.class.x1.ne=cut2(x1[b],threshvec) + risk.class.x2.ne=cut2(x2[b],threshvec) + levels(risk.class.x1.ne)=thresh + levels(risk.class.x2.ne)=thresh + cM.ne=confusionMatrix(risk.class.x2.ne,risk.class.x1.ne) + pup.ne=0 + pdown.ne=0 + for (i in 1:(n.thresh-1)){pup.ne=pup.ev+sum(cM.ne\$table[(i+1):n.thresh,i])} + for (i in 2:n.thresh){pdown.ne=pdown.ne+sum(cM.ne\$table[1:(i-1),i])} + pdown.ne=pdown.ne/nb + pup.ne=pup.ne/nb + nri = pup.ev - pdown.ev - (pup.ne - pdown.ne) + se.nri = sqrt((pup.ev + pdown.ev)/na + (pup.ne + pdown.ne)/nb) + z.nri = nri/se.nri + nri.ev = pup.ev - pdown.ev + se.nri.ev = sqrt((pup.ev + pdown.ev)/na) + z.nri.ev = nri.ev/se.nri.ev + nri.ne = pdown.ne - pup.ne + se.nri.ne = sqrt((pdown.ne + pup.ne)/nb) + z.nri.ne = nri.ne/se.nri.ne + ### Category Free NRI calculations + cfpup.ev = mean(d[a] > 0) + cfpup.ne = mean(d[b] > 0) + cfpdown.ev = mean(d[a] < 0) + cfpdown.ne = mean(d[b] < 0) + cfnri = cfpup.ev - cfpdown.ev - (cfpup.ne - cfpdown.ne) + se.cfnri = sqrt((cfpup.ev + cfpdown.ev)/na + (cfpup.ne + cfpdown.ne)/nb) + z.cfnri = cfnri/se.cfnri + cfnri.ev = cfpup.ev - cfpdown.ev + se.cfnri.ev = sqrt((cfpup.ev + cfpdown.ev)/na) + z.cfnri.ev = cfnri.ev/se.cfnri.ev + cfnri.ne = cfpdown.ne - cfpup.ne + se.cfnri.ne = sqrt((cfpdown.ne + cfpup.ne)/nb) + z.cfnri.ne = cfnri.ne/se.cfnri.ne + ### IDI calculations + improveSens = sum(d[a])/na + improveSpec = -sum(d[b])/nb + idi.ev = mean(improveSens) + idi.ne = mean(improveSpec) + idi = idi.ev - idi.ne + var.ev = var(d[a])/na + se.idi.ev = sqrt(var.ev) + z.idi.ev = idi.ev/se.idi.ev + var.ne = var(d[b])/nb + se.idi.ne = sqrt(var.ne) + z.idi.ne = idi.ne/se.idi.ne + se.idi = sqrt(var.ev + var.ne) + z.idi = idi/se.idi + ### AUC calculations + roc.x1 = roc(y, x1) + auc.x1 = auc(roc.x1) + ci.auc.x1 = ci.auc(roc.x1) + se.auc.x1 = (ci.auc.x1[3] - auc.x1)/qnorm(0.975) + roc.x2 = roc(y, x2) + auc.x2 = auc(roc.x2) + ci.auc.x2 = ci.auc(roc.x2) + se.auc.x2 = (ci.auc.x2[3] - auc.x2)/qnorm(0.975) + roc.test.x1.x2 = roc.test(roc.x1, roc.x2) ###Uses the default Delong method + sens.x1 = roc.x1\$sensitivities + spec.x1 = 1 - roc.x1\$specificities + n.x1 = length(sens.x1) + x1 = roc.x1\$thresholds + x1 = x1[c(-1,-n.x1)] + sens.x1 = sens.x1[c(-1,-n.x1)] + spec.x1 = spec.x1[c(-1,-n.x1)] + sens.x2 = roc.x2\$sensitivities + spec.x2 = 1 - roc.x2\$specificities + n.x2 = length(sens.x2) + x2 = roc.x2\$thresholds + x2 = x2[c(-1,-n.x2)] + sens.x2 = sens.x2[c(-1,-n.x2)] + spec.x2 = spec.x2[c(-1,-n.x2)] + ### Integrated sensitivity and 1-specificity calculations + is.x1 = trapz(x = x1, y = sens.x1) ### area under curves (relates to integrated sens, 1-spec) + is.x2 = trapz(x = x2, y = sens.x2) + ip.x1 = trapz(x = x1, y = spec.x1) + ip.x2 = trapz(x = x2, y = spec.x2) + + ### Output + output = c(n, na, nb, pup.ev, pup.ne, pdown.ev, pdown.ne, nri, se.nri, z.nri, + nri.ev, se.nri.ev, z.nri.ev, nri.ne, se.nri.ne, z.nri.ne, + cfpup.ev, cfpup.ne, cfpdown.ev, cfpdown.ne, cfnri, se.cfnri, z.cfnri, + cfnri.ev, se.cfnri.ev, z.cfnri.ev, cfnri.ne, se.cfnri.ne, z.cfnri.ne, + improveSens, improveSpec, idi.ev, se.idi.ev, z.idi.ev, idi.ne, + se.idi.ne, z.idi.ne, idi, se.idi, z.idi, is.x1, NA, is.x2, NA, + ip.x1, NA, ip.x2, NA, auc.x1, se.auc.x1, auc.x2, se.auc.x2, + roc.test.x1.x2\$p.value,incidence) + names(output) = c("n", "na", "nb", "pup.ev", "pup.ne", "pdown.ev", "pdown.ne", + "nri", "se.nri", "z.nri", "nri.ev", "se.nri.ev", "z.nri.ev", + "nri.ne", "se.nri.ne", "z.nri.ne", + "cfpup.ev", "cfpup.ne", "cfpdown.ev", "cfpdown.ne", + "cfnri", "se.cfnri", "z.cfnri", "cfnri.ev", "se.cfnri.ev", "z.cfnri.ev", + "cfnri.ne", "se.cfnri.ne", "z.cfnri.ne", "improveSens", "improveSpec", + "idi.ev", "se.idi.ev", "z.idi.ev", "idi.ne", "se.idi.ne", + "z.idi.ne", "idi", "se.idi", "z.idi", "is.x1", "se.is.x1", + "is.x2", "se.is.x2", "ip.x1", "se.ip.x1", "ip.x2", "se.ip.x2", + "auc.x1", "se.auc.x1", "auc.x2", "se.auc.x2", + "roc.test.x1.x2.pvalue","incidence") + resdf = data.frame(N=n, Na=na, Nb=nb, pup.ev=pup.ev, pup.ne=pup.ne, pdown.ev=pdown.ev, pdown.ne=pdown.ne, NRI=nri, NRI.se=se.nri, NRI.z=z.nri, + NRI.ev=nri.ev, NRI.ev.se=se.nri.ev, NRI.ev.z=z.nri.ev, NRI.ne=nri.ne, NRI.ne.se=se.nri.ne, NRI.ne.z=z.nri.ne, + cfpup.ev=cfpup.ev, cfpup.ne=cfpup.ne, cfpdown.ev=cfpdown.ev, cfpdown.ne=cfpdown.ne, CFNRI=cfnri, CFNRI.se=se.cfnri, CFNRI.z=z.cfnri, + CFNRI.ev=cfnri.ev, CFNRI.ev.se=se.cfnri.ev, CFNRI.ev.z=z.cfnri.ev, CFNRI.ne=cfnri.ne, CFNRI.ne.se=se.cfnri.ne, CFNRI.ne.z=z.cfnri.ne, + improvSens=improveSens, improvSpec=improveSpec, IDI.ev=idi.ev, IDI.ev.se=se.idi.ev, IDI.ev.z=z.idi.ev, IDI.ne=idi.ne, + IDI.ne.se=se.idi.ne, IDI.ne.z=z.idi.ne, IDI=idi, IDI.se=se.idi, IDI.z=z.idi, isx1=is.x1, isx2=is.x2, + ipxi=ip.x1, ipx2=ip.x2, AUC.x1=auc.x1, AUC.x1.se=se.auc.x1, AUC.x2=auc.x2, AUC.x2.se=se.auc.x2, + roctestpval=roc.test.x1.x2\$p.value,incidence=incidence) + tr = t(resdf) + tresdf = data.frame(measure=colnames(resdf),value=tr[,1]) + return(list(resdf=tresdf,output=output)) +} + + +### More comprehensive summary statistics from a raplot +### Choice of confidence intervals determined through asymptotics or bootstrapping (n.boot = ### of bootstrap resamples) +### dp is number of decimal places for results table + +summary.raplot = function(x1, x2, y, threshvec, cis = "boot", conf.level = 0.95, n.boot = 2000, dp = 4, stat_ra=NA) +{ + results = stat_ra + if (cis == "boot") { + print.noquote("Bootstrap estimates for SE") + results.boot = matrix(NA, n.boot, length(names(results))) + + colnames(results.boot) = names(results) + + for (i in 1:n.boot) { + ###boot.index = sample(length(cc.status), replace = TRUE) + ###risk.model1.boot = risk.model1[boot.index] + ###risk.model2.boot = risk.model2[boot.index] + ###cc.status.boot = cc.status[boot.index] + boot.index = sample(length(y), replace = TRUE) + risk.model1.boot = x1[boot.index] + risk.model2.boot = x2[boot.index] + cc.status.boot = y[boot.index] + r = statistics.raplot(x1 = risk.model1.boot, x2 = risk.model2.boot, y = cc.status.boot) + results.boot[i, ] = r\$output + } + + results.se.boot = apply(results.boot, 2, sd) + print(paste(results.se.boot,collapse=',')) + + + results[grep("se", names(results))] = results.se.boot[grep("se", names(results)) - 1] + + } + + + + ### calculate cis and return + + z = abs(qnorm((1 - conf.level)/2)) + + results.matrix = matrix(NA, 24, 2) + + results.matrix[1, ] = c("Total (n)", results["n"]) + results.matrix[2, ] = c("Events (n)", results["na"]) + results.matrix[3, ] = c("Non-events (n)", results["nb"]) + results.matrix[4, ] = c("Category free NRI and summary statistics","-------------------------") + results.matrix[5, ] = c("cfNRI events (%)", + paste(round(100*results["cfnri.ev"], dp-2), " (", + round(100*results["cfnri.ev"] - z * 100*results["se.cfnri.ev"], dp-2), + " to ", round(100*results["cfnri.ev"] + + z * 100*results["se.cfnri.ev"], dp-2), ")", sep = "")) + results.matrix[6, ] = c("cfNRI non-events (%)", + paste(round(100*results["cfnri.ne"], dp-2), " (", + round(100*results["cfnri.ne"] - z * 100*results["se.cfnri.ne"], dp)-2, + " to ", round(100*results["cfnri.ne"] + z * 100*results["se.cfnri.ne"], + dp-2), ")", sep = "")) + results.matrix[7, ] = c("cfNRI (%)", + paste(round(100*results["cfnri"], dp-2), " (", + round(100*results["cfnri"] - z * 100*results["se.cfnri"], dp-2), + " to ", round(100*results["cfnri"] + z * 100*results["se.cfnri"], + dp-2), ")", sep = "")) + results.matrix[8, ] = c("NRI and summary statistics","-------------------------") + results.matrix[9, ] = c("NRI events (%)", + paste(round(100*results["nri.ev"], dp-2), " (", + round(100*results["nri.ev"] - z * 100*results["se.nri.ev"], dp-2), + " to ", round(100*results["nri.ev"] + + z * 100*results["se.nri.ev"], dp-2), ")", sep = "")) + results.matrix[10, ] = c("NRI non-events (%)", + paste(round(100*results["nri.ne"], dp-2), " (", + round(100*results["nri.ne"] - z * 100*results["se.nri.ne"], dp-2), + " to ", round(100*results["nri.ne"] + z * 100*results["se.nri.ne"], + dp-2), ")", sep = "")) + results.matrix[11, ] = c("NRI (%)", + paste(round(100*results["nri"], dp-2), " (", + round(100*results["nri"] - z * 100*results["se.nri"], dp-2), + " to ", round(100*results["nri"] + z * 100*results["se.nri"], + dp-2), ")", sep = "")) + results.matrix[12, ] = c("IDI and summary statistics","-------------------------") + results.matrix[13, ] = c("IDI events", + paste(round(results["idi.ev"], dp), " (", + round(results["idi.ev"] - z * results["se.idi.ev"], dp), + " to ", round(results["idi.ev"] + z * results["se.idi.ev"], + dp), ")", sep = "")) + results.matrix[14, ] = c("IDI non-events", + paste(round(results["idi.ne"], dp), " (", + round(results["idi.ne"] - z * results["se.idi.ne"], dp), + " to ", round(results["idi.ne"] + z * results["se.idi.ne"], + dp), ")", sep = "")) + results.matrix[15, ] = c("IDI", + paste(round(results["idi"], dp), " (", + round(results["idi"] - z * results["se.idi"], dp), + " to ", round(results["idi"] + z * results["se.idi"], + dp), ")", sep = "")) + results.matrix[16, ] = c("IS (null model)", + paste(round(results["is.x1"], dp), " (", + round(results["is.x1"] - z * results["se.is.x1"], dp), + " to ", round(results["is.x1"] + z * results["se.is.x1"], + dp), ")", sep = "")) + results.matrix[17, ] = c("IS (alt model)", + paste(round(results["is.x2"], dp), " (", + round(results["is.x2"] - z * results["se.is.x2"], dp), + " to ", round(results["is.x2"] + z * results["se.is.x2"], + dp), ")", sep = "")) + results.matrix[18, ] = c("IP (null model)", + paste(round(results["ip.x1"], dp), " (", + round(results["ip.x1"] - z * results["se.ip.x1"], dp), + " to ", round(results["ip.x1"] + z * results["se.ip.x1"], + dp), ")", sep = "")) + results.matrix[19, ] = c("IP (alt model)", + paste(round(results["ip.x2"], dp), " (", + round(results["ip.x2"] - z * results["se.ip.x2"], dp), + " to ", round(results["ip.x2"] + z * results["se.ip.x2"], + dp), ")", sep = "")) + results.matrix[20, ] = c("AUC","-------------------------") + results.matrix[21, ] = c("AUC (null model)", + paste(round(results["auc.x1"], dp), " (", + round(results["auc.x1"] - z * results["se.auc.x1"], dp), + " to ", round(results["auc.x1"] + z * results["se.auc.x1"], + dp), ")", sep = "")) + results.matrix[22, ] = c("AUC (alt model)", + paste(round(results["auc.x2"], dp), " (", + round(results["auc.x2"] - z * results["se.auc.x2"], dp), + " to ", round(results["auc.x2"] + z * results["se.auc.x2"], + dp), ")", sep = "")) + results.matrix[23, ] = c("difference (P)", round(results["roc.test.x1.x2.pvalue"], dp)) + results.matrix[24, ] = c("Incidence", round(results["incidence"], dp)) + + return(results.matrix) +} + + + +]]> + +options(width=120) +options(digits=5) +logf = file("rgNRI.log", open = "a") +sink(logf,type = c("output", "message")) +Out_Dir = "$html_file.files_path" +Input1 = "$input1" +Input2 = "$input2" +myTitle = "$title" +outtab = "$nri_file" +input1_obs = $input1_observed +input1_pred = $input1_predicted +input1_id = $input1_id +input2_obs = $input2_observed +input2_pred = $input2_predicted +input2_id = $input2_id +in1 = read.table(Input1,head=T,sep='\t') +in2 = read.table(Input2,head=T,sep='\t') +id1 = in1[,input1_id] +id2 = in2[,input2_id] +useme1 = in1[which(id1 %in% id2),] +useme2 = in2[which(id2 %in% id1),] +id1 = useme1[,input1_id] +id2 = useme2[,input2_id] +useme1 = useme1[order(id1),] +useme2 = useme2[order(id2),] +x1 = useme1[,input1_pred] +x2 = useme2[,input2_pred] +y1 = useme1[,input1_obs] +y2 = useme2[,input2_obs] +n.boot = $CImeth.nboot +conf.level = 0.95 +cis = "$CImeth.cis" +digits = 4 +nydiff = sum(y1 != y2) +if (nydiff > 0) { + print.noquote(paste('Input error: observed status column has',nydiff,'differences - cannot reliably proceed')) + quit(save="no",status=1) + } +y = y2 +outplot = 'rgNRI_EventRisk.pdf' +res = raplot(x1=x1, x2=x2, y=y, outplot=outplot,title=myTitle) + +stats = statistics.raplot(x1=x1, x2=x2, y=y) +res1 = stats\$resdf +out1 = stats\$output +print.noquote('Results:') +print.noquote(res1,digits=4) +res2 = summary.raplot(x1=x1, x2=x2, y=y, cis = cis, conf.level = conf.level, n.boot = n.boot, dp = digits, stat_ra=out1) +print.noquote('Summary:') +print.noquote(res2,digits=4) +write.table(format(res1,digits=4),outtab,quote=F, col.names=F,sep="\t",row.names=F) +print.noquote('SessionInfo for this R session:') +sessionInfo() +print.noquote('warnings for this R session:') +warnings() +sink() +</configfile> +</configfiles> + <inputs> + <param name="title" type="text" value="NRI test" label="Plot Title" help="Will appear as the title for the comparison plot"/> + <param name="input1" type="data" format="tabular" label="Select a tabular file from the baseline model with predicted and observed outcome column for each subject" + multiple='False' help="Observed and predicted status columns must be selected from this file below - NOTE both models must be in same order with exact matches in all observed outcomes" optional="False"/> + <param name="input1_observed" label="Select column containing observed outcome (0 for no event, 1 for an event)" type="data_column" data_ref="input1" numerical="True" + multiple="False" use_header_names="True" optional="False" help = "Observed outcomes are compared in the two files to check that the datasets are from the same data"/> + <param name="input1_predicted" label="Select column containing predicted event probabilies from baseline model" type="data_column" data_ref="input1" numerical="True" + multiple="False" use_header_names="True" optional="False" help="Must be in range 0-1"/> + <param name="input1_id" label="Select column containing subject ID from baseline model" type="data_column" data_ref="input1" numerical="True" + multiple="False" use_header_names="True" optional="False" help="Subect IDs are needed to match subjects to compare predictions in the two inputs"/> + <param name="input2" type="data" format="tabular" label="Select a tabular file from the new model with predicted and observed outcome columns for each subject" + multiple='False' help="Observed and predicted status columns must be selected from this file below" /> + <param name="input2_observed" label="Select column containing observed outcome (0 for no event, 1 for an event)" type="data_column" data_ref="input2" numerical="True" + multiple="False" use_header_names="True" optional="False" help = "Observed outcomes are compared in the two files to check that the datasets are from the same data"/> + <param name="input2_predicted" label="Select column containing predicted event probabilities from the new model" type="data_column" data_ref="input2" numerical="True" + multiple="False" use_header_names="True" optional="False" help="Must be in range 0-1"/> + <param name="input2_id" label="Select column containing subject ID from the new model" type="data_column" data_ref="input2" numerical="True" + multiple="False" use_header_names="True" optional="False" help="Subect IDs are needed to match subjects to compare predictions in the two inputs"/> + <conditional name="CImeth"> + <param name="cis" type="select" label="CI calculation method" + help="Bootstrap will take time - a long time for thousands - asymptotic is quick and informative"> + <option value="asymptotic" selected="true">Asymptotic estimate</option> + <option value="boot">Bootstrap for empirical CIs</option> + </param> + <when value="boot"> + <param name="nboot" type="integer" value="1000" label="Number of bootstrap replicates"/> + </when> + <when value="asymptotic"> + <param name="nboot" type="hidden" value="1000"/> + </when> + </conditional> + </inputs> + <outputs> + <data format="html" name="html_file" label="${title}.html"/> + <data format="tabular" name="nri_file" label="${title}_nrires.xls"/> + </outputs> + <tests> + <test> + <param name='title' value='nri_test1' /> + <param name='input1' value='nri_test1.xls' ftype='tabular' /> + <param name='input2' value='nri_test1.xls' ftype='tabular' /> + <param name='input1_id' value="1" /> + <param name='input1_observed' value="2" /> + <param name='input1_predicted' value="3" /> + <param name='input2_observed' value="2" /> + <param name='input2_predicted' value="4" /> + <output name='html_file' file='nri_test1_out.html' compare='diff' lines_diff='10' /> + <output name='nri_file' file='nri_test1_out.xls' /> + </test> +</tests> +<help> + +**Before you start** + +This is a simple tool to calculate various measures of improvement in prediction between two models described in pickering_paper_ +It is based on an R script pickering_code_ written by Dr John W Pickering and Dr David Cairns from sunny Otago University which +has been debugged and slightly adjusted to fit a Galaxy tool wrapper. + + +**What it does** + +Copied from the documentation in pickering_code_ :: + + + Functions to create risk assessment plots and associated summary statistics + + + (c) 2012 Dr John W Pickering, john.pickering@otago.ac.nz, and Dr David Cairns + Last modified August 2014 + + Redistribution and use in source and binary forms, with or without + modification, are permitted provided that the following conditions are met: + * Redistributions of source code must retain the above copyright + notice, this list of conditions and the following disclaimer. + * Redistributions in binary form must reproduce the above copyright + notice, this list of conditions and the following disclaimer in + the documentation and/or other materials provided with the distribution + + FUNCTIONS + raplot + Produces a Risk Assessment Plot and outputs the coordinates of the four curves + Based on: Pickering, J. W. and Endre, Z. H. (2012). New Metrics for Assessing Diagnostic Potential of + Candidate Biomarkers. Clinical Journal of the American Society of Nephrology, 7, 1355–1364. doi:10.2215/CJN.09590911 + + statistics.raplot + Produces the NRIs, IDIs, IS, IP, AUCs. + Based on: Pencina, M. J., D'Agostino, R. B. and Steyerberg, E. W. (2011). Extensions of net reclassification improvement calculations to + measure usefulness of new biomarkers. Statistics in Medicine, 30(1), 11–21. doi:10.1002/sim.4085 + Pencina, M. J., D'Agostino, R. B. and Vasan, R. S. (2008). Evaluating the added predictive ability of a new marker: From area under the + ROC curve to reclassification and beyond. + Statistics in Medicine, 27(2), 157–172. doi:10.1002/sim.2929 + DeLong, E., DeLong, D. and Clarke-Pearson, D. (1988). Comparing the areas under 2 or more correlated receiver operating characteristic curves - a nonparametric approach. + Biometrics, 44(3), 837–845. + + summary.raplot + Produces the NRIs, IDIs, IS, IP, AUCs with confidence intervals using a bootstrap or asymptotic procedure. (I prefer bootstrap which is chosed by cis=c("boot")) + + + Required arguments for all functions: + x1 is calculated risk (eg from a glm) for the null model, i.e. predict(,type="response") on a glm object + x2 is calculated risk (eg from a glm) for the alternative model + y is the case-control indicator (0 for controls, 1 for cases) + Optional argument + t are the boundaries of the risks for each group (ie 0, 1 and the thresholds beteween. eg c(0,0,3,0,7,1)). If missing, defaults to c(0, the incidence, 1) + + +**Input** + +The observed and predicted outcomes from two models to be compared. + +**Output** + +Lots'o'measures (TM) see pickering_paper_ for details + +**Attributions** + +pickering_paper_ is the paper the caclulations performed by this tool is based on + +pickering_code_ is the R function from John Pickering exposed by this Galaxy tool with minor modifications and hacks by Ross Lazarus. + +Galaxy_ (that's what you are using right now!) for gluing everything together + +Otherwise, all code and documentation comprising this tool was written by Ross Lazarus and is +licensed to you under the LGPL_ like other rgenetics artefacts + +.. _LGPL: http://www.gnu.org/copyleft/lesser.html +.. _pickering_code: http://www.researchgate.net/publication/264672640_R_function_for_Risk_Assessment_Plot__reclassification_metrics_NRI_IDI_cfNRI +.. _pickering_paper: http://cjasn.asnjournals.org/content/early/2012/05/24/CJN.09590911.full +.. _Galaxy: http://getgalaxy.org + + +</help> + +<citations> + <citation type="doi">doi: 10.2215/CJN.09590911</citation> +</citations> +</tool> + + +
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/rglasso_cox.xml Sat Oct 31 02:26:24 2015 -0400 @@ -0,0 +1,910 @@ +<tool id="rglasso_cox" name="Lasso" version="0.03"> + <description>and cox regression using elastic net</description> + <requirements> + <requirement type="package" version="3.2.2">R_3_2_2</requirement> + <requirement type="package" version="1.3.18">graphicsmagick</requirement> + <requirement type="package" version="9.10">ghostscript</requirement> + <requirement type="package" version="3.2">glmnet_lars_3_2</requirement> + </requirements> + <command interpreter="python"> + rgToolFactory.py --script_path "$runme" --interpreter "Rscript" --tool_name "rglasso" + --output_dir "$html_file.files_path" --output_html "$html_file" --make_HTML "yes" + </command> +<configfiles> +<configfile name="runme"> +<![CDATA[ +library('glmnet') +library('lars') +library('survival') +library('pec') + + +message=function(x) {print.noquote(paste(x,sep=''))} + + +ross.cv.glmnet = function (x, y, weights, offset = NULL, lambda = NULL, type.measure = c("mse", + "deviance", "class", "auc", "mae"), nfolds = 10, foldid, + grouped = TRUE, keep = FALSE, parallel = FALSE, ...) +{ + if (missing(type.measure)) + type.measure = "default" + else type.measure = match.arg(type.measure) + if (!is.null(lambda) && length(lambda) < 2) + stop("Need more than one value of lambda for cv.glmnet") + N = nrow(x) + if (missing(weights)) + weights = rep(1, N) + else weights = as.double(weights) + y = drop(y) + glmnet.call = match.call(expand.dots = TRUE) + sel = match(c("type.measure", "nfolds", "foldid", "grouped", + "keep"), names(glmnet.call), F) + if (any(sel)) + glmnet.call = glmnet.call[-sel] + glmnet.call[[1]] = as.name("glmnet") + glmnet.object = glmnet(x, y, weights = weights, offset = offset, + lambda = lambda, ...) + glmnet.object\$call = glmnet.call + is.offset = glmnet.object\$offset + lambda = glmnet.object\$lambda + if (inherits(glmnet.object, "multnet")) { + nz = predict(glmnet.object, type = "nonzero") + nz = sapply(nz, function(x) sapply(x, length)) + nz = ceiling(apply(nz, 1, median)) + } + else nz = sapply(predict(glmnet.object, type = "nonzero"), + length) + if (missing(foldid)) + foldid = sample(rep(seq(nfolds), length = N)) + else nfolds = max(foldid) + if (nfolds < 3) + stop("nfolds must be bigger than 3; nfolds=10 recommended") + outlist = as.list(seq(nfolds)) + if (parallel && require(foreach)) { + outlist = foreach(i = seq(nfolds), .packages = c("glmnet")) %dopar% + { + sel = foldid == i + if (is.matrix(y)) + y_sub = y[!sel, ] + else y_sub = y[!sel] + if (is.offset) + offset_sub = as.matrix(offset)[!sel, ] + else offset_sub = NULL + glmnet(x[!sel, , drop = FALSE], y_sub, lambda = lambda, + offset = offset_sub, weights = weights[!sel], + ...) + } + } + else { + for (i in seq(nfolds)) { + sel = foldid == i + if (is.matrix(y)) + y_sub = y[!sel, ] + else y_sub = y[!sel] + if (is.offset) + offset_sub = as.matrix(offset)[!sel, ] + else offset_sub = NULL + outlist[[i]] = glmnet(x[!sel, , drop = FALSE], + y_sub, lambda = lambda, offset = offset_sub, + weights = weights[!sel], ...) + } + } + fun = paste("cv", class(glmnet.object)[[1]], sep = ".") + cvstuff = do.call(fun, list(outlist, lambda, x, y, weights, + offset, foldid, type.measure, grouped, keep)) + cvm = cvstuff\$cvm + cvsd = cvstuff\$cvsd + cvname = cvstuff\$name + + out = list(lambda = lambda, cvm = cvm, cvsd = cvsd, cvup = cvm + + cvsd, cvlo = cvm - cvsd, nzero = nz, name = cvname, glmnet.fit = glmnet.object) + if (keep) + out = c(out, list(fit.preval = cvstuff\$fit.preval, foldid = foldid)) + + lamin = if (type.measure == "auc") + getmin(lambda, -cvm, cvsd) + else getmin(lambda, cvm, cvsd) + out = c(out, as.list(lamin)) + hitsse = rep(0,ncol(x)) + hitsmin = rep(0,ncol(x)) + names(hitsse) = colnames(x) + names(hitsmin) = colnames(x) + olmin = lamin\$lambda.min + ol1sd = lamin\$lambda.1se + lambs = c(olmin,ol1sd) + names(lambs) = c('olmin','ol1sd') + for (cvfit in outlist) { + colmin = which(cvfit\$lambda == olmin) + col1se = which(cvfit\$lambda == ol1sd) + nzmin = which(cvfit\$beta[,colmin] != 0) + nz1se = which(cvfit\$beta[,col1se] != 0) + hitsse[nz1se] = hitsse[nz1se] + 1 + hitsmin[nzmin] = hitsmin[nzmin] + 1 + } + obj = c(out,list(cvhits.1se=hitsse,cvhits.min=hitsmin)) + class(obj) = "cv.glmnet" + obj +} + +mdsPlot = function(dm,myTitle,groups=NA,outpdfname,transpose=T) +{ + + samples = colnames(dm) + mt = myTitle + pcols=c('maroon') + if (! is.na(groups)) + { + gu = unique(groups) + colours = rainbow(length(gu),start=0.1,end=0.9) + pcols = colours[match(groups,gu)] + } + mydata = dm + if (transpose==T) + { + mydata = t(dm) + } + npred = ncol(mydata) + d = dist(mydata) + fit = cmdscale(d,eig=TRUE, k=min(10,npred-2)) + xmds = fit\$points[,1] + ymds = fit\$points[,2] + pdf(outpdfname) + plot(xmds, ymds, xlab="Dimension 1", ylab="Dimension 2", + main=paste(mt,"MDS Plot"),type="n", col=pcols, cex=0.35) + text(xmds, ymds, labels = row.names(mydata), cex=0.35, col=pcols) + grid(col="lightgray",lty="dotted") + dev.off() +} + + +getpredp_logistic = function(x,yvec,yvarname,id) +{ + yvals = unique(yvec) + if (length(yvals) != 2) { + message(c('ERROR: y does not have 2 values =',paste(yvals,collapse=','))) + return(NA) + } + cols = colnames(x) + if (length(cols) == 0) { + message('ERROR: No columns in input x? Cannot predict!') + return(NA) + } + cn = paste(cols, collapse = ' + ') + + formstring=paste("y ~",cn) + form = as.formula(formstring) + ok = complete.cases(x) + + if (sum(ok) < length(ok)) { + x = x[ok,] + yvec = yvec[ok] + id = id[ok] + } + nx = data.frame(id=id,x,y=yvec) + print('nx,yvec:') + print(head(nx,n=3)) + print(yvec) + mdl = glm(form, data=nx, family="binomial", na.action=na.omit) + message(c('Model format =',formstring)) + message(paste('Predictive model details used to generate logistic outcome probabilities for',yvarname,':')) + print(summary(md1)) + print(anova(md1)) + predp = predict(md1,nx,type="response") + p1 = data.frame(id=id,pred_response=predp,obs_response=yvec) + return(p1) +} + +getpredp_cox = function(x,time,status,id,predict_at) +{ + cols = colnames(x) + if (length(cols) == 0) { + message('ERROR: No columns in input x? Cannot predict!') + return(NA) + } + cn = paste(colnames(x), collapse = ' + ') + + formstring=paste("Surv(time, status) ~",cn) + + form = as.formula(formstring) + + ok = complete.cases(x) + + if (sum(ok) < length(ok)) { + x = x[ok,] + time = time[ok] + status = status[ok] + id = id[ok] + } + nx = data.frame(x,time=time,status=status) + m1 = coxph(form, data=nx,singular.ok=TRUE) + print.noquote('Predictive model details used to generate survival probabilities:') + print.noquote(m1) + predpq = predictSurvProb(object=m1, newdata=nx, times=predict_at) + predpq = 1-predpq + colnames(predpq) = paste('p_surv_to',predict_at,sep='_') + p1 = data.frame(id=id,predpq,time=time,status=status) + return(p1) +} + + +dolasso_cox = function(x,y,debugOn=F,maxsteps=10000,nfold=10,xcolnames,ycolnames,optLambda='lambda.1se',out_full=F,out_full_file=NA, + out_pred=F,out_pred_file=NA,cox_id=NA, descr='Cox test',do_standard=F,alpha=0.9,penalty,predict_at,mdsplots=F) +{ + logf = file("cox_rglasso.log", open = "a") + sink(logf,type = c("output", "message")) + res = NULL + if (mdsplots==T) { + outpdfname = 'cox_x_in_sample_space_MDS.pdf' + p = try({ mdsPlot(x,'measurements in sample space',groups=NA,outpdfname=outpdfname,transpose=T) },T) + if (class(p) == "try-error") + { + print.noquote(paste('Unable to produce predictors in sample space mds plot',p)) + } + outpdfname = 'cox_samples_in_x_space_MDS.pdf' + p = try({mdsPlot(x,'samples in measurement space',groups=y,outpdfname=outpdfname,transpose=F) },T) + if (class(p) == "try-error") + { + print.noquote(paste('Unable to produce samples in measurement space mds plots',p)) + } + } + if (is.na(predict_at)) { predict_at = quantile(y) } + message(paste('@@@ Cox model will be predicted at times =',paste(predict_at,collapse=','))) + do_standard = do_standard + standardize = do_standard + normalize = do_standard + p = try({larsres = glmnet(x,y,family='cox',standardize=standardize,alpha=alpha,penalty.factor=penalty )},T) + if (class(p) == "try-error") + { + print.noquote('Unable to run cox glmnet on your data') + print.noquote(p) + sink() + return(NA) + } + if (out_full == T) + { + b = as.matrix(larsres\$beta) + nb = length(colnames(b)) + bcoef = b[,nb] + lastl = larsres\$lambda[length(larsres\$lambda)] + allres = data.frame(x=rownames(b),beta=bcoef,lambda=lastl) + write.table(format(allres,digits=5),out_full_file,quote=FALSE, sep="\t",row.names=F) + } + + outpdf = paste('cox',descr,'glmnetdev.pdf',sep='_') + try( + { + pdf(outpdf) + plot(larsres,main='cox glmnet',label=T) + grid() + dev.off() + },T) + + larscv = NA + + p = try({larscv=ross.cv.glmnet(x,y,family=fam,type.measure='deviance',penalty=penalty)},T) + if (class(p) == "try-error") { + print.noquote(paste('Unable to cross validate your data',p)) + sink() + return(NA) + } + lse = larscv\$cvhits.1se + lmin = larscv\$cvhits.min + tot = lse + lmin + allhits = data.frame(hits_lambda_1se = lse,hits_lambda_min = lmin) + nzhits = allhits[which(tot != 0),] + message('Times each predictor was selected in CV models (excluding zero count predictors)') + print.noquote(nzhits) + out_nz_file = 'cox_cross_validation_model_counts.xls' + write.table(nzhits,out_nz_file,quote=FALSE, sep="\t",row.names=F) + + outpdf = paste('cox',descr,'glmnet_cvdeviance.pdf',sep='_') + + p = try( + { + pdf(outpdf) + plot(larscv,main='Deviance',label=T) + grid() + dev.off() + },T) + if (optLambda == 'lambda.min') { + best_lambda = larscv\$lambda.min + bestcoef = as.matrix(coef(larscv, s = "lambda.min")) + } else { + best_lambda = larscv\$lambda.1se + bestcoef = as.matrix(coef(larscv, s = "lambda.1se")) + } + inmodel = which(bestcoef != 0) + coefs = bestcoef[inmodel] + preds = rownames(bestcoef)[inmodel] + + names(coefs) = preds + pen = as.logical( ! penalty[inmodel]) + if (out_pred==T) + { + if (length(inmodel) > 0 ) { + predcols = inmodel + xmat = as.matrix(x[,predcols]) + colnames(xmat) = preds + bestpred = getpredp_cox(x=xmat,time=y[,'time'],status=y[,'status'],id=cox_id, predict_at=predict_at) + pred = data.frame(responsep=bestpred, best_lambda=best_lambda,lamchoice=optLambda,alpha=alpha) + write.table(pred,out_pred_file,quote=FALSE, sep="\t",row.names=F) + } else { print.noquote('WARNING: No coefficients in selected model to predict with - no predictions made') } + } + if (debugOn) { + print.noquote(paste('best_lambda=',best_lambda,'saving cox respreds=',paste(names(coefs),collapse=','),'as predictors of survival. Coefs=',paste(coefs,collapse=','))) + } + p = try({res = data.frame(regulator=names(coefs),partial_likelihood=coefs,forced_in=pen,glmnet_model='cox',best_lambda=best_lambda, + lambdaChoice=optLambda,alpha=alpha)},T) + if (class(p) == "try-error") { + message(paste('@@@ unable to return a dataframe',p)) + sink() + return(NA) + } + print.noquote('@@@ Results preview:') + print.noquote(res,digits=5) + sink() + return(res) + +} + + +do_lasso = function(x=NA,y=NA,do_standard=T,debugOn=T,defaultFam="gaussian",optLambda='minLambda',descr='description', indx=1,target='target',sane=F, + alpha=0.9,nfold=10,penalty=c(),out_pred=F,out_pred_file='outpred',yvarname='yvar',id=c(),mdsplots=F) +{ + logf = file(paste(target,"rglasso.log",sep='_'), open = "a") + sink(logf,type = c("output", "message")) + res = NA + phe_is_bin = (length(unique(y)) == 2) + forcedin = paste(colnames(x)[which(penalty == 0)],collapse=',') + fam = "gaussian" + if (defaultFam %in% c("poisson","binomial","gaussian","multinomial")) fam=defaultFam + if (phe_is_bin == T) { + fam = "binomial" + } + print.noquote(paste('target=',target,'is binary=',phe_is_bin,'dim(x)=',paste(dim(x),collapse=','),'length(y)=',length(y),'force=',forcedin,'fam=',fam)) + standardize = do_standard + p = try({larsres = glmnet(x,y,family=fam,standardize=standardize,maxit=10000,alpha=alpha,penalty.factor=penalty) },T) + if (class(p) == "try-error") + { + print(paste('ERROR: unable to run glmnet for target',target,'error=',p)) + sink() + return(NA) + } + + mt = paste('Glmnet fraction deviance for',target) + outpdf = paste(target,'glmnetPath.pdf',sep='_') + pdf(outpdf) + plot(larsres,main=mt,label=T) + grid() + dev.off() + + outpdf = paste(target,'glmnetDeviance.pdf',sep='_') + + mt2 = paste('Glmnet lambda for',target) + + pdf(outpdf) + plot(larsres,xvar="lambda",main=mt2,label=T) + grid() + dev.off() + + larscv = NA + if (fam=="binomial") { + tmain = paste(target,'AUC') + outpdf = paste(target,'glmnetCV_AUC.pdf',sep='_') + p = try({larscv = ross.cv.glmnet(x=x,y=y,family=fam,type.measure='auc')},T) + } else { + tmain = paste(target,'CV MSE') + outpdf = paste(target,'glmnetCV_MSE.pdf',sep='_') + p = try({larscv = ross.cv.glmnet(x=x,y=y,family=fam,type.measure='mse')},T) + } + if (class(p) == "try-error") + { + print(paste('ERROR: unable to run cross validation for target',target,'error=',p)) + sink() + return(NA) + } + + pdf(outpdf) + plot(larscv,main=tmain) + grid() + dev.off() + + lse = larscv\$cvhits.1se + lmin = larscv\$cvhits.min + tot = lse + lmin + allhits = data.frame(pred=colnames(x),hits_lambda_1se = lse,hits_lambda_min = lmin) + nzhits = allhits[which(tot != 0),] + message('Total hit count for each predictor over all CV models (excluding zero count predictors)') + print.noquote(nzhits) + out_nz_file = paste(target,'cross_validation_model_counts.xls',sep='_') + write.table(nzhits,out_nz_file,quote=FALSE, sep="\t",row.names=F) + + ipenalty = c(0,penalty) + if (optLambda == 'lambda.min') { + best_lambda = larscv\$lambda.min + bestpred = as.matrix(coef(larscv, s = "lambda.min")) + } else { + best_lambda = larscv\$lambda.1se + bestpred = as.matrix(coef(larscv, s = "lambda.1se")) + } + inmodel = which(bestpred != 0) + coefs = bestpred[inmodel,1] + preds = rownames(bestpred)[inmodel] + iforced = ipenalty[inmodel] + forced = ! as.logical(iforced) + names(coefs) = preds + ncoef = length(coefs) - 1 + if (out_pred==T && fam=="binomial") + { + print.noquote(paste('Predicting',target,'probabilities from binomial glmnet at alpha',alpha,'and lambda',best_lambda)) + bestpred = predict.glmnet(larsres,s=best_lambda,newx=x,type="response") + bestpred = exp(bestpred)/(1+exp(bestpred)) + pred = data.frame(id=id,y=y,predp=as.vector(bestpred), best_lambda=best_lambda) + write.table(pred,out_pred_file,quote=FALSE, sep="\t",row.names=F) + } + if (debugOn) {cat(indx,'best_lambda=',best_lambda,'saving',fam,'respreds=',names(coefs),'as predictors of',target,'coefs=',coefs,'\n')} + res = try(data.frame(i=indx,pred=target,regulator=names(coefs),coef=coefs,forced_in=forced,glmnet_model=fam,ncoef=ncoef, + best_lambda=best_lambda,lambdaChoice=optLambda,alpha=alpha),T) + if (class(res) == "try-error") { + sink() + return(NA) } + print.noquote(res) + sink() + return(res) +} + + +dolasso_generic = function(predvars=NA,depvars=NA,debugOn=T,maxsteps=100, alpha=0.9,nfold=10,xcolnames=c(),ycolnames=c(),optLambda='minLambda', out_pred_file=NA, + descr="describe me",do_standard=F,defaultFam="gaussian",penalty=c(),out_pred=F,cox_id=c(),mdsplots=F,xfilt=0.95) +{ + logf = file("rglasso.log", open = "a") + sink(logf,type = c("output", "message")) + xdat = predvars + xm = data.matrix(xdat) + res = NULL + id = cox_id + depnames = ycolnames + ndep = length(depnames) + if (mdsplots==T) { + outpdfname = 'rglasso_x_in_sample_space_MDS.pdf' + p = try({ mdsPlot(xm,'measurements in sample space',groups=NA,outpdfname=outpdfname,transpose=T) },T) + if (class(p) == "try-error") + { + print.noquote(paste('Unable to produce predictors in sample space mds plot. Error:',p)) + } + outpdfname = 'rglasso_samples_in_x_space_MDS.pdf' + p = try({mdsPlot(xm,'samples in measurement space',groups=NA,outpdfname=outpdfname,transpose=F) },T) + if (class(p) == "try-error") + { + print.noquote(paste('Unable to produce samples in measurement space mds plot. Error:',p)) + } + } + ndat = nrow(xm) + cfracs = colSums(! is.na(xm))/ndat + keepme = (cfracs >= xfilt) + print.noquote(paste('Removing', sum(! keepme), 'xvars with more than',xfilt,'fraction missing')) + vars = apply(xm,2,var,na.rm=T) + xm = xm[,keepme] + for (i in c(1:max(1,ndep))) { + target = depnames[i] + if (length(target) < 1) { target='y' } + if (i %% 100 == 0) { cat(i,target,'\n') } + if (ndep <= 1) { + y=depvars + } else { + y = depvars[,i] + } + if (fam == "binomial") {y = as.factor(y)} + x = xm + id = cox_id + if (fam != "cox") { + ok = complete.cases(x,y) + if (sum(! ok) > 0) { + message(paste('@@@ Removing',sum(! ok),'cases with missing y of',length(y),'@@@')) + y = y[(ok)] + x = x[(ok),] + id = id[(ok)] + } + } + ok = complete.cases(y) + if (sum(ok) == 0 ) { + print(paste("No complete cases found for",target,"in input x dim =",paste(dim(xm),collapse=','),"length y=",length(y))) + } else { + if (i == 1) { outpred = out_pred_file + } else { + outpred = paste(target,'predicted_output.xls') + } + regres = do_lasso(x=x,y=y,do_standard=do_standard,debugOn=debugOn,defaultFam=defaultFam,optLambda=optLambda,out_pred_file=outpred, + descr=descr,indx=i,target=target,alpha=alpha,nfold=nfold,penalty=penalty,out_pred=out_pred,yvarname=target,id=id,mdsplots=mdsplots) + if (! is.na(regres)) { res = rbind(res,regres) } + } + } + print.noquote('@@@ Results preview:') + print.noquote(res,digits=5) + sink() + return(res) +} + + +corPlot=function(xdat=c(),main='main title',is_raw=T) +{ + library(pheatmap) + library(gplots) + if (is_raw) { + cxdat = cor(xdat,method="spearman",use="pairwise.complete.obs") + } else { + cxdat=xdat + } + xro = nrow(cxdat) + if (xro > 1000) stop("Too many rows for heatmap, who can read?!") + fontsize_col = 5.0 + pheatmap(cxdat, main=main, show_colnames = F, width=30, height=30, + fontsize_row=fontsize_col, border_color=NA) +} + + +runTest = function(n=10) +{ + set.seed (NULL) + Y = data.frame(y1=runif (n),y2=runif(n)) + Xv <- runif(n*n) + X <- matrix(Xv, nrow = n, ncol = n) + + mydf <- data.frame(Y, X) + + regres_out = dolasso_generic(predvars=X,depvars=Y,debugOn=T,p.cutoff = 0.05,maxsteps=10000,nfold=10, + descr='randomdata',do_standard=do_standard,defaultFam="gaussian",alpha=0.05) + return(regres_out) +} +]]> +options(width=512) +options(digits=5) +alpha = $alpha +nfold = $nfold +optLambda = "$optLambda" +Out_Dir = "$html_file.files_path" +Input = "$input1" +indat = read.table(Input,head=T,sep='\t') +datcols = colnames(indat) +myTitle = "$title" +outtab = "$model_file" +do_standard = as.logical("$do_standard") +mdsplots = as.logical("$mdsplots") +fam = "$model.fam" +xvar_cols_in = "$xvar_cols" +force_xvar_cols_in = "$force_xvar_cols" +xvar_cols = as.numeric(strsplit(xvar_cols_in,",")[[1]]) +force_xvar_cols = c() +penalties = rep(1,length(datcols)) +forced_in = NA + +logxform = "$logxform_cols" +if (logxform != "None") { + logxform_cols = as.numeric(strsplit(logxform,",")[[1]]) + if (length(logxform_cols) > 0) { + small = 1e-10 + sset = indat[,logxform_cols] + zeros = which(sset==0,arr.ind=T) + nz = nrow(zeros) + if (nz > 0) { + message(paste('Log transforming encountered',nz,'zeros - added 1e-10')) + sset[zeros] = sset[zeros] + small + lset = log(sset) + indat[,logxform_cols] = lset + } + } +} +if (force_xvar_cols_in != "None") +{ + force_xvar_cols = as.numeric(strsplit(force_xvar_cols_in,",")[[1]]) + allx = c(xvar_cols,force_xvar_cols) + xvar_cols = unique(allx) + xvar_cols = xvar_cols[order(xvar_cols)] + penalties[force_xvar_cols] = 0 +} +penalty = penalties[xvar_cols] +forcedin = paste(datcols[which(penalties == 0)],collapse=',') +cox_id_col = NA +cox_id = NA + +message(paste('@@@ Using alpha =',alpha,'for all models')) +x = indat[,xvar_cols] +nx = nrow(x) +cx = ncol(x) +message(paste('@@@@ Input has',nx,'samples and',cx,'predictors')) +if (cx > nx) { +message('@@@ WARNING: Models will have more variables than cases so glmnet will likely return one of many possible solutions! Please DO NOT expect reliable results - glmnet is clever but not magical @@@') +} + +xcolnames = datcols[xvar_cols] + +if (file.exists(Out_Dir) == F) dir.create(Out_Dir) +out_full = F +out_full_file = NA +out_pred_file = "" +out_pred = as.logical("$model.output_pred") + +#if $model.fam == "binomial" or $model.fam == "cox": + cox_id_col = $model.cox_id + cox_id = indat[,cox_id_col] + if (out_pred == T) { + out_pred_file="$output_pred_file" + rownames(x) = cox_id + } +#end if +#if $model.fam == "cox": + cox_time = $model.cox_time + cox_status = $model.cox_status + out_full = as.logical("$model.output_full") + if (out_full == T) { out_full_file="$output_full_file" } + yvar_cols = c(cox_time,cox_status) + ycolnames = c('time','status') + istat = as.double(indat[,cox_status]) + itime = as.double(indat[,cox_time]) + predict_at = quantile(itime) + if ("$model.predict_at" > "") + { + pa = "$model.predict_at" + predict_at = as.numeric(strsplit(pa,",")[[1]]) + } + y = data.frame(time = itime, status = istat) + ustat = unique(istat) + if ((length(ustat) != 2) | (! 1 %in% ustat ) | (! 0 %in% ustat)) + { + print.noquote(paste('INPUT ERROR: status must have 0 (censored) or 1 (event) but found',paste(ustat,collapse=',') )) + quit(save='no',status=1) + } + y = as.matrix(y) + x = as.matrix(x) + print.noquote(paste('@@@ Cox model will predict yvar=',datcols[cox_status],'using cols=',paste(xcolnames,collapse=','),'n preds=',length(xcolnames), + 'forced in=',forcedin)) + regres_out = dolasso_cox(x=x,y=y,debugOn=F,maxsteps=10000,nfold=nfold,xcolnames=xcolnames,ycolnames=ycolnames,optLambda=optLambda,out_full=out_full,out_full_file=out_full_file, + out_pred=out_pred,out_pred_file=out_pred_file,cox_id=cox_id,descr=myTitle,do_standard=do_standard,alpha=alpha,penalty=penalty,predict_at=predict_at,mdsplots=mdsplots) +#else: + yvar_cols = "$model.yvar_cols" + yvar_cols = as.numeric(strsplit(yvar_cols,",")[[1]]) + ycolnames = datcols[yvar_cols] + print.noquote(paste('@@@',fam,'model will predict yvar=',paste(ycolnames,collapse=','),'using cols=',paste(xcolnames,collapse=','),'n preds=',length(xcolnames), + 'forced in=',forcedin)) + y = data.matrix(indat[,yvar_cols]) + print.noquote(paste('Model will use',fam,'link function to predict yvar=',paste(ycolnames,collapse=','),'n preds=',length(xcolnames),'forced in=',forcedin)) + regres_out = dolasso_generic(predvars=x,depvars=y,debugOn=F, maxsteps=10000,nfold=nfold,xcolnames=xcolnames,ycolnames=ycolnames,optLambda=optLambda,out_pred_file=out_pred_file, + descr=myTitle,do_standard=do_standard,defaultFam=fam,alpha=alpha,penalty=penalty,out_pred=out_pred,cox_id=cox_id,mdsplots=mdsplots) +#end if + +write.table(format(regres_out,digits=5),outtab,quote=FALSE, sep="\t",row.names=F) +print.noquote('@@@ SessionInfo for this R session:') +sessionInfo() +warnings() + +</configfile> +</configfiles> + <inputs> + <param name="title" type="text" value="lasso test" label="Title for job outputs" help="Typing a short, meaningful text here will help remind you (and explain to others) what the outputs represent"> + <sanitizer invalid_char=""> + <valid initial="string.letters,string.digits"><add value="_" /> </valid> + </sanitizer> + </param> + <param name="input1" type="data" format="tabular" label="Select an input tabular text file from your history. Rows represent samples; Columns are measured phenotypes" + multiple='False' optional="False" help="Tabular text data with samples as rows, phenotypes as columns with a header row of column identifiers" /> + <param name="xvar_cols" label="Select columns containing numeric variables to use as predictor (x) variables" type="data_column" data_ref="input1" numerical="False" + multiple="True" use_header_names="True" force_select="True" /> + <param name="force_xvar_cols" label="Select numeric columns containing variables ALWAYS included as predictors in cross validation" type="data_column" data_ref="input1" numerical="False" + multiple="True" use_header_names="True" force_select="False"/> + <conditional name="model"> + <param name="fam" type="select" label="GLM Link function for models" + help="Binary dependant variables will automatically be set to Binomial no matter what this is set to"> + <option value="gaussian" selected="true">Gaussian - continuous dependent (y)</option> + <option value="binomial">Binomial dependent variables</option> + <option value="poisson">Poisson (eg counts)</option> + <option value="cox">Cox models - require special setup for y variables - see below</option> + </param> + <when value="gaussian"> + <param name="yvar_cols" label="Select numeric columns containing variables to use as the dependent (y) in elasticnet" type="data_column" data_ref="input1" numerical="False" + multiple="True" use_header_names="True" help = "If multiple, each will be modelled against all the x variables and reported separately." force_select="True"/> + <param name="output_full" type="hidden" value='F' /> + <param name="output_pred" type="hidden" value='F' /> + <param name="cox_id" label="Select column containing a unique sample identifier" + help = "Only really needed for output sample specific predicted values downstream." + type="data_column" data_ref="input1" numerical="False" force_select="True" + multiple="False" use_header_names="True" /> + </when> + <when value="binomial"> + <param name="yvar_cols" label="Select numeric columns containing variables to use as the dependent (y) in elasticnet" type="data_column" data_ref="input1" numerical="False" + multiple="True" use_header_names="True" help = "If multiple, each will be modelled against all the x variables and reported separately." force_select="True"/> + <param name="output_full" type="hidden" value='F' /> + <param name="output_pred" type="select" label="Create a tabular output with predicted values for each subject from the optimal model for (eg) NRI estimates" > + <option value="F" selected="true">No predicted value output file</option> + <option value="T">Create a predicted value output file</option> + </param> + <param name="cox_id" label="Select column containing a unique sample identifier" + help = "Only really needed for output sample specific predicted values downstream." + type="data_column" data_ref="input1" numerical="False" force_select="True" + multiple="False" use_header_names="True" /> + <param name="predict_at" type="hidden" value='' /> + + </when> + <when value="poisson"> + <param name="yvar_cols" label="Select columns containing variables to use as the dependent (y) in elasticnet" type="data_column" data_ref="input1" numerical="True" + multiple="True" use_header_names="True" help = "If multiple, each will be modelled against all the x variables and reported separately." force_select="True"/> + <param name="output_full" type="hidden" value='F' /> + <param name="output_pred" type="hidden" value='F' /> + <param name="predict_at" type="hidden" value='' /> + <param name="cox_id" label="Select column containing a unique sample identifier" + help = "Optional. Only really needed for output sample specific predicted values downstream. Free - enjoy" + type="data_column" data_ref="input1" numerical="True" force_select="False" + multiple="False" use_header_names="True" /> + </when> + <when value="cox"> + <param name="cox_time" label="Select column containing time under observation for Cox regression" + type="data_column" data_ref="input1" numerical="True" force_select="True" + multiple="False" use_header_names="True" help = "This MUST contain a time period - eg continuous years or days to failure or right censoring"/> + <param name="cox_status" label="Select column containing status = 1 for outcome of interest at the end of the time under observation or 0 for right censoring" + type="data_column" data_ref="input1" numerical="True" force_select="True" + multiple="False" use_header_names="True" help = "This MUST contain 1 for subjects who had an event at that time or 0 for a right censored observation"/> + <param name="cox_id" label="Select column containing a unique sample identifier" + help = "Optional. Only really needed for output sample specific predicted values downstream. Free - enjoy" + type="data_column" data_ref="input1" numerical="False" force_select="False" + multiple="False" use_header_names="True" /> + <param name="output_full" type="select" label="Create a tabular output with coefficients for all predictors" > + <option value="F" selected="true">No full model output file</option> + <option value="T">Create a full model output file</option> + </param> + <param name="output_pred" type="select" label="Create a tabular output with predicted values for each subject from the optimal model for (eg) NRI estimates" > + <option value="F" selected="true">No predicted value output file</option> + <option value="T">Create a predicted value output file</option> + </param> + <param name="predict_at" type="text" value='' label="Provide a comma separated list of times to make a prediction for each subject" + optional="True" help="Default (blank) will return predictions at 0%,25%,50%,75%,100% of the observed times which should be informative" /> + + </when> + </conditional> + <param name="optLambda" type="select" label="Value to use when reporting optimal model and coefficients" help="minLambda will have more predictors - 1SDLambda will be more parsimonious"> + <option value="lambda.1se" selected="true">Lambda + 1 SE of min MSE or AUC (fewer coefficients - more false negatives)</option> + <option value="lambda.min">Lambda at min MSE or max AUC (more coefficients - more false positives)</option> + </param> + <param name="logxform_cols" optional="True" label="Select numeric columns to be log transformed before use as predictors or dependent variables" type="data_column" + data_ref="input1" numerical="True" multiple="True" use_header_names="True" help = "The wisdom of doing this depends entirely on your predictors - eg can help diminish long-tailed outlier influence" + force_select="False"/> + <param name="do_standard" type="select" label="Standardise x vars" + help="If all measurements on same scale, may not be needed. Coefficients are always returned on the original scale."> + <option value="False" selected="true">No standardisation of predictors</option>l + <option value="True">Standardise predictors before model</option> + </param> + <param name="mdsplots" type="select" label="Generate MDS plots of samples in measurement space and measurements in sample space" > + <option value="False" selected="true">No MDS plots</option>l + <option value="True">Yes create MDS plots</option> + </param> + <param name="alpha" type="float" value="0.95" min="0.01" max="1.0" label="Alpha - see glmnet docs. 1 for pure lasso. 0.0 for pure ridge regression" + help="Default 0.95 allows lasso to cope better with expected predictor collinearity. Use (eg) 0.5 for hybrid regularised regression or (eg) 0.025 for ridge regression"/> + <param name="nfold" type="integer" value="10" label="Number of folds for internal cross validation" + help="Default of 10 is usually ok"/> + </inputs> + <outputs> + <data format="html" name="html_file" label="${title}.html"/> + <data format="tabular" name="model_file" label="${title}_modelres.xls"/> + <data format="tabular" name="output_full_file" label="${title}_full_cox_model.xls"> + <filter>model['output_full'] == 'T'</filter> + </data> + <data format="tabular" name="output_pred_file" label="${title}_predicted_from_model.xls"> + <filter>model['output_pred'] == 'T'</filter> + </data> + </outputs> + <tests> + <test> + <param name='input1' value='cox_test.xls' ftype='tabular' /> + <param name='treatment_name' value='case' /> + <param name='title' value='Cox glmnet test' /> + <param name='nfold' value='10' /> + <param name='logxform_cols' value='' /> + <param name='alpha' value='0.95' /> + <param name='do_standard' value="True" /> + <param name='cox_time' value='1' /> + <param name='cox_status' value='2' /> + <param name='cox_id' value='1' /> + <param name='predict_at' value='' /> + <param name='fam' value='cox' /> + <param name='yvar_cols' value='' /> + <param name='xvar_cols' value='3,4,5' /> + <param name='force_xvar_cols' value='3' /> + <param name='output_full' value='F' /> + <param name='output_pred' value='F' /> + <output name='model_file' file='coxlassotest_modelres.xls'> + <assert_contents> + <has_text text="rhubarb" /> + <has_text text="TRUE" /> + <!-- 	 is XML escape code for tab --> + <!-- has_line line="regulator	partial_likelihood	forced_in	glmnet_model	best_lambda" / --> + <has_line line="regulator	partial_likelihood	forced_in	glmnet_model	best_lambda	lambdaChoice	alpha" /> + <has_n_columns n="7" /> + </assert_contents> + </output> + <output name='html_file' file='coxlassotest.html' compare='diff' lines_diff='16' /> + </test> +</tests> +<help> + +**Before you start** + +Please read the glmnet documentation @ glmnet_ + +This Galaxy wrapper merely exposes that code and the glmnet_ documentation is essential reading +before getting useful results here. + +**What it does** + +From documentation at glmnet_ :: + + Glmnet is a package that fits a generalized linear model via penalized maximum likelihood. + The regularization path is computed for the lasso or elasticnet penalty at a grid of values for the regularization parameter lambda. + The algorithm is extremely fast, and can exploit sparsity in the input matrix x. + It fits linear, logistic and multinomial, poisson, and Cox regression models. + A variety of predictions can be made from the fitted models. + +Internal cross validation is used to optimise the choice of lambda based on CV AUC for logistic (binomial outcome) models, or CV mse for gaussian. + +**Warning about the tyrany of dimensionality** + +Yes, this package will select 'optimal' models even when you (optimistically) supply more predictors than you have cases. +The model returned is unlikely to represent the only informative regularisation path through your data - if you run repeatedly with +exactly the same settings, you will probably see many different models being selected. +This is not a software bug - the real problem is that you just don't have enough information in your data. + +Sufficiently big jobs will take a while (eg each lasso regression with 20k features on 1k samples takes about 2-3 minutes on our aged cluster) + +**Input** + +Assuming you have more measurements than samples, you supply data as a tabular text file where each row is a sample and columns +are variables. You specify which columns are dependent (predictors) and which are observations for each sample. Each of multiple +dependent variable columns will be run and reported independently. Predictors can be forced in to the model. + +**Output** + +For each selected dependent regression variable, a brief report of the model coefficients predicted at the +'optimal' nfold CV value of lambda. + +**Predicted event probabilities for Cox and Logistic models** + +If you want to compare (eg) two competing clinical predictions, there's a companion generic NRI tool +for predicted event probabilities. Estimates dozens of measures of improvement in prediction. Currently only works for identical id subjects +but can probably be extended to independent sample predictions. + +Given a model, we can generate a predicted p (for status 1) in binomial or cox frameworks so models can be evaluated in terms of NRI. +Of course, estimates are likely substantially inflated over 'real world' performance by being estimated from the same sample - but you probably +already knew that since you were smart enough to reach this far down into the on screen help. The author salutes you, intrepid reader! + +It may seem an odd thing to do, but we can predict p for an event for each subject from our original data, given a parsimonious model. Doing +this for two separate models (eg, forcing in an additional known explanatory measurement to the new model) allows comparison of the two models +predicted status for each subject, or the same model in independent populations to see how badly it does + +**Attributions** + +glmnet_ is the R package exposed by this Galaxy tool. + +Galaxy_ (that's what you are using right now!) for gluing everything together + +Otherwise, all code and documentation comprising this tool was written by Ross Lazarus and is +licensed to you under the LGPL_ like other rgenetics artefacts + +.. _LGPL: http://www.gnu.org/copyleft/lesser.html +.. _glmnet: http://web.stanford.edu/~hastie/glmnet/glmnet_alpha.html +.. _Galaxy: http://getgalaxy.org +</help> + +<citations> + <citation type="bibtex"> +@Article{Friedman2010, title = {Regularization Paths for Generalized Linear Models via Coordinate Descent}, + author = {Jerome Friedman and Trevor Hastie and Robert Tibshirani}, + journal = {Journal of Statistical Software}, + year = {2010}, + volume = {33}, + number = {1}, + pages = {1--22}, + url = {http://www.jstatsoft.org/v33/i01/} + } + </citation> + <citation type="doi"> +10.1093/bioinformatics/bts573 + </citation> +</citations> +</tool>
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/test-data/cox_test.xls Sat Oct 31 02:26:24 2015 -0400 @@ -0,0 +1,502 @@ +time status rhubarb vegemite apple +575.966708 0 452.405468 30.339584 32286.089057 +1539.245319 1 329.689929 30.603839 15735.863202 +2072.798422 0 534.379263 32.474983 22639.136685 +1638.450154 0 362.522161 30.925996 31108.364370 +1771.625630 0 417.167751 31.444652 31317.491931 +4413.484706 1 459.947526 30.069033 11232.037882 +1519.376431 1 576.868414 30.102709 16539.948172 +4805.620290 0 554.756164 30.685550 9979.125726 +657.497124 0 439.481951 32.305768 22301.057122 +181.498205 1 581.498091 32.549519 10499.737992 +1630.552909 0 580.670671 30.825057 22321.821037 +2285.043061 1 588.034529 32.799742 8067.710675 +3697.405336 0 400.195630 31.826099 20975.730587 +3821.226480 0 561.466202 32.333230 17887.547850 +4810.837196 0 447.004938 32.547717 25383.925946 +3195.877681 1 482.169788 31.420740 27659.409774 +4556.750833 1 415.144072 31.960725 18401.138992 +3797.192061 1 484.492359 30.789882 7431.214452 +3734.484446 1 587.339066 32.052766 29443.923798 +567.665985 1 377.599777 31.203617 30148.076752 +2551.367122 0 447.623689 31.795708 15922.972458 +3.724234 1 534.080187 31.175131 8259.786167 +3528.906130 0 471.293653 31.225006 18307.543198 +3649.459162 1 443.804109 32.390564 9720.974163 +1490.640766 1 456.837379 32.099378 15292.782705 +4023.629633 1 443.869307 30.022985 23145.513801 +1389.288557 0 497.850711 30.361120 28757.227200 +535.427950 1 437.634304 32.307966 18289.733136 +3115.819953 0 587.593182 31.317353 21533.733397 +1561.895566 1 472.250365 31.550546 8724.067134 +2579.785003 0 331.955830 32.992176 12105.988070 +1592.819387 0 361.242026 31.011367 12823.339723 +4473.715814 0 334.539193 30.448049 9621.600487 +4328.198133 0 537.782739 32.683122 19090.292579 +3761.244687 1 392.916387 32.665889 22119.375302 +4638.730675 0 514.168146 30.548414 20651.266975 +1531.289211 0 451.404672 32.437627 21748.461852 +3844.140593 0 427.404840 32.894070 25824.198908 +4542.766339 0 366.992470 30.814372 6090.765585 +3296.236419 0 575.910564 32.509193 4050.864668 +710.915095 0 395.062092 32.940952 9086.464980 +3211.769396 0 428.608771 30.918481 29996.819084 +4448.442217 0 300.153203 31.160941 8735.458385 +4215.329607 0 504.064526 32.040819 29036.000397 +3249.725224 0 489.538009 32.863547 27520.919952 +3388.033079 0 485.938597 30.920035 20395.348731 +1465.239279 0 435.641977 30.193093 20469.748911 +2022.519649 0 362.011214 32.756787 18367.094394 +1474.662644 0 448.727985 31.194883 18456.157990 +2680.178585 0 553.521680 32.339108 17386.094437 +254.953131 1 467.930054 32.705913 28823.817908 +2080.749522 0 387.588403 30.520279 3389.418182 +4609.398123 1 471.568096 31.101077 27574.777543 +2506.254687 1 466.233702 32.928613 24693.847389 +1418.686322 0 519.402535 32.685931 9055.300824 +1716.501599 0 305.671024 31.453871 25382.505028 +2080.668376 0 589.514396 31.143248 24546.670311 +4051.441999 1 431.028785 30.854820 32291.736106 +221.739141 1 435.123612 30.465214 19077.175842 +3945.564918 0 559.566964 32.871624 18320.616591 +3543.190633 1 548.960006 30.269711 6259.113819 +4141.263169 1 418.013926 30.159280 12806.181141 +4613.643735 1 467.593153 30.584036 28998.377681 +2526.696757 1 442.602416 31.075723 24392.589804 +661.056610 0 306.186743 30.071900 4931.492293 +420.059232 0 557.043309 30.889234 30410.871000 +3777.621646 1 564.397202 30.810439 25359.485010 +4210.010862 0 541.943504 31.109403 18240.759045 +2251.752128 1 442.108369 31.308836 14136.141189 +4037.391365 1 467.892619 32.089289 14540.468274 +2478.012769 0 305.600832 30.653583 7117.580065 +2435.693015 0 551.918277 30.355254 29969.743529 +3816.690227 0 314.516061 31.511033 9880.067119 +531.400047 0 565.141705 30.722177 11252.478408 +4367.861436 0 421.500739 31.696812 16879.486148 +634.857316 1 476.587731 30.741441 31351.695315 +2757.953233 1 484.171370 31.510266 7767.961089 +676.813083 0 412.948756 30.174301 11123.425270 +3672.369695 0 374.771725 30.531392 32349.940422 +2572.842389 0 554.137199 32.777709 12727.962812 +4960.699040 0 560.164589 30.206776 3453.492367 +4316.475328 1 484.168551 32.069106 31984.273406 +4015.881070 1 400.274241 30.658886 11783.352436 +3052.451190 1 398.216589 30.358316 27928.803920 +1398.823206 0 446.521178 30.185715 31193.500034 +3922.909840 1 575.609185 32.724842 10233.363838 +182.803149 0 370.875903 30.856525 23954.116536 +2632.176915 1 426.948756 30.641145 30860.455084 +1060.390433 0 300.347626 31.143117 23615.898639 +3178.728283 1 514.261332 31.356224 26944.147061 +763.218442 1 463.367562 31.274306 7509.665226 +953.708474 0 384.301310 31.065789 15102.331493 +4434.527995 1 431.967331 32.369140 7177.785265 +623.914857 0 495.381873 32.915097 22739.985454 +3157.386533 1 515.293513 32.210421 31073.191185 +521.293142 0 593.449301 30.668076 29606.157493 +2845.587324 1 367.180408 32.424788 7052.350189 +4832.119270 1 465.679571 30.019029 5249.534096 +2856.445302 1 591.090965 31.133381 32498.377584 +2572.527257 1 409.021421 31.885768 13743.016360 +1413.194580 1 330.844512 32.740409 5441.606816 +2918.271731 0 304.169439 31.256918 22384.990257 +1674.024061 1 589.100300 32.763039 13305.351548 +2808.138084 0 491.357574 31.150814 17202.472869 +1424.663181 0 359.990029 31.343898 28564.024266 +2718.446177 1 580.778432 31.842978 29231.311831 +3735.730841 0 573.549637 30.079704 31952.172593 +4641.585855 0 365.172295 32.894223 17289.384320 +4071.637795 0 505.251583 32.472570 9633.518042 +471.788564 1 322.257982 32.546866 30076.426907 +1266.346159 1 335.749836 31.646670 32477.736291 +2695.370600 1 509.475837 31.947731 22351.940960 +3631.475811 0 584.575140 30.366829 13700.704821 +488.174184 0 522.383785 30.637536 29035.941989 +2645.805895 1 390.848623 31.639290 3004.701971 +4037.703358 1 425.444485 32.775358 28720.752294 +2196.610999 1 547.619331 32.640598 29108.134173 +1081.264261 1 350.044271 32.021783 19844.288986 +1074.305164 1 497.881588 32.048148 16325.089356 +1643.851420 0 372.156064 30.815571 12924.845649 +4919.287039 0 470.345027 31.043515 6953.198088 +3961.891671 1 353.750322 32.386195 10074.185995 +1643.427847 0 521.758887 31.784944 10303.308946 +2714.177818 0 363.519932 31.954710 4089.970927 +1958.737870 0 372.752456 32.778447 25735.088968 +72.430537 1 374.015311 30.406160 3519.910494 +666.824775 0 443.814245 31.382885 32222.429084 +179.203304 0 339.335011 30.911475 29294.393947 +3193.179459 1 329.239965 32.270432 7539.157039 +1405.955780 1 575.062041 31.769778 26056.208652 +4378.622335 1 583.779711 30.280724 12958.062814 +4933.985803 0 474.998284 30.835755 8298.116085 +479.970673 0 311.469610 31.110600 16498.771798 +194.917302 1 533.735116 32.682195 17363.293404 +3126.230367 1 568.043474 32.209623 22525.398240 +3605.775657 0 519.068951 30.750083 21588.128190 +68.236883 1 316.956971 32.735304 5303.699685 +2956.632516 1 359.835070 31.735527 26701.676218 +1273.924283 1 346.804469 31.185254 29977.854210 +1530.707021 1 568.442488 31.510433 20716.251168 +3847.696838 1 487.304670 32.242322 5276.023316 +3121.510941 0 518.612631 30.448810 10702.449428 +2693.401144 1 533.008711 32.964097 31390.744850 +1774.742268 1 588.145493 30.487710 21800.531169 +2216.272824 0 467.093275 31.533770 7991.454914 +1526.314833 0 433.703579 32.809884 11209.353005 +2089.304275 0 321.084472 31.021678 27856.425967 +3046.990363 1 460.307711 32.631130 7182.074497 +1317.038946 0 593.108863 31.275244 26560.416891 +1352.101840 0 306.529955 30.979817 21571.506540 +3163.866962 0 304.262526 30.303302 24717.381911 +3417.487105 1 367.150160 31.345893 13987.483107 +4582.950578 0 478.091198 32.280854 23336.546777 +1183.214096 1 435.761798 31.182880 32631.545653 +549.532217 0 436.107384 32.534391 23409.562965 +4018.079504 1 514.150908 32.582211 11007.862101 +4537.862661 0 474.720396 31.380999 6160.620917 +2023.882772 1 576.659734 31.873967 24564.400593 +301.339657 0 525.192309 30.612934 19844.142010 +989.521636 0 323.610331 30.471812 6456.387005 +3647.414661 1 353.160773 30.687803 20893.860619 +184.390450 0 367.556798 32.894773 7004.911101 +2851.709060 1 420.820133 32.238140 4715.183972 +3988.318894 0 352.860472 31.705620 18426.123706 +2244.796525 0 332.927674 31.617697 20768.384056 +726.779082 0 370.641547 32.668547 5629.086601 +2379.683675 1 566.970531 30.287987 26943.017249 +2039.968848 0 550.042019 31.203815 6686.226792 +949.545757 1 302.786734 30.845201 7393.079915 +2509.519553 1 378.264518 31.980310 16762.542739 +2949.998835 0 472.701800 32.896473 12708.457939 +3958.575946 1 402.750276 30.318068 10667.350272 +550.167142 0 394.589505 30.431805 32262.087079 +1462.278890 1 534.727448 31.514861 27759.667272 +2696.454522 1 535.656182 32.096399 32431.581354 +3756.029506 1 332.505291 32.486460 26887.712404 +1948.023625 1 470.099353 31.100519 31437.137850 +2008.758594 1 489.846790 31.656842 18513.611970 +889.052994 0 378.045215 30.360086 30999.351259 +956.860624 0 498.874188 31.238091 32986.760129 +1596.180150 1 468.477647 32.876509 10708.862352 +2044.623310 0 341.364451 30.831209 4377.873160 +1670.382802 0 344.655922 31.528413 25861.700033 +2751.419144 0 525.374048 31.427506 23893.998314 +1139.794772 0 510.822662 32.192052 27853.949410 +1951.275598 1 355.760660 30.362940 26312.111688 +3652.432308 0 330.985338 32.585428 25167.952966 +500.427055 1 453.806271 30.774837 16661.345691 +2884.409660 0 529.815414 31.683479 3317.970560 +2365.322179 0 540.382433 31.999929 29498.364669 +3853.672492 1 367.171163 32.671321 7512.370365 +938.719173 0 569.038125 31.387703 23444.121607 +3058.327259 1 411.675057 31.700811 28738.835254 +4215.442937 0 425.942827 31.617841 19505.555179 +1639.143130 1 310.802062 31.996376 24378.789768 +3015.816015 1 397.865213 30.774914 32208.055822 +1361.137094 1 314.093792 30.698064 21667.130306 +1900.651819 1 582.449408 30.643572 26755.743826 +2038.834435 1 359.535511 31.758170 26389.582514 +4217.395801 1 532.389233 32.367523 28128.901119 +4048.824843 1 415.408543 30.661216 6695.688224 +3278.940683 0 590.282168 31.148497 17706.150830 +2590.370816 1 467.475790 32.894899 31984.209283 +2961.953200 1 409.350826 30.412768 26800.587858 +757.305077 0 531.036597 32.692111 23056.009036 +4823.679979 0 547.031467 31.700254 4679.396986 +2159.887871 1 461.722894 31.554845 23118.215328 +1990.158592 0 436.583574 32.276433 22516.752446 +4010.528172 0 315.169873 31.546055 4521.472077 +527.527936 0 429.812828 32.161079 26933.467576 +2517.395024 1 442.643281 32.903345 11566.976782 +3973.530406 0 444.298671 30.860264 30738.467879 +909.115735 1 548.278268 31.969666 3061.067051 +1003.138129 0 599.492409 32.731826 5268.191640 +496.214281 1 525.539783 32.288968 18314.877509 +3321.884314 1 518.880235 31.593901 27585.403619 +531.060305 0 311.564795 30.660864 26630.390327 +1985.865186 0 468.954849 30.022861 7126.223629 +1075.080978 0 566.983564 32.548614 4189.576013 +2899.914447 0 501.981481 31.522821 3699.352956 +1209.686792 1 397.096875 30.589148 26885.709293 +1061.484419 0 399.802412 31.376845 21631.836258 +2341.148192 0 599.041273 32.705976 7273.381892 +2350.754927 1 578.118121 32.505017 16776.353449 +3715.290307 1 374.328239 31.642496 24018.734911 +353.266473 1 306.976330 30.935877 17228.695180 +739.123546 1 529.022431 30.659845 18290.776722 +465.537407 0 308.998940 31.824504 3898.463033 +3524.840397 1 385.275512 30.348328 6065.321632 +751.616003 0 322.281940 30.576017 20828.524508 +4989.405682 1 381.074571 30.437314 16844.697669 +488.563077 0 318.125000 31.138994 28854.306278 +4558.130890 0 479.058472 31.076337 3808.209751 +3644.754920 0 437.984557 30.038643 24302.584238 +4595.447449 0 323.287824 31.904699 24911.538572 +1064.168831 1 541.625023 30.168168 12380.733046 +4197.068049 0 418.764938 31.133507 14481.902057 +3593.303853 0 545.685419 30.904012 27860.611155 +3842.381434 0 596.665096 31.650039 32071.065698 +781.763222 0 544.786508 32.895055 4905.803177 +3941.386557 1 362.501808 32.489031 8386.976863 +4285.917626 1 373.105515 31.909823 18733.487831 +4152.435642 1 440.036083 30.439913 4669.690889 +2326.656610 0 409.061810 31.219523 17268.378337 +887.243489 1 458.991888 32.904742 32560.558550 +2383.154926 0 402.990252 31.501777 32165.195285 +1575.020117 1 567.576879 31.223480 11851.574021 +2681.793474 0 382.275765 30.672297 31505.220294 +3745.755933 1 524.831384 31.755858 20277.860327 +2310.372497 0 381.887619 32.549019 23884.208743 +1934.629896 0 430.992059 30.090286 16064.349955 +228.815046 1 510.344577 30.272399 13314.791074 +106.494754 0 421.895617 30.249703 6022.257615 +107.758807 1 352.048161 31.015999 27010.353099 +2366.497161 1 589.257754 32.599191 29060.870687 +3652.314844 1 325.392214 30.378902 9562.454775 +4628.382437 0 389.680904 30.488690 21333.408550 +3992.522719 1 460.945057 30.259246 21474.935828 +3116.879118 0 470.551054 30.384387 23511.934422 +3669.891947 1 468.454570 30.792238 6696.682439 +2455.615625 1 536.995869 30.130133 27385.537567 +334.254171 0 510.509005 31.195189 15166.910603 +2392.466005 1 503.464732 31.609503 28757.028622 +563.741633 0 535.281334 31.242527 11455.193241 +2683.204132 0 530.400221 30.122120 21995.611523 +116.891969 0 458.723324 30.381523 5852.247596 +2817.127744 1 547.927894 32.877328 22278.789600 +323.704077 0 432.630710 31.643437 8793.259513 +2223.194666 1 584.593612 31.798636 4232.455807 +4836.806252 0 328.278741 31.551587 28407.563886 +2302.391657 1 304.186017 30.683099 21760.214902 +1597.778913 1 593.734794 32.275770 26311.401980 +3231.228313 0 469.955959 31.155902 21789.637085 +4329.078171 1 326.299894 30.629070 28598.928501 +3004.531335 0 345.517768 31.874809 20075.187439 +3426.658667 0 407.558711 32.882949 32961.004750 +2942.293151 1 474.436479 30.681717 10297.918976 +2790.908728 1 372.634696 31.883438 18617.426989 +3645.736003 0 338.768900 30.271779 13340.959790 +4954.230936 0 475.204218 32.305344 31368.726270 +1374.169666 0 322.794825 32.054615 29507.251607 +952.311143 0 466.066676 31.283949 31527.207980 +1582.935333 1 430.497522 32.636945 6365.677867 +1574.488078 1 512.547475 32.000036 3214.800133 +3675.314941 0 420.357568 30.365188 24034.010725 +3737.900636 0 310.463467 32.944674 19212.098849 +453.453905 1 307.998144 31.734006 6696.312537 +180.430188 0 500.121289 31.647589 25439.918801 +2223.386978 0 342.168201 31.413275 10321.858985 +3347.152707 0 450.255234 32.419409 21979.139997 +2032.988552 1 359.268392 32.444386 4151.429949 +1855.874006 1 486.902129 32.194183 28734.066535 +576.472090 0 491.131111 30.250436 21183.624322 +2629.739140 1 490.880496 30.706030 20741.966223 +3234.449151 0 591.702238 31.565719 9467.501732 +4568.779861 0 391.821980 31.778473 12171.711242 +1369.135728 1 352.424832 32.501343 24914.054783 +1203.818730 1 453.408627 30.021814 9344.375802 +1328.193901 1 326.714975 31.957412 22506.010391 +3143.392191 0 488.722719 30.906428 20404.329854 +1516.300551 1 524.551137 30.783381 26784.582992 +4401.773609 1 490.129039 31.474209 27673.989380 +471.342147 0 345.653763 30.714091 27614.068357 +2846.556659 0 387.773914 32.101213 25758.358083 +3113.332960 0 340.017571 31.250862 8517.797718 +569.863929 1 507.285752 32.723247 7350.729977 +1351.433366 1 355.255156 31.629009 10035.728199 +2057.479962 0 504.044203 31.725231 14699.764489 +2029.207221 1 326.501058 31.883646 25956.067323 +4987.734533 0 391.140321 30.918262 15723.785267 +2728.957211 1 498.936107 32.017269 21531.521853 +2885.926679 0 458.513681 31.221876 23111.352960 +1214.796671 0 531.593343 32.409879 26260.051485 +4073.228494 1 563.045005 31.629048 5742.795393 +632.393708 1 428.087852 32.360714 20236.813182 +4548.671698 0 387.645357 30.048297 28095.140888 +735.589017 0 510.876342 30.202972 28887.002260 +2683.928468 1 352.751328 31.385745 3968.105767 +3761.692766 1 470.650685 31.713693 21645.943688 +3949.678238 0 572.913454 32.230970 18919.434641 +1695.346973 0 387.076598 32.045632 8265.498458 +4805.534106 0 322.040787 30.256378 20922.671974 +3834.130476 0 548.536590 30.581470 12530.815196 +3700.095266 0 434.154537 30.131024 28914.087338 +2083.132166 1 435.468758 31.961937 24155.150686 +2566.436697 1 455.368285 30.113525 23903.526223 +3372.416477 1 439.390127 32.380970 8312.312507 +3341.243804 1 496.573292 32.055384 9484.908468 +1289.638171 1 346.709851 31.116968 13853.905474 +441.798476 1 395.115019 32.820316 27178.509065 +3876.357092 1 501.568592 30.100469 3872.161975 +3082.204963 0 445.140351 31.403349 21141.060076 +1222.633709 1 510.019589 32.832816 12522.872061 +3257.092925 1 567.760441 30.935411 17780.698060 +4978.318869 0 488.434181 31.349966 20399.317688 +268.617967 0 563.410236 30.677835 13353.114481 +3668.370880 1 432.604476 31.336602 18298.862821 +1616.920156 1 528.362539 30.664065 18290.365790 +1582.470140 1 566.928554 31.829775 31526.762426 +1894.952274 0 576.611659 32.930078 5451.111476 +130.599064 1 357.672426 30.192175 5127.046275 +1658.869611 0 518.418202 30.222626 21787.150970 +1930.327624 0 539.744396 31.409771 11998.451350 +1360.478477 1 357.741711 32.358872 4894.722814 +2762.405397 1 498.202228 32.110196 20926.903902 +66.089100 0 525.721147 30.742496 28773.910408 +4617.260796 0 388.807122 31.052491 19010.672351 +4046.674710 1 363.050835 32.600458 8208.205191 +103.101768 0 401.874765 30.885852 29276.307064 +2870.170146 1 479.961929 31.997551 19841.555525 +1069.648823 1 598.655109 31.936390 10655.006702 +1988.650561 0 483.709224 31.937744 11014.379454 +2983.772162 0 328.380916 32.023817 8974.076899 +4534.002381 0 311.260565 30.136494 25947.551507 +3064.855784 1 369.080218 31.164070 29129.542680 +1746.878916 0 313.512424 30.992758 20731.078669 +4018.774695 1 497.866125 31.834090 3971.821700 +3742.871097 1 353.541863 30.115952 20815.885581 +2828.315773 1 454.655490 32.424599 26891.927439 +637.033003 1 424.553960 31.428681 20585.183759 +3087.856465 1 536.195443 31.597533 12661.823700 +604.873967 1 378.515493 30.925004 22439.100233 +3748.531033 1 515.115463 32.625710 6408.313067 +1558.535013 1 433.937532 32.380082 8383.657584 +2365.992186 0 409.216710 31.115409 32511.614098 +1183.022513 0 362.544664 30.959448 5065.906749 +3789.375438 0 595.724663 30.914010 20244.524890 +850.273991 1 387.136972 31.618840 5094.567581 +1231.544018 1 322.683917 31.457545 20767.599384 +1581.181089 1 305.296620 31.856051 17539.248499 +4880.147078 1 316.598007 31.317619 11086.097072 +4013.027211 0 485.648876 31.953258 30138.577997 +31.538821 0 486.966298 30.424936 3940.149283 +1735.735564 1 453.551097 30.972695 29273.785203 +4682.372724 1 426.064369 30.915967 29848.714453 +1419.673303 1 410.849400 31.957769 20423.937410 +4735.693777 0 331.994083 30.314568 21510.813276 +796.562868 0 502.136561 31.584336 14765.778268 +2637.389588 1 550.253910 31.759164 29206.835372 +1596.336472 0 367.679269 32.352327 5586.166324 +4673.788569 0 380.575280 32.897426 12723.911545 +2339.522676 0 405.444884 30.501805 20468.692308 +2935.072522 0 458.797029 30.585909 9955.177397 +2865.796979 0 497.765242 31.556959 27770.325805 +1032.502528 1 390.204545 31.050812 26330.543896 +1571.491134 1 454.897576 32.799730 4600.369582 +453.244362 0 513.834562 31.693168 18448.030313 +2566.490430 1 370.231825 31.371622 3530.831525 +2470.800433 1 506.755104 32.707049 15365.573438 +1674.578908 1 484.850818 31.657126 31034.118967 +1569.707807 0 380.339060 31.625216 26765.694375 +2048.318765 1 502.814703 32.061583 26794.094870 +4878.355942 1 308.718223 31.785469 9019.516491 +3394.948104 0 354.869800 31.553766 18461.635011 +3784.598905 0 506.145921 32.828899 22372.102070 +3648.947834 0 503.931222 32.468882 27102.757907 +1722.690965 1 389.086536 31.876640 5309.945334 +1266.660624 0 452.533007 30.468222 28572.903114 +794.701939 1 459.408945 32.075307 21311.040821 +1461.841041 0 370.527899 32.775555 10211.935095 +1976.852418 0 492.898069 30.817372 17854.725954 +859.962977 1 304.267257 30.348074 3439.856318 +2900.907710 1 396.919281 31.272958 4858.314485 +3693.615130 1 597.052502 30.903533 25716.531228 +1281.197465 1 470.854828 30.164368 9315.105349 +126.181622 1 381.237663 30.660067 28411.425041 +2517.244075 0 559.007071 30.373940 11135.538898 +2069.720704 1 544.708646 31.610539 16441.322232 +578.677623 1 591.867328 32.146388 21419.935158 +1918.809273 1 308.991680 31.457030 31410.129144 +1074.094412 0 336.860228 30.308287 17215.627286 +3034.994121 0 560.991229 31.862778 17867.418001 +4794.119298 1 480.292221 30.685158 11676.367480 +3307.911419 0 389.291016 30.558138 31108.361766 +457.891802 1 588.926068 32.330020 21119.054157 +1258.403059 0 332.016717 31.309245 16741.712001 +3312.429446 0 568.695418 31.010073 14155.107714 +3839.897879 0 370.810835 31.722406 16960.037012 +2068.214455 0 392.107828 30.058746 32853.457910 +2526.811696 0 510.287902 30.913056 22333.505644 +914.090787 1 323.087744 32.697177 14094.308393 +1890.199376 0 597.011219 31.259415 19326.502842 +1234.225123 0 536.105198 30.519191 9869.430841 +350.596239 1 320.162343 30.373749 5065.307636 +2846.841176 0 403.168597 31.624641 25354.766351 +1853.530932 1 564.965357 30.095729 11551.597590 +1183.876207 1 572.411649 32.013383 17167.971369 +3953.600477 0 404.253211 30.364392 17515.400202 +3289.696058 1 576.771217 31.015263 14654.467190 +2975.113909 1 587.011738 31.223106 7437.225919 +1404.141852 0 372.793625 30.325109 13500.525950 +1634.182114 0 569.705673 32.338410 27233.776345 +462.646125 1 467.892571 30.013450 27556.195432 +4138.443866 0 318.632678 30.348871 21795.067964 +4139.367443 1 304.173328 32.473154 29126.297266 +969.741373 1 599.886523 32.409973 16299.878199 +1494.463287 1 429.995117 30.484126 17461.996880 +4336.770459 1 416.241963 30.152207 32032.419174 +3064.782012 1 428.313584 31.688383 31731.608033 +2581.013362 0 490.480997 32.990983 17803.475008 +413.230864 1 430.932119 30.358681 18338.089371 +4469.602591 0 544.949864 30.599950 17579.676286 +2299.774314 1 532.263729 32.087371 25491.390739 +3570.957635 0 429.164688 30.798996 21783.084247 +2730.725224 0 434.512962 32.353024 8584.173078 +1859.613968 0 465.494240 30.853468 24199.875583 +3239.212454 1 585.752947 31.392592 17429.754381 +1626.243010 0 548.642457 32.586064 28320.002213 +4552.584886 0 411.908494 32.541147 18556.711124 +719.687465 1 437.155107 32.439905 18592.734838 +4528.629292 1 598.029535 30.261328 28073.717443 +1445.996134 1 407.243137 30.545202 28549.910302 +2143.373187 1 581.134422 30.728743 26291.822745 +4543.013656 0 406.670256 30.716994 24689.785142 +3479.794580 0 585.305633 31.528087 18146.824396 +2206.562888 0 526.292818 32.250928 6586.281804 +1796.678722 0 441.391846 32.716749 16294.189483 +2705.347154 0 551.894682 31.578422 11390.827042 +172.103511 1 320.919515 31.863888 29876.119334 +315.943273 1 512.197902 32.302445 17856.311999 +2750.695088 1 565.301177 31.514433 19939.221077 +4398.626158 0 418.977196 32.337330 4024.240756 +1522.050792 0 554.491924 32.249187 30588.435646 +3385.739992 1 496.030916 31.480018 17815.393181 +955.825566 1 393.353810 31.434487 26026.721071 +3971.619403 1 552.518065 30.319670 5068.002389 +2010.848206 1 569.817785 31.925723 26117.233780 +2181.752401 1 535.925058 30.952605 3509.328188 +905.194864 1 588.512406 30.229715 15643.423750 +3274.240917 0 590.044688 32.642379 13352.913630 +2061.196669 0 457.883033 30.745459 19655.487967 +4056.062516 0 354.337157 30.418176 24440.653038 +611.999600 0 432.598681 30.201091 19387.208118 +45.061679 1 437.997979 31.939322 21882.265038 +2091.842258 1 445.122275 32.963802 4496.459612 +2785.821696 1 535.953261 32.295237 5166.117813 +1272.363761 0 498.289820 31.101707 13536.602174 +2161.870929 0 487.379529 30.927239 22791.304718 +1230.121506 0 326.397921 31.953971 21880.979367 +2065.316190 0 405.291555 30.660933 28025.816633 +3230.173916 1 365.768097 30.055539 23327.051800 +4572.451029 1 323.399890 31.386897 21097.082477 +4916.172714 0 384.428209 32.754127 10981.900793 +4882.524497 0 399.434058 30.795923 18117.477126 +1457.322605 1 412.012101 30.834341 20230.063544 +1988.228627 0 566.906278 31.371606 15413.227181 +4337.262787 1 316.056609 32.556132 11299.235579 +2762.610907 1 454.177332 32.779747 14436.150973 +4756.268480 1 562.493863 32.727774 26209.569931 +3031.941096 0 313.505196 31.386484 9197.983181 +4586.200757 0 329.983728 30.115395 27039.250998 +300.892434 0 520.668472 32.713615 14579.471742 +392.962829 0 513.293270 30.568538 10311.182381 +2116.838847 1 355.318371 30.883049 18887.246656 +779.846218 1 407.952315 30.145534 22063.619725 +4236.655978 0 326.344600 32.797447 16004.751246 +4008.051049 1 319.382840 30.262423 20092.936132 +1550.891721 1 570.189663 32.238898 11686.528747 +2616.214276 1 448.726264 31.633295 26368.520969 +2775.847491 1 532.436984 30.938190 21751.684582 +
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/test-data/coxlassotest.html Sat Oct 31 02:26:24 2015 -0400 @@ -0,0 +1,167 @@ +<!DOCTYPE html PUBLIC "-//W3C//DTD XHTML 1.0 Transitional//EN" "http://www.w3.org/TR/xhtml1/DTD/xhtml1-transitional.dtd"> + <html xmlns="http://www.w3.org/1999/xhtml" xml:lang="en" lang="en"> + <head> <meta http-equiv="Content-Type" content="text/html; charset=utf-8" /> + <meta name="generator" content="Galaxy rgToolFactory.py tool output - see http://g2.trac.bx.psu.edu/" /> + <title></title> + <link rel="stylesheet" href="/static/style/base.css" type="text/css" /> + </head> + <body> + <div class="toolFormBody"> + +<div class="infomessage">Galaxy Tool "rglasso" run at 18/02/2015 22:06:08</div><br/> +<div class="toolFormTitle">cox images and outputs</div> +(Click on a thumbnail image to download the corresponding original PDF image)<br/> +<div><table class="simple" cellpadding="2" cellspacing="2"> +<tr> +<td><a href="cox_Coxglmnettest_glmnet_cvdeviance.pdf"><img src="cox_Coxglmnettest_glmnet_cvdeviance.png" title="Click to download a PDF of cox_Coxglmnettest_glmnet_cvdeviance.pdf" hspace="5" width="400" + alt="Image called cox_Coxglmnettest_glmnet_cvdeviance.pdf"/></a></td> + +<td><a href="cox_Coxglmnettest_glmnetdev.pdf"><img src="cox_Coxglmnettest_glmnetdev.png" title="Click to download a PDF of cox_Coxglmnettest_glmnetdev.pdf" hspace="5" width="400" + alt="Image called cox_Coxglmnettest_glmnetdev.pdf"/></a></td> +</tr> + +</table></div> + +<div class="toolFormTitle">cox log output</div> + +<pre> + +[1] @@@ Cox model will be predicted at times = 3.724234,1325.40516225,2373.090418,3650.1730825,4989.405682 + +[1] Times each predictor was selected in CV models (excluding zero count predictors) + + hits_lambda_1se hits_lambda_min + +rhubarb 10 10 + +vegemite 6 6 + +[1] @@@ Results preview: + + regulator partial_likelihood forced_in glmnet_model best_lambda lambdaChoice alpha + +rhubarb rhubarb 0.0012215 TRUE cox 0.04197 lambda.1se 0.95 + + +</pre> + +<div class="toolFormTitle">rglasso log output</div> + +<pre> + +Loading required package: Matrix + +Loading required package: methods + +Loaded glmnet 1.9-8 + +Loaded lars 1.2 + +Loading required package: splines + +Warning messages: + +1: In if (is.na(predict_at)) { : + + the condition has length > 1 and only the first element will be used + +2: In if (class(p) == "try-error") { : + + the condition has length > 1 and only the first element will be used + +3: In plot.window(...) : "label" is not a graphical parameter + +4: In plot.xy(xy, type, ...) : "label" is not a graphical parameter + +5: In axis(side = side, at = at, labels = labels, ...) : + + "label" is not a graphical parameter + +6: In axis(side = side, at = at, labels = labels, ...) : + + "label" is not a graphical parameter + +7: In box(...) : "label" is not a graphical parameter + +8: In title(...) : "label" is not a graphical parameter + + +</pre> + +<div class="toolFormTitle">Other log output</div> + +<pre> + +## Toolfactory running rglasso as Rscript script + +[1] @@@ Using alpha = 0.95 for all models + +[1] @@@@ Input has 500 samples and 3 predictors + +[1] @@@ Cox model will predict yvar= status using cols= rhubarb,vegemite,apple n preds= 3 forced in= rhubarb + +[1] @@@ SessionInfo for this R session: + +R version 3.1.0 (2014-04-10) + +Platform: x86_64-unknown-linux-gnu (64-bit) + +locale: + + [1] LC_CTYPE=en_AU.UTF-8 LC_NUMERIC=C LC_TIME=en_AU.UTF-8 LC_COLLATE=en_AU.UTF-8 LC_MONETARY=en_AU.UTF-8 LC_MESSAGES=en_AU.UTF-8 LC_PAPER=en_AU.UTF-8 LC_NAME=C LC_ADDRESS=C LC_TELEPHONE=C LC_MEASUREMENT=en_AU.UTF-8 LC_IDENTIFICATION=C + +attached base packages: + +[1] splines methods stats graphics grDevices utils datasets base + +other attached packages: + +[1] pec_2.4.4 survival_2.37-7 lars_1.2 glmnet_1.9-8 Matrix_1.1-5 + +loaded via a namespace (and not attached): + +[1] codetools_0.2-10 foreach_1.4.2 grid_3.1.0 iterators_1.0.7 lattice_0.20-29 lava_1.3 prodlim_1.5.1 + +Warning messages: + +1: In if (is.na(predict_at)) { ... : + + the condition has length > 1 and only the first element will be used + +2: In if (class(p) == "try-error") { ... : + + the condition has length > 1 and only the first element will be used + +3: In plot.window(...) : "label" is not a graphical parameter + +4: In plot.xy(xy, type, ...) : "label" is not a graphical parameter + +5: In axis(side = side, at = at, labels = labels, ...) : + + "label" is not a graphical parameter + +6: In axis(side = side, at = at, labels = labels, ...) : + + "label" is not a graphical parameter + +7: In box(...) : "label" is not a graphical parameter + +8: In title(...) : "label" is not a graphical parameter + + +</pre> + +<div class="toolFormTitle">All output files available for downloading</div> + +<div><table class="colored" cellpadding="3" cellspacing="3"><tr><th>Output File Name (click to view)</th><th>Size</th></tr> + +<tr><td><a href="cox_Coxglmnettest_glmnet_cvdeviance.pdf">cox_Coxglmnettest_glmnet_cvdeviance.pdf</a></td><td>6.4 KB</td></tr> +<tr class="odd_row"><td><a href="cox_Coxglmnettest_glmnetdev.pdf">cox_Coxglmnettest_glmnetdev.pdf</a></td><td>5.3 KB</td></tr> +<tr><td><a href="cox_cross_validation_model_counts.xls">cox_cross_validation_model_counts.xls</a></td><td>42 B</td></tr> +<tr class="odd_row"><td><a href="cox_rglasso.log">cox_rglasso.log</a></td><td>522 B</td></tr> +<tr><td><a href="rglasso.Rscript">rglasso.Rscript</a></td><td>21.5 KB</td></tr> +<tr class="odd_row"><td><a href="rglasso_error.log">rglasso_error.log</a></td><td>803 B</td></tr> +<tr><td><a href="rglasso_runner.log">rglasso_runner.log</a></td><td>1.7 KB</td></tr> +</table></div><br/> +</div></body></html> +
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/test-data/coxlassotest_modelres.xls Sat Oct 31 02:26:24 2015 -0400 @@ -0,0 +1,2 @@ +regulator partial_likelihood forced_in glmnet_model best_lambda lambdaChoice alpha +rhubarb 0.0012215 TRUE cox 0.04197 lambda.1se 0.95
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/test-data/genTest.R Sat Oct 31 02:26:24 2015 -0400 @@ -0,0 +1,10 @@ +ids=c(1:50) +io1 = rep(c(0,0,0,0,1),10) +ip2 = runif(50)+0.1 +ip2[which(ip2>1.0)] = 1.0 +ip1 = runif(50)+0.05 +ip1[which(ip1>1.0)] = 1.0 +df=data.frame(id=ids,input1_observed=io1,input1_predicted=ip1,input2_predicted=ip2) +fout='test-data/nri_test1.xls' +write.table(df,file=fout, quote=FALSE, sep="\t",row.names=F) +# planemo test --job_output_files /home/rlazarus/tmp --test_output /home/rlazarus/tmp/startest/foo.html --update_test_data --galaxy_root /home/rlazarus/galaxy rg_nri.xml
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/test-data/nri_test1.xls Sat Oct 31 02:26:24 2015 -0400 @@ -0,0 +1,51 @@ +id input1_observed input1_predicted input2_predicted +1 0 0.431491391919553 0.786344331596047 +2 0 0.235500678373501 0.37771512279287 +3 0 0.458084875205532 0.343597211781889 +4 0 0.400509147858247 0.50323077281937 +5 1 0.463376078708097 0.28904068809934 +6 0 0.638350193109363 0.745497886231169 +7 0 0.909975615050644 1 +8 0 0.0712537909392267 0.613464191928506 +9 0 0.475019342219457 0.757945845043287 +10 1 0.146289554610848 0.948540200153366 +11 0 0.947521302569658 0.452743433788419 +12 0 0.512624199455604 0.107415618747473 +13 0 0.260238706320524 0.165854234388098 +14 0 0.252348658815026 0.866530107706785 +15 1 0.552841167151928 0.781692814920098 +16 0 0.203356911847368 0.397802447108552 +17 0 0.138653858751059 0.966900746012107 +18 0 0.237485485337675 0.126699626818299 +19 0 0.317045996850356 0.646952681383118 +20 1 0.271802965737879 0.828731852956116 +21 0 0.522744563594461 0.421856845868751 +22 0 0.819749742234126 0.628745510149747 +23 0 0.27876661689952 0.826974717108533 +24 0 0.077877387823537 0.12339704008773 +25 1 0.911161796143279 1 +26 0 0.567042255960405 1 +27 0 0.223712538788095 0.738010874623433 +28 0 0.306199585506693 0.388159935269505 +29 0 0.305726604955271 0.618934138352051 +30 1 0.919456032756716 1 +31 0 0.493467168603092 0.470883537689224 +32 0 0.626042458415031 0.753118774993345 +33 0 0.739080456737429 0.812446624180302 +34 0 0.786149174161255 0.768807894224301 +35 1 0.65491428244859 0.855509085487574 +36 0 0.62289561489597 0.585996081028134 +37 0 0.297479379642755 0.357051227381453 +38 0 0.796507286094129 0.775113077275455 +39 0 0.581927774194628 0.569900369411334 +40 1 0.76229884494096 0.877794308261946 +41 0 0.845583727071062 0.804134764615446 +42 0 0.79128703516908 1 +43 0 0.0507813768461347 0.97111975508742 +44 0 0.49284805515781 0.101141295302659 +45 1 0.665627157501876 0.121852198988199 +46 0 0.059045681450516 0.996109608234838 +47 0 0.0814651515800506 0.852028644317761 +48 0 0.162336849560961 1 +49 0 0.577456943178549 0.343083152920008 +50 1 0.905623372131959 0.651731353392825
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/test-data/nri_test1_out.html Sat Oct 31 02:26:24 2015 -0400 @@ -0,0 +1,34 @@ +<!DOCTYPE html PUBLIC "-//W3C//DTD XHTML 1.0 Transitional//EN" "http://www.w3.org/TR/xhtml1/DTD/xhtml1-transitional.dtd"> + <html xmlns="http://www.w3.org/1999/xhtml" xml:lang="en" lang="en"> + <head> <meta http-equiv="Content-Type" content="text/html; charset=utf-8" /> + <meta name="generator" content="Galaxy rgToolFactory.py tool output - see http://g2.trac.bx.psu.edu/" /> + <title></title> + <link rel="stylesheet" href="/static/style/base.css" type="text/css" /> + </head> + <body> + <div class="toolFormBody"> + +<div class="infomessage">Galaxy Tool "rg_NRI" run at 08/01/2015 16:12:57</div><br/> +<div class="toolFormTitle">rg log output</div> + +<pre> + +Error in library("e1071") : there is no package called ‘e1071’ + +Execution halted + + +</pre> + +<div class="toolFormTitle">Other log output</div> +/tmp/tmpq72Dni/job_working_directory/000/2/dataset_2_files/rg_NRI_runner.log is empty<br/> +<div class="toolFormTitle">All output files available for downloading</div> + +<div><table class="colored" cellpadding="3" cellspacing="3"><tr><th>Output File Name (click to view)</th><th>Size</th></tr> + +<tr><td><a href="rg_NRI.Rscript">rg_NRI.Rscript</a></td><td>18.1 KB</td></tr> +<tr class="odd_row"><td><a href="rg_NRI_error.log">rg_NRI_error.log</a></td><td>84 B</td></tr> +<tr><td><a href="rg_NRI_runner.log">rg_NRI_runner.log</a></td><td>48 B</td></tr> +</table></div><br/> +</div></body></html> +
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/tool_dependencies.xml Sat Oct 31 02:26:24 2015 -0400 @@ -0,0 +1,120 @@ +<?xml version="1.0"?> +<tool_dependency> + <package name="R_3_2_2" version="3.2.2"> + <repository changeset_revision="883acf7a3ddb" name="package_r_3_2_2" owner="mvdbeek" prior_installation_required="True" toolshed="https://testtoolshed.g2.bx.psu.edu" /> + </package> + <package name="graphicsmagick" version="1.3.18"> + <repository changeset_revision="bff3f66adff2" name="package_graphicsmagick_1_3" owner="iuc" prior_installation_required="True" toolshed="https://testtoolshed.g2.bx.psu.edu" /> + </package> + <package name="ghostscript" version="9.10"> + <repository changeset_revision="9345d2740f0c" name="package_ghostscript_9_10" owner="devteam" prior_installation_required="True" toolshed="https://testtoolshed.g2.bx.psu.edu" /> + </package> + <package name="glmnet_lars_3_2" version="3.2"> + <install version="1.0"> + <actions> + <action type="setup_r_environment"> + <repository changeset_revision="883acf7a3ddb" name="package_r_3_2_2" owner="mvdbeek" toolshed="https://testtoolshed.g2.bx.psu.edu"> + <package name="R_3_2_2" version="3.2.2" /> + </repository> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/RColorBrewer_1.0-5.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/RColorBrewer_1.1-2.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/Rcpp_0.11.3.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/Rcpp_0.12.1.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/RcppArmadillo_0.4.450.1.0.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/RcppArmadillo_0.4.500.0.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/RcppArmadillo_0.4.550.1.0.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/RcppArmadillo_0.4.600.0.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/RcppEigen_0.3.2.5.1.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/e1071_1.6-4.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/e1071_1.6-7.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/plyr_1.8.1.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/plyr_1.8.3.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/digest_0.6.4.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/digest_0.6.8.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/gtable_0.1.2.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/stringi_1.0-1.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/magrittr_1.5.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/stringr_0.6.2.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/stringr_1.0.0.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/reshape2_1.4.1.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/reshape2_1.4.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/dichromat_2.0-0.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/colorspace_1.2-4.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/colorspace_1.2-6.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/munsell_0.4.2.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/labeling_0.3.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/scales_0.2.4.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/scales_0.3.0.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/proto_0.3-10.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/ggplot2_1.0.0.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/ggplot2_1.0.1.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/minqa_1.2.4.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/nloptr_1.0.4.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/lme4_1.1-10.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/lme4_1.1-7.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/pbkrtest_0.4-2.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/SparseM_1.05.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/SparseM_1.6.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/SparseM_1.7.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/MatrixModels_0.4-1.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/quantreg_5.05.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/quantreg_5.19.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/car_2.0-22.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/car_2.1-0.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/caret_6.0-35.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/caret_6.0-37.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/caret_6.0-41.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/caret_6.0-58.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/iterators_1.0.7.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/iterators_1.0.8.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/foreach_1.4.2.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/foreach_1.4.3.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/pROC_1.7.3.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/pROC_1.8.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/Formula_1.1-2.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/Formula_1.2-1.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/latticeExtra_0.6-26.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/acepack_1.3-3.3.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/gridExtra_2.0.0.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/Hmisc_3.14-5.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/Hmisc_3.14-6.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/Hmisc_3.17-0.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/pracma_1.7.3.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/pracma_1.7.9.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/pracma_1.8.6.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/survival_2.37-7.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/survival_2.38-3.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/lars_1.2.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/glmnet_1.9-8.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/glmnet_2.0-2.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/numDeriv_2012.9-1.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/numDeriv_2014.2-1.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/lava_1.3.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/lava_1.4.1.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/prodlim_1.5.1.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/prodlim_1.5.5.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/polspline_1.1.12.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/polspline_1.1.9.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/mvtnorm_1.0-2.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/mvtnorm_1.0-3.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/TH.data_1.0-5.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/TH.data_1.0-6.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/zoo_1.7-11.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/zoo_1.7-12.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/sandwich_2.3-2.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/sandwich_2.3-4.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/multcomp_1.3-8.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/multcomp_1.4-1.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/rms_4.2-1.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/rms_4.4-0.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/pec_2.3.7.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/pec_2.4.4.tar.gz?raw=true</package> + <package>https://github.com/fubar2/galaxy_tool_source/blob/master/RELEASE_2_14/pec_2.4.7.tar.gz?raw=true</package> + </action> + </actions> + </install> + <readme> + Yee-Haw! Lasso for Galaxy! + </readme> + </package> +</tool_dependency>