|
16
|
1 ########################################################
|
|
|
2 #
|
|
|
3 # creation date : 07/01/16
|
|
|
4 # last modification : 03/07/16
|
|
|
5 # author : Dr Nicolas Beaume
|
|
|
6 # owner : IRRI
|
|
|
7 #
|
|
|
8 ########################################################
|
|
|
9 log <- file(paste(getwd(), "log_SVM.txt", sep="/"), open = "wt")
|
|
|
10 sink(file = log, type="message")
|
|
|
11 library("e1071")
|
|
|
12 ############################ helper functions #######################
|
|
|
13 svmModel <- function(train, target, kernel="radial", g=NULL, c=NULL, coef=NULL, d=NULL) {
|
|
|
14 # tuning parameters then train
|
|
|
15 model <- NULL
|
|
|
16 if(is.null(g)){g <- 10^(-6:0)}
|
|
|
17 if(is.null(c)){c <- 10^(0:2)}
|
|
|
18 switch(kernel,
|
|
|
19 sigmoid={
|
|
|
20 tune <- tune.svm(train, target, gamma = , cost = 10^(0:2), kernel="sigmoid");
|
|
|
21 g <- tune$best.parameters[[1]];
|
|
|
22 c <- tune$best.parameters[[2]];
|
|
|
23 model <- svm(x=train, y=target, gamma = g, cost = c, kernel = "sigmoid")},
|
|
|
24 linear={
|
|
|
25 tune <- tune.svm(train, target, cost = c, kernel="linear");
|
|
|
26 c <- tune$best.parameters[[2]];
|
|
|
27 model <- svm(x=train, y=target, gamma = g, cost = c, kernel = "linear")},
|
|
|
28 polynomial={
|
|
|
29 if(is.null(coef)){coef <- 0:3};
|
|
|
30 if(is.null(d)){d <- 0:4};
|
|
|
31 tune <- tune.svm(train, target, gamma = g, cost = c, degree = d, coef0 = coef, kernel="polynomial");
|
|
|
32 d <- tune$best.parameters[[1]];
|
|
|
33 g <- tune$best.parameters[[2]];
|
|
|
34 coef <- tune$best.parameters[[3]];
|
|
|
35 c <- tune$best.parameters[[4]];
|
|
|
36 model <- svm(x=train, y=target, gamma = g, cost = c, kernel = "polynomial", degree = d, coef0 = coef)},
|
|
|
37 {
|
|
|
38 tune <- tune.svm(train, target, gamma = g, cost = c, kernel="radial");
|
|
|
39 g <- tune$best.parameters[[1]];
|
|
|
40 c <- tune$best.parameters[[2]];
|
|
|
41 model <- svm(x=train, y=target, gamma = g, cost = c, kernel = "radial")}
|
|
|
42 )
|
|
|
43 return(model)
|
|
|
44 }
|
|
|
45 ################################## main function ###########################
|
|
|
46 svmSelection <- function(genotype, evaluation = T, outFile, folds, kernel="radial", g=NULL, c=NULL, coef=NULL, d=NULL) {
|
|
|
47 # build model
|
|
|
48 labelIndex <- match("label", colnames(genotype))
|
|
|
49 if(evaluation) {
|
|
|
50 prediction <- NULL
|
|
|
51 for(i in 1:length(folds)) {
|
|
|
52 test <- folds[[i]]
|
|
|
53 train <- unlist(folds[-i])
|
|
|
54 svm.fit <- svmModel(train = genotype[train,-labelIndex], target = genotype[train,labelIndex], kernel=kernel, g=g, c=c, coef=coef, d=d)
|
|
|
55 prediction <- c(prediction, list(predict(svm.fit, genotype[test,-labelIndex])))
|
|
|
56 }
|
|
|
57 saveRDS(prediction, file=paste(outFile, ".rds", sep = ""))
|
|
|
58 } else {
|
|
|
59 model <- svmModel(train = genotype[,-labelIndex], target = genotype[,labelIndex], kernel=kernel, g=g, c=c, coef=coef, d=d)
|
|
|
60 saveRDS(model, file=paste(outFile, ".rds", sep = ""))
|
|
|
61 }
|
|
|
62 }
|
|
|
63
|
|
|
64 ############################ main #############################
|
|
|
65
|
|
|
66 cmd <- commandArgs(T)
|
|
|
67 source(cmd[1])
|
|
|
68 if(as.numeric(g) == -1) {g <- NULL}
|
|
|
69 if(as.numeric(c) == -1) {c <- NULL}
|
|
|
70 if(as.numeric(coef) == -1) {coef <- NULL}
|
|
|
71 if(as.numeric(d) == -1) {d <- NULL}
|
|
|
72 # check if evaluation is required
|
|
|
73 evaluation <- F
|
|
|
74 if(as.integer(doEvaluation) == 1) {
|
|
|
75 evaluation <- T
|
|
|
76 con = file(folds)
|
|
|
77 folds <- readLines(con = con, n = 1, ok=T)
|
|
|
78 close(con)
|
|
|
79 folds <- readRDS(folds)
|
|
|
80 }
|
|
|
81 # load genotype and phenotype
|
|
|
82 con = file(genotype)
|
|
|
83 genotype <- readLines(con = con, n = 1, ok=T)
|
|
|
84 close(con)
|
|
|
85 genotype <- read.table(genotype, sep="\t", h=T)
|
|
|
86 # phenotype is written as a table (in columns) but it must be sent as a vector for mixed.solve
|
|
|
87 phenotype <- read.table(phenotype, sep="\t", h=T)[,1]
|
|
|
88 # run !
|
|
|
89 svmSelection(genotype = data.frame(genotype, label=phenotype, check.names = F, stringsAsFactors = F),
|
|
|
90 evaluation = evaluation, outFile = out, folds = folds, g=g, c=c, coef=coef, d=d, kernel=kernel)
|
|
|
91 cat(paste(paste(out, ".rds", sep = ""), "\n", sep=""))
|
