process output

R/analysis/jags.run.R

@@ -4,7 +4,6 @@
 source('R/analysis/lmer.run.R')
 #source('R/analysis/stan.run.R')
 
-### TODO CHANGE TO FIT WITH NEW LMER DATA AND
 
 ### FUNCTION TO GENERATE THE GOOD FORMAT OF DATA
 
@@ -22,36 +21,36 @@ if(grepl('sigma',names(pars)[i])){
  }else{
   if(grepl('param',names(pars)[i])){
       if(grepl('logD', names(pars)[i])){
-    a.t <- 2/3+(chain-2.5)*0.15
+    a.t <- 2/3+(chain-sum(chains.v)/length(chains.v))*0.15
     a <- rnorm(list.jags[[pars[i]]],mean = a.t, sd = 0.1)
       }else{
         if(grepl('intercept', names(pars)[i])){
-    a.t <- 0.2+(chain-2.5)*0.1
+    a.t <- 0.2+(chain-sum(chains.v)/length(chains.v))*0.1
     a <- rnorm(list.jags[[pars[i]]],mean = a.t, sd = 0.1)
         }else{
          if(grepl('sumBn', names(pars)[i])){
-    a.t <- -0.25+(chain-2.5)*0.1
+    a.t <- -0.25+(chain-sum(chains.v)/length(chains.v))*0.1
     a <- rnorm(list.jags[[pars[i]]],mean = a.t, sd = 0.1)
-         }else{          
-          a.t <- 0+(chain-2.5)*0.15
+         }else{
+          a.t <- 0+(chain-sum(chains.v)/length(chains.v))*0.15
        a <- rnorm(list.jags[[pars[i]]],mean = a.t, sd = 0.1)
-        } 
+        }
        }
-      }    
+      }
   }else{
      if(grepl('logD', names(pars)[i])){
-    a <- 2/3+(chain-2.5)*0.15
+    a <- 2/3+(chain-sum(chains.v)/length(chains.v))*0.15
       }else{
         if(grepl('intercept', names(pars)[i])){
-    a <- 0.2+(chain-2.5)*0.1
+    a <- 0.2+(chain-sum(chains.v)/length(chains.v))*0.1
         }else{
           if(grepl('sumBn', names(pars)[i])){
-    a <- -0.25+(chain-2.5)*0.1
-          }else{                    
-    a <- 0+(chain-2.5)*0.15
+    a <- -0.25+(chain-sum(chains.v)/length(chains.v))*0.1
+          }else{
+    a <- 0+(chain-sum(chains.v)/length(chains.v))*0.15
          }
       }
-     }    
+     }
   }
  }
     return(a)
@@ -80,9 +79,9 @@ fun.call.jags.and.save <- function(jags.model,
                                    path.out,
                                    var.sample,
                                    iter = 5000,
-                                   warmup = 500,
+                                   warmup = 2000,
                                    chains = 3,
-                                   thin = 50,
+                                   thin = 5,
                                    method = 'rjags'){
 require(runjags)
  start <-  Sys.time()
@@ -151,7 +150,6 @@ print('OK')
                                select.set. = select.set,
                                sample.vec.TF = TRUE,
                                merge.biomes.TF = merge.biomes.TF)
-    browser()
     cat("Ok data with Nobs", jags.list[[1]]$N_indiv, 'run real',
          "\n")
     res <- fun.call.jags.and.save(jags.model = model,

launch.cluster/launch_all_jags.bash

@@ -5,9 +5,9 @@ export LD_LIBRARY_PATH=/usr/lib64/R/library
 
 mkdir -p Rscript_temp
 
-iter=4000
-warmup=1000
-thin=4
+iter=6000
+warmup=2000
+thin=18
 # test parameter
 nbargs=$#
 echo "number of arguments="$nbargs

launch.cluster/launch_all_jags.init.bash

@@ -21,10 +21,10 @@ fi
 
 
 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(c(model.files.jags.Tf.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 = TRUE);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(c(model.files.jags.Tf.1), run.jags.b,trait = $trait,data.type='simple', sample.size = $samplesize, var.sample = 'wwf',iter = $iter, warmup = $warmup, chains = 3, thin = $thin, init.TF = TRUE);print('done')\"" > Rscript_temp/speciesjags${trait}.sh
 	qsub Rscript_temp/speciesjags${trait}.sh -d ~/trait.competition.workshop -l nodes=1:ppn=1,mem=8gb -N "jags${trait}" -q opt32G -j oe
 
 

scripts/analysis/lmer.all.no.log.output.R

@@ -6,17 +6,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.0.rds',
-                   models = c(model.files.lmer.Tf.0),
-                   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.simple.set.0.rds',
+##                    models = c(model.files.lmer.Tf.0),
+##                    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.simple.set.rds',
-                   models = c(model.files.lmer.Tf.1),
-                   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.simple.set.rds',
+##                    models = c(model.files.lmer.Tf.1),
+##                    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.simple.set.TP.rds',