annotate prediction.R @ 100:09827a9f75a6 draft

Uploaded

########################################################
#
# creation date : 26/01/16
# last modification : 02/06/16
# author : Dr Nicolas Beaume
# owner : IRRI
#
########################################################

library(rrBLUP)
suppressWarnings(suppressMessages(library(randomForest)))
library(e1071)
suppressWarnings(suppressMessages(library(glmnet)))
library(methods)


############################ main #############################
classifierNames <- c("list", "randomForest", "svm", "glmnet")
# load argument
cmd <- commandArgs(trailingOnly = T)
source(cmd[1])
# load data
con = file(genotype)
genotype <- readLines(con = con, n = 1, ok=T)
close(con)
genotype <- read.table(genotype, sep="\t", h=T)
con = file(model)
model <- readLines(con = con, n = 1, ok=T)
close(con)
model <- readRDS(model)
# check if the classifier name is valid
if(all(is.na(match(class(model), classifierNames)))) {
  stop(paste(class(model), "is not recognized as a valid model. Valid models are : ", classifierNames))
}
# run prediction according to the classifier
# rrBLUP prediction
if(any(class(model) == "list")) {
  predictions <- as.matrix(genotype) %*% as.matrix(model$u)
  predictions <- predictions[,1]+model$beta
  predictions <- data.frame(lines=rownames(genotype), predictions=predictions)
# LASSO prediction  
} else if(any(class(model) == "glmnet")) {
  predictions <- predict(model, as.matrix(genotype), type = "response")
  predictions <- data.frame(lines=rownames(genotype), predictions=predictions)
# SVM or RandomForest prediction (predict is a wrapper for many machine learning function)
} else {
  predictions <- predict(model, genotype)
  predictions <- data.frame(lines=names(predictions), predictions=predictions)
}
# save results
write.table(predictions, file=out, sep="\t", row.names = F)