Commit 8593643c authored by kunstler's avatar kunstler
Browse files

launch jags not in parallel

parent a77de9cc
......@@ -6,76 +6,28 @@ read.jags.output <- function(file.name){
tryCatch(readRDS(file.name), error=function(cond)return(NULL)) # Choose a return value in case of error
}
summarise.jags.output <- function(x, probs = c(0.05, 0.1, 0.5, 0.9, 0.95)){
summarise.jags.output <- function(x,
param.to.save = c('intercept', 'mean_logD', 'mean_Tf',
'mean_sumBn', 'mean_sumTfBn', 'mean_sumTnBn',
'mean_sumTnTfBn_abs')){
list(niter = x$n.iter,
Rhat = max(x$BUGSoutput$summary[,"Rhat"]),
DIC = x$BUGSoutput$DIC,
deviance = quantile(x$BUGSoutput$sims.matrix[, "deviance"], probs = probs),
slope.trait = if( sum(colnames(x$BUGSoutput$sims.matrix)=="slope.trait")>0) {
quantile(x$BUGSoutput$sims.matrix[, "slope.trait"],probs = probs)}else{
rep(NA,length(probs))},
slope.trait.r = if( sum(colnames(x$BUGSoutput$sims.matrix)=="slope.trait.r")>0) {
quantile(x$BUGSoutput$sims.matrix[, "slope.trait.r"],probs = probs)}else{
rep(NA,length(probs))},
slope.trait.e = if( sum(colnames(x$BUGSoutput$sims.matrix)=="slope.trait.e")>0) {
quantile(x$BUGSoutput$sims.matrix[, "slope.trait.e"],probs = probs)}else{
rep(NA,length(probs))},
slope.trait.re = if( sum(colnames(x$BUGSoutput$sims.matrix)=="slope.trait.re")>0) {
quantile(x$BUGSoutput$sims.matrix[, "slope.trait.re"],probs = probs)}else{
rep(NA,length(probs))},
intercept.trait = if( sum(colnames(x$BUGSoutput$sims.matrix)=="intercept.trait")>0) {
quantile(x$BUGSoutput$sims.matrix[, "intercept.trait"],probs = probs)}else{
rep(NA,length(probs))})
summary = x$BUGSoutput$summary[param.to.save,],
DIC = x$BUGSoutput$DIC,
deviance = quantile(x$BUGSoutput$sims.matrix[, "deviance"], probs = c(0.025, 0.5, 0.975)))
}
summarise.jags.output.table <- function(f, probs = c(0.05, 0.1, 0.5, 0.9, 0.95) ){
data.frame(probs = probs,
files.details(f),
as.data.frame(summarise.jags.output(read.jags.output(f), probs = probs)),
row.names = NULL)
}
summarise.jags.output.table.old <- function(f, probs = c(0.05, 0.1, 0.5, 0.9, 0.95) ){
data.frame(probs = probs,
files.details.old(f),
as.data.frame(summarise.jags.output(read.jags.output(f), probs = probs)),
row.names = NULL)
list(file.names= files.details(f),
list.summary= summarise.jags.output(read.jags.output(f)))
}
files.details <- function(x){
s <- data.frame(t(strsplit(x, "/", fixed=TRUE)[[1]]), stringsAsFactors= FALSE)
names(s) <- c("d1", "d2", "model", "trait", "set", "ecocode", "file" )
s
}
files.details.old <- function(x){
s <- t(strsplit(x, "/", fixed=TRUE)[[1]])
if(grepl("SLA",s[,5])==1){
z <- "SLA"}
if(grepl("Leaf.N",s[,5])==1){
z <- "Leaf.N"}
if(grepl("Seed.mass",s[,5])==1){
z <- "Seed.mass"}
if(grepl("Wood.density",s[,5])==1){
z <- "Wood.density"}
if(grepl("Max.height",s[,5])==1){
z <- "Max.height"}
if(grepl("NULL",s[,5])==1){
u <- "NULL"}
if(grepl("AD",s[,5])==1){
u <- "AD"}
if(grepl("HD",s[,5])==1){
u <- "HD"}
if(grepl("RE",s[,5])==1){
u <- "RE"}
if(grepl("RE2",s[,5])==1){
u <- "RE2"}
s <- data.frame(s,trait=z,model=u)
names(s) <- c("d1", "d2", "set", "ecocode", "file", "traits", "model" )
s
}
plot.param.quant <- function(x,y.mat,param.name,ecocode){
......
......@@ -25,10 +25,10 @@ samplesize=$1
for trait in "'SLA'" "'Wood.density'" "'Max.height'"; do
echo "/usr/local/R/R-3.1.1/bin/Rscript -e \"source('R/analysis/jags.run.R'); run.multiple.model.for.set.one.trait(model.files.jags.Tf.1[1], run.jags.b, trait = $trait,data.type='simple', sample.size = NA, var.sample = 'wwf',iter = $iter, warmup = $warmup, chains = 3, thin = $thin, init.TF = FALSE, method = 'parallel');print('done')\"" > Rscript_temp/speciesjags${trait}.sh
echo "/usr/local/R/R-3.1.1/bin/Rscript -e \"source('R/analysis/jags.run.R'); run.multiple.model.for.set.one.trait(model.files.jags.Tf.1[1], run.jags.b, trait = $trait,data.type='simple', sample.size = NA, var.sample = 'wwf',iter = $iter, warmup = $warmup, chains = 3, thin = $thin, init.TF = FALSE);print('done')\"" > Rscript_temp/speciesjags${trait}.sh
qsub Rscript_temp/speciesjags${trait}.sh -d ~/trait.competition.workshop -l nodes=1:ppn=4,mem=16gb -N "jags${trait}" -q opt32G -j oe
echo "/usr/local/R/R-3.1.1/bin/Rscript -e \"source('R/analysis/jags.run.R'); run.multiple.model.for.set.one.trait(model.files.jags.Tf.1[2], run.jags.b, trait = $trait,data.type='simple', sample.size = NA, var.sample = 'wwf',iter = $iter, warmup = $warmup, chains = 3, thin = $thin, init.TF = FALSE, merge.biomes.TF = TRUE, method = 'parallel');print('done')\"" > Rscript_temp/speciesjags2${trait}.sh
echo "/usr/local/R/R-3.1.1/bin/Rscript -e \"source('R/analysis/jags.run.R'); run.multiple.model.for.set.one.trait(model.files.jags.Tf.1[2], run.jags.b, trait = $trait,data.type='simple', sample.size = NA, var.sample = 'wwf',iter = $iter, warmup = $warmup, chains = 3, thin = $thin, init.TF = FALSE, merge.biomes.TF = TRUE);print('done')\"" > Rscript_temp/speciesjags2${trait}.sh
qsub Rscript_temp/speciesjags2${trait}.sh -d ~/trait.competition.workshop -l nodes=1:ppn=4,mem=16gb -N "jags2${trait}" -q opt32G -j oe
done
......@@ -7,17 +7,17 @@ source("R/analysis/lmer.run.R")
## ## save list of all output for NA
## format.all.output.lmer(file.name = "NA.wwf.results.nolog.all.rds",
## list.file.name = 'list.lmer.out.all.NA.simple.set.rds',
## models = c(model.files.lmer.Tf.0, model.files.lmer.Tf.1,
## model.files.lmer.Tf.2, model.files.lmer.Tf.3),
## traits = c("SLA", "Wood.density", "Max.height"))
format.all.output.lmer(file.name = "NA.wwf.results.nolog.all.rds",
list.file.name = 'list.lmer.out.all.NA.intra.set.rds',
models = c(model.files.lmer.Tf.intra.1,
model.files.lmer.Tf.intra.2),
traits = c("SLA", "Wood.density", "Max.height"),
data.type = 'intra')
list.file.name = 'list.lmer.out.all.NA.simple.set.rds',
models = c(model.files.lmer.Tf.0, model.files.lmer.Tf.1,
model.files.lmer.Tf.2, model.files.lmer.Tf.3),
traits = c("SLA", "Wood.density", "Max.height"))
## format.all.output.lmer(file.name = "NA.wwf.results.nolog.all.rds",
## list.file.name = 'list.lmer.out.all.NA.intra.set.rds',
## models = c(model.files.lmer.Tf.intra.1,
## model.files.lmer.Tf.intra.2),
## traits = c("SLA", "Wood.density", "Max.height"),
## data.type = 'intra')
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment