| Previous changeset 60:50894fc3c231 (2016-10-26) Next changeset 62:613c81a4f7c0 (2016-10-26) |
|
Commit message:
Deleted selected files |
|
removed:
aggregation.R aggregation.xml |
| b |
| diff -r 50894fc3c231 -r 339ae95a3d46 aggregation.R --- a/aggregation.R Wed Oct 26 18:22:03 2016 -0400 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 |
| [ |
| b'@@ -1,320 +0,0 @@\n-########################################################\n-#\n-# creation date : 25/10/16\n-# last modification : 25/10/16\n-# author : Dr Nicolas Beaume\n-#\n-########################################################\n-\n-suppressWarnings(suppressMessages(library(GA)))\n-library("miscTools")\n-library(rpart)\n-suppressWarnings(suppressMessages(library(randomForest)))\n-library(e1071)\n-suppressWarnings(suppressMessages(library(glmnet)))\n-options(warn=-1)\n-############################ helper functions #######################\n-\n-##### Genetic algorithm\n-\n-# compute r2 by computing the classic formula\n-# compare the sum of square difference from target to prediciton\n-# to the sum of square difference from target to the mean of the target\n-r2 <- function(target, prediction) {\n- sst <- sum((target-mean(target))^2)\n- ssr <- sum((target-prediction)^2)\n- return(1-ssr/sst)\n-}\n-\n-optimizeOneIndividual <- function(values, trueValue) {\n- # change the value into a function\n- f <- function(w) {sum(values * w/sum(w))}\n- fitness <- function(x) {1/abs(trueValue-f(x))}\n- resp <- ga(type = "real-valued", fitness = fitness, min = rep(0, length(values)), max = rep(1, length(values)), \n- maxiter = 1000, monitor = NULL, keepBest = T)\n- resp@solution <- resp@solution/sum(resp@solution)\n- return(resp)\n-}\n-\n-optimizeWeight <- function(values, trueValue, n=1000) {\n- fitnessAll <- function(w) {\n- predicted <- apply(values, 1, weightedPrediction.vec, w)\n- return(mean(r2(trueValue, predicted)))\n- #return(mean(1/abs(trueValue-predicted)))\n- }\n- resp <- ga(type = "real-valued", fitness = fitnessAll, min = rep(0, ncol(values)), max = rep(1, ncol(values)), \n- maxiter = n, monitor = NULL, keepBest = T)\n- resp@solution <- resp@solution/sum(resp@solution)\n- return(resp)\n-}\n-\n-weightedPrediction <- function(classifiers, w) {\n- if(length(w) > ncol(classifiers)) {\n- warning("more weights than classifiers, extra weigths are ignored")\n- w <- w[1:ncol(classifiers)]\n- } else if(length(w) < ncol(classifiers)) {\n- warning("less weights than classifiers, extra classifiers are ignored")\n- classifiers <- classifiers[,1:length(w)]\n- }\n- prediction <- NULL\n- prediction <- c(prediction, apply(classifiers, 1, weightedPrediction.vec, w))\n- return(prediction)\n-}\n-\n-weightedPrediction.vec <- function(values, w) {\n- return(sum(values * w/sum(w)))\n-}\n-\n-##### meta-decision tree\n-\n-tuneTree <- function(data, target) {\n- data <- data.frame(data, target=target)\n- size <- nrow(data)\n- xerror <- NULL\n- split <- 1:ceiling(size/5)\n- leafSize <- 1:ceiling(size/10)\n- xerror <- matrix(rep(-1, length(split)*length(leafSize)), ncol=length(leafSize))\n- cp <- matrix(rep(-1, length(split)*length(leafSize)), ncol=length(leafSize))\n- for(i in 1:length(split)) {\n- for(j in 1:length(leafSize)) {\n- op <- list(minsplit=split[i], minbucket=leafSize[j])\n- tree <- rpart(target ~., data=data, control=op, method="anova")\n- xerror[i,j] <- tree$cptable[which.min(tree$cptable[,"xerror"]),"xerror"]\n- cp[i,j] <- tree$cptable[which.min(tree$cptable[,"xerror"]),"CP"]\n- }\n- }\n- index <- which(xerror==min(xerror), arr.ind = T)\n- op <- list(minsplit=split[index[1]], minbucket=leafSize[index[2]], cp=cp[index[1], index[2]])\n- return(op)\n-}\n-\n-###### meta-LASSO\n-# create fold by picking at random row indexes\n-createFolds <- function(nbObs, n) {\n- # pick indexes\n- index <- sample(1:n, size=nbObs, replace = T)\n- # populate folds\n- folds <- NULL\n- for(i in 1:n) {\n- folds <- c(folds, list(which(index==i)))\n- }\n- return(folds)\n-}\n-\n-searchParamLASSO <- function(genotype, phenotype, alpha=seq(0,1,0.1), n=7) {\n- folds <- createFolds(nrow(genotype), n = n)\n- acc <- NULL\n- indexAlpha <- 1\n- for(a in alpha) {\n- curAcc <- NULL\n- for(i in 1:n) {\n- train <- genotype[-folds[[i]],]\n- test <- genotype[folds[[i]],]\n- phenoTrain <- phenotype[-folds[[i]]]\n- phenoTest <- phenotype[folds[[i'..b'arget)\n- cv <- cv.glmnet(x=as.matrix(classifiers), y=target, alpha=alpha)\n- model <- glmnet(x=as.matrix(classifiers), y=target, alpha=alpha, lambda = cv$lambda.1se)\n- out <- paste(out, ".rds", sep = "")\n- saveRDS(model, out)\n- return(out)\n- } else {\n- return(predict(model, classifiers), out)\n- }\n-}\n-\n-aggregateRF <- function(classifiers, target=NULL, model=NULL, ntree=NULL, prediction=F, out) {\n- if(!prediction) {\n- ntree <- searchParamRF(genotype = classifiers, phenotype = target,\n- rangeNtree = seq(100, 1000, 100))\n- model <- randomForest(x=classifiers, y=target, ntree = ntree, mtry = ncol(classifiers))\n- out <- paste(out, ".rds", sep = "")\n- saveRDS(model, out)\n- return(out)\n- } else {\n- return(predict(model, classifiers), out)\n- }\n-}\n-\n-aggregateSVM <- function(classifiers, target=NULL, prediction=F, \n- model=NULL, c=NULL, g=NULL, d=NULL, coef=NULL, kernel="radial", out) {\n- if(!prediction) {\n- model <- searchParamSVM(train = classifiers, target = target, kernel = kernel)\n- out <- paste(out, ".rds", sep = "")\n- saveRDS(model, out)\n- return(out)\n- } else {\n- return(predict(model, classifiers), out)\n- }\n-}\n-\n-################################### main #############################\n-# # load argument\n-cmd <- commandArgs(T)\n-source(cmd[1])\n-# check if evaluation is required\n-evaluation <- F\n-if(as.integer(doEvaluation) == 1) {\n- evaluation <- T\n- con = file(folds)\n- folds <- readLines(con = con, n = 1, ok=T)\n- close(con)\n- folds <- readRDS(folds)\n-}\n-# check for model\n-if(model == "None") {\n- model <- NULL\n- prediction <- F\n-} else {\n- prediction <- T\n- con = file(model)\n- model <- readLines(con = con, n = 1, ok=T)\n- close(con)\n- model <- readRDS(model)\n-}\n-# load classifiers and phenotype\n-classifiers <- NULL\n-classifNames <- NULL\n-if(lassoPred !="None"){\n- classifiers <- c(classifiers, lassoPred)\n- classifNames <- c(classifNames, "lasso")\n-}\n-if(rrBLUPPred !="None"){\n- classifiers <- c(classifiers, rrBLUPPred)\n- classifNames <- c(classifNames, "rrBLUP")\n-}\n-if(rfPred !="None"){\n- classifiers <- c(classifiers, rfPred)\n- classifNames <- c(classifNames, "rf")\n-}\n-if(svmPred !="None"){\n- classifiers <- c(classifiers, svmPred)\n- classifNames <- c(classifNames, "svm")\n-}\n-classifPrediction <- NULL\n-for(classif in classifiers) {\n- classifPrediction <- c(classifPrediction, list(read.table(classif, sep="\\t", h=T)))\n-}\n-classifPrediction <- data.frame(classifPrediction)\n-colnames(classifPrediction) <- classifNames\n-# phenotype is written as a table (in columns) but it must be sent as a vector for mixed.solve\n-phenotype <- read.table(phenotype, sep="\\t", h=T)[,1]\n-# aggregate !\n-switch(method,\n- geneticMean={\n- res <- aggregateGeneticMean(classifiers = classifPrediction, target = phenotype,\n- out = out, prediction = prediction, model=model)\n- },\n- dt={\n- res <- aggregateDT(classifiers = classifPrediction, target = phenotype,\n- out = out, prediction = prediction, model=model)\n- },\n- lasso={\n- res <- aggregateLASSO(classifiers = data.matrix(classifPrediction), target = phenotype,\n- out = out, prediction = prediction, model=model)\n- },\n- rf={\n- res <- aggregateRF(classifiers = classifPrediction, target = phenotype,\n- out = out, prediction = prediction, model=model)\n- },\n- # svm\n- {res <- aggregateSVM(classifiers = classifPrediction, target = phenotype, kernel = kernel,\n- out = out, prediction = prediction, model = model)}\n- )\n-if(prediction) {\n- write.table(data.frame(lines=rownames(classifPrediction), res), paste(out,".csv", sep=""),\n- sep="\\t", row.names = F)\n-} else {\n- cat(paste(paste(out, ".rds", sep = ""), "\\n", sep=""))\n-}\n\\ No newline at end of file\n' |
| b |
| diff -r 50894fc3c231 -r 339ae95a3d46 aggregation.xml --- a/aggregation.xml Wed Oct 26 18:22:03 2016 -0400 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 |
| b |
| @@ -1,78 +0,0 @@ -<tool id="aggreg" name="aggregation" version="1.0.0"> - <description>predict phenotype by combining multiple classifiers</description> - <command interpreter="Rscript"> - aggregation.R $config > ${output1} - </command> - - <inputs> - <param name="lassoPred" type="data" optional="true" - label="lasso prediction" help="path to rds containing LASSO prediction" - /> - - <param name="rfPred" type="data" optional="true" - label="rf prediction" help="path to rds containing Random Forest prediction" - /> - - <param name="rrBLUPPred" type="data" optional="true" - label="rrBLUP prediction" help="path to rds containing rrBLUP prediction" - /> - - <param name="svmPred" type="data" optional="true" - label="SVM prediction" help="path to rds containing SVM prediction" - /> - - <param name="phenotype" type="data" - label="phenotype data" help=" a tabular datatype containing the phenotypes " - /> - - <param name="eval" type="integer" value="0" - label="do evaluation" help=" whether to produce a model or to use folds to evaluate the tool. 1 means the tool will be evaluate (and a folds argument is required) any other value produces a model " - /> - - <param name="folds" type="data" optional="true" - label="folds" help=" OPTIONAL ARGUMENT path to a folds file containing folds indexes in a R list called /folds/ such as produced by the folds tools in OGHMA suite. " - /> - - <param name="model" type="data" optional="true" - label="model" help= " a path to a file where the results (depending on the chosen mode) will be writen" - /> - <param name="out" type="text" - label="output path" help= " a path to a rds file" - /> - <param name="method" type="text" value="svm" - label="aggregation method" help= "choose among geneticMean, rrBLUP, lasso, rf or svm" - /> - <param name="kernel" type="text" value="linear" - label="kernel for SVM" help= "choose among linear, polynomial, radial, sigmoid" - /> - - </inputs> - - <configfiles> - <configfile name="config"> -## Desc: this file is sourced in encode wrapper script -## as means to pass all galaxy params to R -"${lassoPred}" -> lassoPred -"${rfPred}" -> rfPred -"${rrBLUPPred}" -> rrBLUPPred -"${svmPred}" -> svmPred -"${phenotype}" -> phenotype -"${model}" -> model -"${out}" -> out -"${eval}" -> evaluation -"${folds}" -> folds -"${method}" -> method -"${kernel}" -> kernel -"${eval}" -> doEvaluation - - </configfile> -</configfiles> - -<outputs> - <data format="tabular" name = "output1" label="aggregation output" /> -</outputs> - - <help> - - </help> - </tool> \ No newline at end of file |