Added a new tool set for performing multivariate statistical analyses

This commit is contained in:
Guruprasad Anada
2010-03-04 11:48:10 -05:00
parent a6a86d0085
commit 45cffbb8ed
10 changed files with 1020 additions and 3 deletions
+6
View File
@@ -150,6 +150,12 @@
<tool file="regVariation/best_regression_subsets.xml" />
<tool file="regVariation/rcve.xml" />
</section>
<section name="Multivariate Analysis" id="multVar">
<tool file="multivariate_stats/pca.xml" />
<tool file="multivariate_stats/cca.xml" />
<tool file="multivariate_stats/kpca.xml" />
<tool file="multivariate_stats/kcca.xml" />
</section>
<section name="Evolution" id="hyphy">
<tool file="hyphy/hyphy_branch_lengths_wrapper.xml" />
<tool file="hyphy/hyphy_nj_tree_wrapper.xml" />
+159
View File
@@ -0,0 +1,159 @@
#!/usr/bin/env python
from galaxy import eggs
import sys, string
from rpy import *
import numpy
def stop_err(msg):
sys.stderr.write(msg)
sys.exit()
infile = sys.argv[1]
x_cols = sys.argv[2].split(',')
y_cols = sys.argv[3].split(',')
x_scale = x_center = "FALSE"
if sys.argv[4] == 'both':
x_scale = x_center = "TRUE"
elif sys.argv[4] == 'center':
x_center = "TRUE"
elif sys.argv[4] == 'scale':
x_scale = "TRUE"
y_scale = y_center = "FALSE"
if sys.argv[5] == 'both':
y_scale = y_center = "TRUE"
elif sys.argv[5] == 'center':
y_center = "TRUE"
elif sys.argv[5] == 'scale':
y_scale = "TRUE"
std_scores = "FALSE"
if sys.argv[6] == "yes":
std_scores = "TRUE"
outfile = sys.argv[7]
outfile2 = sys.argv[8]
fout = open(outfile,'w')
elems = []
for i, line in enumerate( file ( infile )):
line = line.rstrip('\r\n')
if len( line )>0 and not line.startswith( '#' ):
elems = line.split( '\t' )
break
if i == 30:
break # Hopefully we'll never get here...
if len( elems )<1:
stop_err( "The data in your input dataset is either missing or not formatted properly." )
x_vals = []
for k,col in enumerate(x_cols):
x_cols[k] = int(col)-1
x_vals.append([])
y_vals = []
for k,col in enumerate(y_cols):
y_cols[k] = int(col)-1
y_vals.append([])
skipped = 0
for ind,line in enumerate( file( infile )):
if line and not line.startswith( '#' ):
try:
fields = line.strip().split("\t")
valid_line = True
for col in x_cols+y_cols:
try:
assert float(fields[col])
except:
skipped += 1
valid_line = False
break
if valid_line:
for k,col in enumerate(x_cols):
try:
xval = float(fields[col])
except:
xval = NaN#
x_vals[k].append(xval)
for k,col in enumerate(y_cols):
try:
yval = float(fields[col])
except:
yval = NaN#
y_vals[k].append(yval)
except:
skipped += 1
x_vals1 = numpy.asarray(x_vals).transpose()
y_vals1 = numpy.asarray(y_vals).transpose()
x_dat= r.list(array(x_vals1))
y_dat= r.list(array(y_vals1))
try:
r.suppressWarnings(r.library("yacca"))
except:
stop_err("Missing R library yacca.")
set_default_mode(NO_CONVERSION)
try:
xcolnames = ["c%d" %(el+1) for el in x_cols]
ycolnames = ["c%d" %(el+1) for el in y_cols]
cc = r.cca(x=x_dat, y=y_dat, xlab=xcolnames, ylab=ycolnames, xcenter=r(x_center), ycenter=r(y_center), xscale=r(x_scale), yscale=r(y_scale), standardize_scores=r(std_scores))
ftest = r.F_test_cca(cc)
except RException, rex:
stop_err("Encountered error while performing CCA on the input data: %s" %(rex))
set_default_mode(BASIC_CONVERSION)
summary = r.summary(cc)
ncomps = len(summary['corr'])
comps = summary['corr'].keys()
corr = summary['corr'].values()
xlab = summary['xlab']
ylab = summary['ylab']
for i in range(ncomps):
corr[comps.index('CV %s' %(i+1))] = summary['corr'].values()[i]
ftest=ftest.as_py()
print >>fout, "#Component\t%s" %("\t".join(["%s" % el for el in range(1,ncomps+1)]))
print >>fout, "#Correlation\t%s" %("\t".join(["%s" % el for el in corr]))
print >>fout, "#F-statistic\t%s" %("\t".join(["%s" % el for el in ftest['statistic']]))
print >>fout, "#p-value\t%s" %("\t".join(["%s" % el for el in ftest['p.value']]))
print >>fout, "#X-Coefficients\t%s" %("\t".join(["%s" % el for el in range(1,ncomps+1)]))
for i,val in enumerate(summary['xcoef']):
print >>fout, "%s\t%s" %(xlab[i], "\t".join(["%s" % el for el in val]))
print >>fout, "#Y-Coefficients\t%s" %("\t".join(["%s" % el for el in range(1,ncomps+1)]))
for i,val in enumerate(summary['ycoef']):
print >>fout, "%s\t%s" %(ylab[i], "\t".join(["%s" % el for el in val]))
print >>fout, "#X-Loadings\t%s" %("\t".join(["%s" % el for el in range(1,ncomps+1)]))
for i,val in enumerate(summary['xstructcorr']):
print >>fout, "%s\t%s" %(xlab[i], "\t".join(["%s" % el for el in val]))
print >>fout, "#Y-Loadings\t%s" %("\t".join(["%s" % el for el in range(1,ncomps+1)]))
for i,val in enumerate(summary['ystructcorr']):
print >>fout, "%s\t%s" %(ylab[i], "\t".join(["%s" % el for el in val]))
print >>fout, "#X-CrossLoadings\t%s" %("\t".join(["%s" % el for el in range(1,ncomps+1)]))
for i,val in enumerate(summary['xcrosscorr']):
print >>fout, "%s\t%s" %(xlab[i], "\t".join(["%s" % el for el in val]))
print >>fout, "#Y-CrossLoadings\t%s" %("\t".join(["%s" % el for el in range(1,ncomps+1)]))
for i,val in enumerate(summary['ycrosscorr']):
print >>fout, "%s\t%s" %(ylab[i], "\t".join(["%s" % el for el in val]))
r.pdf( outfile2, 8, 8 )
#r.plot(cc)
for i in range(ncomps):
r.helio_plot(cc, cv = i+1, main = r.paste("Explained Variance for CV",i+1), type = "variance")
r.dev_off()
+95
View File
@@ -0,0 +1,95 @@
<tool id="cca1" name="Canonical Correlation Analysis" version="1.0.0">
<description> </description>
<command interpreter="python">
cca.py
$input1
$x_cols
$y_cols
$x_scale
$y_scale
$std_scores
$out_file1
$out_file2
</command>
<inputs>
<param format="tabular" name="input1" type="data" label="Select data" help="Query missing? See TIP below."/>
<param name="x_cols" label="Select columns containing X variables " type="data_column" data_ref="input1" numerical="True" multiple="true" >
<validator type="no_options" message="Please select at least one column."/>
</param>
<param name="y_cols" label="Select columns containing Y variables " type="data_column" data_ref="input1" numerical="True" multiple="true" >
<validator type="no_options" message="Please select at least one column."/>
</param>
<param name="x_scale" type="select" label="Type of Scaling for X variables" help="Can be used to center and/or scale variables">
<option value="none" selected="true">None</option>
<option value="center">Center only</option>
<option value="scale">Scale only</option>
<option value="both">Center and Scale</option>
</param>
<param name="y_scale" type="select" label="Type of Scaling for Y variables" help="Can be used to center and/or scale variables">
<option value="none" selected="true">None</option>
<option value="center">Center only</option>
<option value="scale">Scale only</option>
<option value="both">Center and Scale</option>
</param>
<param name="std_scores" type="select" label="Report standardized scores?" help="Selecting 'Yes' will rescale scores (and coefficients) to produce scores of unit variance">
<option value="no" selected="true">No</option>
<option value="yes">Yes</option>
</param>
</inputs>
<outputs>
<data format="input" name="out_file1" metadata_source="input1" />
<data format="pdf" name="out_file2" />
</outputs>
<requirements>
<requirement type="python-module">rpy</requirement>
</requirements>
<tests>
<test>
<param name="input1" value="iris.tabular"/>
<param name="x_cols" value="3,4"/>
<param name="y_cols" value="1,2"/>
<param name="x_scale" value="both"/>
<param name="y_scale" value="scale"/>
<param name="std_scores" value="yes"/>
<output name="out_file1" file="cca_out1.tabular"/>
<output name="out_file2" file="cca_out2.pdf"/>
</test>
</tests>
<help>
.. class:: infomark
**TIP:** If your data is not TAB delimited, use *Edit Queries-&gt;Convert characters*
-----
.. class:: infomark
**What it does**
This tool uses functions from 'yacca' library from R statistical package to perform Canonical Correlation Analysis (CCA) on the input data. It outputs two files, one containing the summary statistics of the performed CCA, and the other containing helioplots, which display structural loadings of X and Y variables on different canonical components.
*Carter T. Butts (2009). yacca: Yet Another Canonical Correlation Analysis Package. R package version 1.1.*
-----
.. class:: warningmark
**Note**
- This tool currently treats all predictor and response variables as continuous numeric variables. Running the tool on categorical variables might result in incorrect results.
- Rows containing non-numeric (or missing) data in any of the chosen columns will be skipped from the analysis.
- The summary statistics in the output are described below:
- correlation: Canonical correlation between the canonical variates (i.e. transformed variables)
- F-statistic: F-value obtained from F Test for Canonical Correlations Using Rao's Approximation
- p-value: denotes significance of canonical correlations
- Coefficients: represent the coefficients of X and Y variables on each canonical variate
- Loadings: represent the correlations between the original variables in each set and their respective canonical variates
- CrossLoadings: represent the correlations between the original variables in each set and the opposite canonical variates
</help>
</tool>
+146
View File
@@ -0,0 +1,146 @@
#!/usr/bin/env python
"""
Run kernel CCA using kcca() from R 'kernlab' package
usage: %prog [options]
-i, --input=i: Input file
-o, --output1=o: Summary output
-x, --x_cols=x: X-Variable columns
-y, --y_cols=y: Y-Variable columns
-k, --kernel=k: Kernel function
-f, --features=f: Number of canonical components to return
-s, --sigma=s: sigma
-d, --degree=d: degree
-l, --scale=l: scale
-t, --offset=t: offset
-r, --order=r: order
usage: %prog input output1 x_cols y_cols kernel features sigma(or_None) degree(or_None) scale(or_None) offset(or_None) order(or_None)
"""
from galaxy import eggs
import sys, string
from rpy import *
import numpy
import pkg_resources; pkg_resources.require( "bx-python" )
from bx.cookbook import doc_optparse
def stop_err(msg):
sys.stderr.write(msg)
sys.exit()
#Parse Command Line
options, args = doc_optparse.parse( __doc__ )
#{'options= kernel': 'rbfdot', 'var_cols': '1,2,3,4', 'degree': 'None', 'output2': '/afs/bx.psu.edu/home/gua110/workspace/galaxy_bitbucket/database/files/000/dataset_260.dat', 'output1': '/afs/bx.psu.edu/home/gua110/workspace/galaxy_bitbucket/database/files/000/dataset_259.dat', 'scale': 'None', 'offset': 'None', 'input': '/afs/bx.psu.edu/home/gua110/workspace/galaxy_bitbucket/database/files/000/dataset_256.dat', 'sigma': '1.0', 'order': 'None'}
infile = options.input
x_cols = options.x_cols.split(',')
y_cols = options.y_cols.split(',')
kernel = options.kernel
outfile = options.output1
ncomps = int(options.features)
fout = open(outfile,'w')
if ncomps < 1:
print "You chose to return '0' canonical components. Please try rerunning the tool with number of components = 1 or more."
sys.exit()
elems = []
for i, line in enumerate( file ( infile )):
line = line.rstrip('\r\n')
if len( line )>0 and not line.startswith( '#' ):
elems = line.split( '\t' )
break
if i == 30:
break # Hopefully we'll never get here...
if len( elems )<1:
stop_err( "The data in your input dataset is either missing or not formatted properly." )
x_vals = []
for k,col in enumerate(x_cols):
x_cols[k] = int(col)-1
x_vals.append([])
y_vals = []
for k,col in enumerate(y_cols):
y_cols[k] = int(col)-1
y_vals.append([])
NA = 'NA'
skipped = 0
for ind,line in enumerate( file( infile )):
if line and not line.startswith( '#' ):
try:
fields = line.strip().split("\t")
valid_line = True
for col in x_cols+y_cols:
try:
assert float(fields[col])
except:
skipped += 1
valid_line = False
break
if valid_line:
for k,col in enumerate(x_cols):
try:
xval = float(fields[col])
except:
xval = NaN#
x_vals[k].append(xval)
for k,col in enumerate(y_cols):
try:
yval = float(fields[col])
except:
yval = NaN#
y_vals[k].append(yval)
except:
skipped += 1
x_vals1 = numpy.asarray(x_vals).transpose()
y_vals1 = numpy.asarray(y_vals).transpose()
x_dat= r.list(array(x_vals1))
y_dat= r.list(array(y_vals1))
try:
r.suppressWarnings(r.library('kernlab'))
except:
stop_err('Missing R library kernlab')
set_default_mode(NO_CONVERSION)
if kernel=="rbfdot" or kernel=="anovadot":
pars = r.list(sigma=float(options.sigma))
elif kernel=="polydot":
pars = r.list(degree=float(options.degree),scale=float(options.scale),offset=float(options.offset))
elif kernel=="tanhdot":
pars = r.list(scale=float(options.scale),offset=float(options.offset))
elif kernel=="besseldot":
pars = r.list(degree=float(options.degree),sigma=float(options.sigma),order=float(options.order))
elif kernel=="anovadot":
pars = r.list(degree=float(options.degree),sigma=float(options.sigma))
else:
pars = rlist()
try:
kcc = r.kcca(x=x_dat, y=y_dat, kernel=kernel, kpar=pars, ncomps=ncomps)
except RException, rex:
stop_err("Encountered error while performing kCCA on the input data: %s" %(rex))
set_default_mode(BASIC_CONVERSION)
kcor = r.kcor(kcc)
if ncomps == 1:
kcor = [kcor]
xcoef = r.xcoef(kcc)
ycoef = r.ycoef(kcc)
print >>fout, "#Component\t%s" %("\t".join(["%s" % el for el in range(1,ncomps+1)]))
print >>fout, "#Correlation\t%s" %("\t".join(["%s" % el for el in kcor]))
print >>fout, "#Estimated X-coefficients\t%s" %("\t".join(["%s" % el for el in range(1,ncomps+1)]))
for obs,val in enumerate(xcoef):
print >>fout, "%s\t%s" %(obs+1, "\t".join(["%s" % el for el in val]))
print >>fout, "#Estimated Y-coefficients\t%s" %("\t".join(["%s" % el for el in range(1,ncomps+1)]))
for obs,val in enumerate(ycoef):
print >>fout, "%s\t%s" %(obs+1, "\t".join(["%s" % el for el in val]))
+150
View File
@@ -0,0 +1,150 @@
<tool id="kcca1" name="Kernel Canonical Correlation Analysis" version="1.0.0">
<description> </description>
<command interpreter="python">
kcca.py
--input=$input1
--output1=$out_file1
--x_cols=$x_cols
--y_cols=$y_cols
--kernel=$kernelChoice.kernel
--features=$features
#if $kernelChoice.kernel == "rbfdot" or $kernelChoice.kernel == "anovadot":
--sigma=$kernelChoice.sigma
--degree="None"
--scale="None"
--offset="None"
--order="None"
#elif $kernelChoice.kernel == "polydot":
--sigma="None"
--degree=$kernelChoice.degree
--scale=$kernelChoice.scale
--offset=$kernelChoice.offset
--order="None"
#elif $kernelChoice.kernel == "tanhdot":
--sigma="None"
--degree="None"
--scale=$kernelChoice.scale
--offset=$kernelChoice.offset
--order="None"
#elif $kernelChoice.kernel == "besseldot":
--sigma=$kernelChoice.sigma
--degree=$kernelChoice.degree
--scale="None"
--offset="None"
--order=$kernelChoice.order
#elif $kernelChoice.kernel == "anovadot":
--sigma=$kernelChoice.sigma
--degree=$kernelChoice.degree
--scale="None"
--offset="None"
--order="None"
#else:
--sigma="None"
--degree="None"
--scale="None"
--offset="None"
--order="None"
#end if
</command>
<inputs>
<param format="tabular" name="input1" type="data" label="Select data" help="Query missing? See TIP below."/>
<param name="x_cols" label="Select columns containing X variables " type="data_column" data_ref="input1" numerical="True" multiple="true" >
<validator type="no_options" message="Please select at least one column."/>
</param>
<param name="y_cols" label="Select columns containing Y variables " type="data_column" data_ref="input1" numerical="True" multiple="true" >
<validator type="no_options" message="Please select at least one column."/>
</param>
<param name="features" size="10" type="integer" value="2" label="Number of canonical components to return" help="Enter an integer value greater than 0"/>
<conditional name="kernelChoice">
<param name="kernel" type="select" label="Kernel function">
<option value="rbfdot" selected="true">Gaussian Radial Basis Function</option>
<option value="polydot">Polynomial</option>
<option value="vanilladot">Linear</option>
<option value="tanhdot">Hyperbolic</option>
<option value="laplacedot">Laplacian</option>
<option value="besseldot">Bessel</option>
<option value="anovadot">ANOVA Radial Basis Function</option>
<option value="splinedot">Spline</option>
</param>
<when value="vanilladot" />
<when value="splinedot" />
<when value="rbfdot">
<param name="sigma" size="10" type="float" value="1" label="sigma (inverse kernel width)" />
</when>
<when value="laplacedot">
<param name="sigma" size="10" type="float" value="1" label="sigma (inverse kernel width)" />
</when>
<when value="polydot">
<param name="degree" size="10" type="float" value="1" label="degree" />
<param name="scale" size="10" type="float" value="1" label="scale" />
<param name="offset" size="10" type="float" value="1" label="offset" />
</when>
<when value="tanhdot">
<param name="scale" size="10" type="float" value="1" label="scale" />
<param name="offset" size="10" type="float" value="1" label="offset" />
</when>
<when value="besseldot">
<param name="sigma" size="10" type="float" value="1" label="sigma" />
<param name="order" size="10" type="float" value="1" label="order" />
<param name="degree" size="10" type="float" value="1" label="degree" />
</when>
<when value="anovadot">
<param name="sigma" size="10" type="float" value="1" label="sigma" />
<param name="degree" size="10" type="float" value="1" label="degree" />
</when>
</conditional>
</inputs>
<outputs>
<data format="input" name="out_file1" metadata_source="input1" />
</outputs>
<requirements>
<requirement type="python-module">rpy</requirement>
</requirements>
<tests>
<test>
<param name="input1" value="iris.tabular"/>
<param name="x_cols" value="1,2"/>
<param name="y_cols" value="3,4"/>
<param name="kernel" value="anovadot"/>
<param name="features" value="4"/>
<param name="sigma" value="0.1"/>
<param name="degree" value="2"/>
<output name="out_file1" file="kcca_out1.tabular"/>
</test>
<test>
<param name="input1" value="iris.tabular"/>
<param name="x_cols" value="3,4"/>
<param name="y_cols" value="1,2"/>
<param name="kernel" value="rbfdot"/>
<param name="features" value="2"/>
<param name="sigma" value="0.5"/>
<output name="out_file1" file="kcca_out2.tabular"/>
</test>
</tests>
<help>
.. class:: infomark
**TIP:** If your data is not TAB delimited, use *Edit Queries-&gt;Convert characters*
-----
.. class:: infomark
**What it does**
This tool uses functions from 'kernlab' library from R statistical package to perform Kernel Canonical Correlation Analysis (kCCA) on the input data.
*Alexandros Karatzoglou, Alex Smola, Kurt Hornik, Achim Zeileis (2004). kernlab - An S4 Package for Kernel Methods in R. Journal of Statistical Software 11(9), 1-20. URL http://www.jstatsoft.org/v11/i09/*
-----
.. class:: warningmark
**Note**
This tool currently treats all variables as continuous numeric variables. Running the tool on categorical variables might result in incorrect results. Rows containing non-numeric (or missing) data in any of the chosen columns will be skipped from the analysis.
</help>
</tool>
+134
View File
@@ -0,0 +1,134 @@
#!/usr/bin/env python
"""
Run kernel PCA using kpca() from R 'kernlab' package
usage: %prog [options]
-i, --input=i: Input file
-o, --output1=o: Summary output
-p, --output2=p: Figures output
-c, --var_cols=c: Variable columns
-k, --kernel=k: Kernel function
-f, --features=f: Number of principal components to return
-s, --sigma=s: sigma
-d, --degree=d: degree
-l, --scale=l: scale
-t, --offset=t: offset
-r, --order=r: order
usage: %prog input output1 output2 var_cols kernel features sigma(or_None) degree(or_None) scale(or_None) offset(or_None) order(or_None)
"""
from galaxy import eggs
import sys, string
from rpy import *
import numpy
import pkg_resources; pkg_resources.require( "bx-python" )
from bx.cookbook import doc_optparse
def stop_err(msg):
sys.stderr.write(msg)
sys.exit()
#Parse Command Line
options, args = doc_optparse.parse( __doc__ )
#{'options= kernel': 'rbfdot', 'var_cols': '1,2,3,4', 'degree': 'None', 'output2': '/afs/bx.psu.edu/home/gua110/workspace/galaxy_bitbucket/database/files/000/dataset_260.dat', 'output1': '/afs/bx.psu.edu/home/gua110/workspace/galaxy_bitbucket/database/files/000/dataset_259.dat', 'scale': 'None', 'offset': 'None', 'input': '/afs/bx.psu.edu/home/gua110/workspace/galaxy_bitbucket/database/files/000/dataset_256.dat', 'sigma': '1.0', 'order': 'None'}
infile = options.input
x_cols = options.var_cols.split(',')
kernel = options.kernel
outfile = options.output1
outfile2 = options.output2
ncomps = int(options.features)
fout = open(outfile,'w')
elems = []
for i, line in enumerate( file ( infile )):
line = line.rstrip('\r\n')
if len( line )>0 and not line.startswith( '#' ):
elems = line.split( '\t' )
break
if i == 30:
break # Hopefully we'll never get here...
if len( elems )<1:
stop_err( "The data in your input dataset is either missing or not formatted properly." )
x_vals = []
for k,col in enumerate(x_cols):
x_cols[k] = int(col)-1
x_vals.append([])
NA = 'NA'
skipped = 0
for ind,line in enumerate( file( infile )):
if line and not line.startswith( '#' ):
try:
fields = line.strip().split("\t")
for k,col in enumerate(x_cols):
try:
xval = float(fields[col])
except:
#xval = r('NA')
xval = NaN#
x_vals[k].append(xval)
except:
skipped += 1
x_vals1 = numpy.asarray(x_vals).transpose()
dat= r.list(array(x_vals1))
try:
r.suppressWarnings(r.library('kernlab'))
except:
stop_err('Missing R library kernlab')
set_default_mode(NO_CONVERSION)
if kernel=="rbfdot" or kernel=="anovadot":
pars = r.list(sigma=float(options.sigma))
elif kernel=="polydot":
pars = r.list(degree=float(options.degree),scale=float(options.scale),offset=float(options.offset))
elif kernel=="tanhdot":
pars = r.list(scale=float(options.scale),offset=float(options.offset))
elif kernel=="besseldot":
pars = r.list(degree=float(options.degree),sigma=float(options.sigma),order=float(options.order))
elif kernel=="anovadot":
pars = r.list(degree=float(options.degree),sigma=float(options.sigma))
else:
pars = rlist()
try:
kpc = r.kpca(x=r.na_exclude(dat), kernel=kernel, kpar=pars, features=ncomps)
except RException, rex:
stop_err("Encountered error while performing kPCA on the input data: %s" %(rex))
set_default_mode(BASIC_CONVERSION)
eig = r.eig(kpc)
pcv = r.pcv(kpc)
rotated = r.rotated(kpc)
comps = eig.keys()
eigv = eig.values()
for i in range(ncomps):
eigv[comps.index('Comp.%s' %(i+1))] = eig.values()[i]
print >>fout, "#Component\t%s" %("\t".join(["%s" % el for el in range(1,ncomps+1)]))
print >>fout, "#Eigenvalue\t%s" %("\t".join(["%s" % el for el in eig.values()]))
print >>fout, "#Principal component vectors\t%s" %("\t".join(["%s" % el for el in range(1,ncomps+1)]))
for obs,val in enumerate(pcv):
print >>fout, "%s\t%s" %(obs+1, "\t".join(["%s" % el for el in val]))
print >>fout, "#Rotated values\t%s" %("\t".join(["%s" % el for el in range(1,ncomps+1)]))
for obs,val in enumerate(rotated):
print >>fout, "%s\t%s" %(obs+1, "\t".join(["%s" % el for el in val]))
r.pdf( outfile2, 8, 8 )
if ncomps != 1:
r.pairs(rotated,labels=r.list(range(1,ncomps+1)),main="Scatterplot of rotated values")
else:
r.plot(rotated, ylab='Comp.1', main="Scatterplot of rotated values")
r.dev_off()
+151
View File
@@ -0,0 +1,151 @@
<tool id="kpca1" name="Kernel Principal Component Analysis" version="1.0.0">
<description> </description>
<command interpreter="python">
kpca.py
--input=$input1
--output1=$out_file1
--output2=$out_file2
--var_cols=$var_cols
--kernel=$kernelChoice.kernel
--features=$features
#if $kernelChoice.kernel == "rbfdot" or $kernelChoice.kernel == "anovadot":
--sigma=$kernelChoice.sigma
--degree="None"
--scale="None"
--offset="None"
--order="None"
#elif $kernelChoice.kernel == "polydot":
--sigma="None"
--degree=$kernelChoice.degree
--scale=$kernelChoice.scale
--offset=$kernelChoice.offset
--order="None"
#elif $kernelChoice.kernel == "tanhdot":
--sigma="None"
--degree="None"
--scale=$kernelChoice.scale
--offset=$kernelChoice.offset
--order="None"
#elif $kernelChoice.kernel == "besseldot":
--sigma=$kernelChoice.sigma
--degree=$kernelChoice.degree
--scale="None"
--offset="None"
--order=$kernelChoice.order
#elif $kernelChoice.kernel == "anovadot":
--sigma=$kernelChoice.sigma
--degree=$kernelChoice.degree
--scale="None"
--offset="None"
--order="None"
#else:
--sigma="None"
--degree="None"
--scale="None"
--offset="None"
--order="None"
#end if
</command>
<inputs>
<param format="tabular" name="input1" type="data" label="Select data" help="Query missing? See TIP below."/>
<param name="var_cols" label="Select columns containing input variables " type="data_column" data_ref="input1" numerical="True" multiple="true" >
<validator type="no_options" message="Please select at least one column."/>
</param>
<param name="features" size="10" type="integer" value="2" label="Number of principal components to return" help="To return all, enter 0"/>
<conditional name="kernelChoice">
<param name="kernel" type="select" label="Kernel function">
<option value="rbfdot" selected="true">Gaussian Radial Basis Function</option>
<option value="polydot">Polynomial</option>
<option value="vanilladot">Linear</option>
<option value="tanhdot">Hyperbolic</option>
<option value="laplacedot">Laplacian</option>
<option value="besseldot">Bessel</option>
<option value="anovadot">ANOVA Radial Basis Function</option>
<option value="splinedot">Spline</option>
</param>
<when value="vanilladot" />
<when value="splinedot" />
<when value="rbfdot">
<param name="sigma" size="10" type="float" value="1" label="sigma (inverse kernel width)" />
</when>
<when value="laplacedot">
<param name="sigma" size="10" type="float" value="1" label="sigma (inverse kernel width)" />
</when>
<when value="polydot">
<param name="degree" size="10" type="integer" value="1" label="degree" />
<param name="scale" size="10" type="integer" value="1" label="scale" />
<param name="offset" size="10" type="integer" value="1" label="offset" />
</when>
<when value="tanhdot">
<param name="scale" size="10" type="integer" value="1" label="scale" />
<param name="offset" size="10" type="integer" value="1" label="offset" />
</when>
<when value="besseldot">
<param name="sigma" size="10" type="integer" value="1" label="sigma" />
<param name="order" size="10" type="integer" value="1" label="order" />
<param name="degree" size="10" type="integer" value="1" label="degree" />
</when>
<when value="anovadot">
<param name="sigma" size="10" type="integer" value="1" label="sigma" />
<param name="degree" size="10" type="integer" value="1" label="degree" />
</when>
</conditional>
</inputs>
<outputs>
<data format="input" name="out_file1" metadata_source="input1" />
<data format="pdf" name="out_file2" />
</outputs>
<requirements>
<requirement type="python-module">rpy</requirement>
</requirements>
<tests>
<test>
<param name="input1" value="iris.tabular"/>
<param name="var_cols" value="1,2,3,4"/>
<param name="kernel" value="polydot"/>
<param name="features" value="2"/>
<param name="offset" value="0"/>
<param name="scale" value="1"/>
<param name="degree" value="2"/>
<output name="out_file1" file="kpca_out1.tabular"/>
<output name="out_file2" file="kpca_out2.pdf"/>
</test>
<test>
<param name="input1" value="iris.tabular"/>
<param name="var_cols" value="2,3,4"/>
<param name="kernel" value="besseldot"/>
<param name="features" value="1"/>
<param name="sigma" value="1"/>
<param name="order" value="1"/>
<param name="degree" value="1"/>
<output name="out_file1" file="kpca_out3.tabular"/>
<output name="out_file2" file="kpca_out4.pdf"/>
</test>
</tests>
<help>
.. class:: infomark
**TIP:** If your data is not TAB delimited, use *Edit Queries-&gt;Convert characters*
-----
.. class:: infomark
**What it does**
This tool uses functions from 'kernlab' library from R statistical package to perform Kernel Principal Component Analysis (kPCA) on the input data. It outputs two files, one containing the summary statistics of the performed kPCA, and the other containing a scatterplot matrix of rotated values reported by kPCA.
*Alexandros Karatzoglou, Alex Smola, Kurt Hornik, Achim Zeileis (2004). kernlab - An S4 Package for Kernel Methods in R. Journal of Statistical Software 11(9), 1-20. URL http://www.jstatsoft.org/v11/i09/*
-----
.. class:: warningmark
**Note**
This tool currently treats all variables as continuous numeric variables. Running the tool on categorical variables might result in incorrect results. Rows containing non-numeric (or missing) data in any of the chosen columns will be skipped from the analysis.
</help>
</tool>
+98
View File
@@ -0,0 +1,98 @@
#!/usr/bin/env python
from galaxy import eggs
import sys, string
from rpy import *
import numpy
def stop_err(msg):
sys.stderr.write(msg)
sys.exit()
infile = sys.argv[1]
x_cols = sys.argv[2].split(',')
method = sys.argv[3]
outfile = sys.argv[4]
outfile2 = sys.argv[5]
fout = open(outfile,'w')
elems = []
for i, line in enumerate( file ( infile )):
line = line.rstrip('\r\n')
if len( line )>0 and not line.startswith( '#' ):
elems = line.split( '\t' )
break
if i == 30:
break # Hopefully we'll never get here...
if len( elems )<1:
stop_err( "The data in your input dataset is either missing or not formatted properly." )
x_vals = []
for k,col in enumerate(x_cols):
x_cols[k] = int(col)-1
x_vals.append([])
NA = 'NA'
skipped = 0
for ind,line in enumerate( file( infile )):
if line and not line.startswith( '#' ):
try:
fields = line.strip().split("\t")
for k,col in enumerate(x_cols):
try:
xval = float(fields[col])
except:
#xval = r('NA')
xval = NaN#
x_vals[k].append(xval)
except:
skipped += 1
x_vals1 = numpy.asarray(x_vals).transpose()
dat= r.list(array(x_vals1))
set_default_mode(NO_CONVERSION)
try:
if method == "cor":
pc = r.princomp(r.na_exclude(dat), cor = r("TRUE"))
else:
pc = r.princomp(r.na_exclude(dat), cor = r("FALSE"))
except RException, rex:
stop_err("Encountered error while performing PCA on the input data: %s" %(rex))
set_default_mode(BASIC_CONVERSION)
summary = r.summary(pc, loadings="TRUE")
ncomps = len(summary['sdev'])
comps = summary['sdev'].keys()
sd = summary['sdev'].values()
for i in range(ncomps):
sd[comps.index('Comp.%s' %(i+1))] = summary['sdev'].values()[i]
print >>fout, "#Component\t%s" %("\t".join(["%s" % el for el in range(1,ncomps+1)]))
print >>fout, "#Std. deviation\t%s" %("\t".join(["%s" % el for el in sd]))
total_var = 0
vars = []
for s in sd:
var = s*s
total_var += var
vars.append(var)
for i,var in enumerate(vars):
vars[i] = vars[i]/total_var
print >>fout, "#Proportion of variance explained\t%s" %("\t".join(["%s" % el for el in vars]))
print >>fout, "#Loadings\t%s" %("\t".join(["%s" % el for el in range(1,ncomps+1)]))
xcolnames = ["c%d" %(el+1) for el in x_cols]
for i,val in enumerate(summary['loadings']):
print >>fout, "%s\t%s" %(xcolnames[i], "\t".join(["%s" % el for el in val]))
print >>fout, "#Scores\t%s" %("\t".join(["%s" % el for el in range(1,ncomps+1)]))
for obs,sc in enumerate(summary['scores']):
print >>fout, "%s\t%s" %(obs+1, "\t".join(["%s" % el for el in sc]))
r.pdf( outfile2, 8, 8 )
r.biplot(pc)
r.dev_off()
+76
View File
@@ -0,0 +1,76 @@
<tool id="pca1" name="Principal Component Analysis" version="1.0.0">
<description> </description>
<command interpreter="python">
pca.py
$input1
$var_cols
$method
$out_file1
$out_file2
</command>
<inputs>
<param format="tabular" name="input1" type="data" label="Select data" help="Query missing? See TIP below."/>
<param name="var_cols" label="Select columns containing input variables " type="data_column" data_ref="input1" numerical="True" multiple="true" >
<validator type="no_options" message="Please select at least one column."/>
</param>
<param name="method" type="select" label="Method" help="The correlation matrix can only be used if there are no constant variables">
<option value="cor" selected="true">Correlation</option>
<option value="cov">Covariance</option>
</param>
</inputs>
<outputs>
<data format="input" name="out_file1" metadata_source="input1" />
<data format="pdf" name="out_file2" />
</outputs>
<requirements>
<requirement type="python-module">rpy</requirement>
</requirements>
<tests>
<test>
<param name="input1" value="iris.tabular"/>
<param name="var_cols" value="1,2,3,4"/>
<param name="method" value="cor"/>
<output name="out_file1" file="pca_out1.tabular"/>
<output name="out_file2" file="pca_out2.pdf"/>
</test>
<test>
<param name="input1" value="iris.tabular"/>
<param name="var_cols" value="1,2,3,4"/>
<param name="method" value="cov"/>
<output name="out_file1" file="pca_out3.tabular"/>
<output name="out_file2" file="pca_out4.pdf"/>
</test>
</tests>
<help>
.. class:: infomark
**TIP:** If your data is not TAB delimited, use *Edit Queries-&gt;Convert characters*
-----
.. class:: infomark
**What it does**
This tool uses the 'princomp' function from R statistical package to perform Principal Component Analysis (PCA) on the given numeric input data. It outputs two files, one containing the summary statistics of PCA, and the other containing biplots of the observations and principal components.
*R Development Core Team (2009). R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. ISBN 3-900051-07-0, URL http://www.R-project.org.*
-----
.. class:: warningmark
**Note**
- This tool currently treats all variables as continuous numeric variables. Running the tool on categorical variables might result in incorrect results. Rows containing non-numeric (or missing) data in any of the chosen columns will be skipped from the analysis.
- The summary statistics in the output are described below:
- Std. deviation: Standard deviations of the principal components
- Loadings: a list of eigen-vectors
- Scores: Scores of the input data on the principal components
</help>
</tool>
+5 -3
View File
@@ -11,8 +11,8 @@
</command>
<inputs>
<param format="tabular" name="input1" type="data" label="Select data" help="Query missing? See TIP below."/>
<param name="response_col" label="Response column (Y)" type="data_column" data_ref="input1" />
<param name="predictor_cols" label="Predictor columns (X)" type="data_column" data_ref="input1" multiple="true" >
<param name="response_col" label="Response column (Y)" type="data_column" data_ref="input1" numerical="True"/>
<param name="predictor_cols" label="Predictor columns (X)" type="data_column" data_ref="input1" numerical="True" multiple="true" >
<validator type="no_options" message="Please select at least one column."/>
</param>
</inputs>
@@ -47,13 +47,15 @@
This tool uses the 'lm' function from R statistical package to perform linear regression on the input data. It outputs two files, one containing the summary statistics of the performed regression, and the other containing diagnostic plots to check whether model assumptions are satisfied.
*R Development Core Team (2009). R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. ISBN 3-900051-07-0, URL http://www.R-project.org.*
-----
.. class:: warningmark
**Note**
- This tool currently treats all predictor and response variables as continuous variables. Running the tool on categorical variables might result in incorrect results.
- This tool currently treats all predictor and response variables as continuous numeric variables. Running the tool on categorical variables might result in incorrect results.
- Rows containing non-numeric (or missing) data in any of the chosen columns will be skipped from the analysis.