From fb51848436a0d6a9f76c0160a927d274c60e084c Mon Sep 17 00:00:00 2001
From: Daniel Falster <daniel.falster@mq.edu.au>
Date: Tue, 17 Sep 2013 16:43:08 +1000
Subject: [PATCH] Tweak files for running jags locally
---
.gitignore | 1 +
R/analysis/{fun.analysis.R => JAGS.fun.R} | 34 +++++++++++++++---
R/analysis/run.local.R | 42 +++++++++++++----------
readme.md | 1 +
4 files changed, 54 insertions(+), 24 deletions(-)
rename R/analysis/{fun.analysis.R => JAGS.fun.R} (93%)
diff --git a/.gitignore b/.gitignore
index 4213df5..d652326 100644
--- a/.gitignore
+++ b/.gitignore
@@ -14,3 +14,4 @@ figs
*.rds
docs/metadata/sites
.Rproj.user
+GROWTH.MODEL.*
diff --git a/R/analysis/fun.analysis.R b/R/analysis/JAGS.fun.R
similarity index 93%
rename from R/analysis/fun.analysis.R
rename to R/analysis/JAGS.fun.R
index 892cc18..82ae982 100644
--- a/R/analysis/fun.analysis.R
+++ b/R/analysis/JAGS.fun.R
@@ -1,7 +1,8 @@
###########################
###########################
##### FUNCTIONS USED IN FORMATTING DATA FOR JAGS AND RUN ANALYSIS
-
+source("R/packages.R")
+check_packages(c("rjags", "R2jags"))
######
######
@@ -165,7 +166,8 @@ return(jags.data)
#########################################
### WRITE JAGS MODELS
-
+write.jags.model.HD <- function(file = "GROWTH.MODEL.LOGLIN.HD"){
+
### LOG LIN HD
GROWTH.MODEL.LOGLIN.HD <-
"################################################################################
@@ -233,7 +235,10 @@ sigma.alpha ~ dunif(0.000001,5)
} # End of the jags model
"
-cat(GROWTH.MODEL.LOGLIN.HD , file = "GROWTH.MODEL.LOGLIN.HD", sep=" ", fill = FALSE, labels = NULL, append = FALSE)
+cat(GROWTH.MODEL.LOGLIN.HD , file = file, sep=" ", fill = FALSE, labels = NULL, append = FALSE)
+}
+
+write.jags.model.AD <- function(file = "GROWTH.MODEL.LOGLIN.AD"){
### LOG LIN HD
GROWTH.MODEL.LOGLIN.AD <-
@@ -302,8 +307,8 @@ sigma.alpha ~ dunif(0.000001,5)
} # End of the jags model
"
-cat(GROWTH.MODEL.LOGLIN.AD , file = "GROWTH.MODEL.LOGLIN.AD", sep=" ", fill = FALSE, labels = NULL, append = FALSE)
-
+cat(GROWTH.MODEL.LOGLIN.AD , file = file, sep=" ", fill = FALSE, labels = NULL, append = FALSE)
+}
################################################
###########################
@@ -455,3 +460,22 @@ n.chain <- dim(jags.output$BUGSoutput$sims.array)[2]
length.chain <- dim(jags.output$BUGSoutput$sims.array)[1]
plot(x,y,col=rep(1:n.chain,each=length.chain),xlab=param.name.1,ylab=param.name.2)
}
+
+
+get.ecoregions <- function(site){
+ sub(paste("./output/processed/",site,"/",sep=""),"",list.dirs(paste("./output/processed/",site,sep="")))[-1]
+}
+
+run.jags.for.trait <- function(site, trait, A1, A2, ecoregions = get.ecoregions(site)){
+ for (e in ecoregions)
+ try(fun.run.4models.jags(site=site,ecoregion.name=e,min.obs=10,
+ name.trait=trait,
+ A1=A1,
+ A2=A2,
+ n.chains=3,
+ n.iter=3000,
+ n.burnin=1000,
+ n.thin=1,
+ sample.size=2000,
+ na.action.traits="fill with global mean"))
+}
diff --git a/R/analysis/run.local.R b/R/analysis/run.local.R
index a808bb9..2c3da7f 100644
--- a/R/analysis/run.local.R
+++ b/R/analysis/run.local.R
@@ -1,20 +1,24 @@
#####
-## run jags code
-
-source("./R/analysis/fun.analysis.R")
-site <- "NVS"
-
-ecoregion.code <- sub(paste("./output/processed/",site,"/",sep=""),"",list.dirs(paste("./output/processed/",site,sep="")))[-1]
-
-for (i in ecoregion.code[1]){
-fun.run.4models.jags(site=site,ecoregion.name=i,min.obs=10,
- name.trait="SLA",
- A1="MAT",
- A2="MAP",
- n.chains=3,
- n.iter=3000,
- n.burnin=1000,
- n.thin=1,
- sample.size=NA,
- na.action.traits="fill with global mean")
-}
+
+source("R/analysis/JAGS.fun.R")
+
+write.jags.model.HD()
+write.jags.model.AD()
+
+traits <- c("SLA", "Leaf.N","Seed.mass","SLA","Wood.density","Max.height")
+
+for(this.trait in traits)
+ run.jags.for.trait("France", trait = this.trait, A1 = "sgdd", A2 = "WB.s")
+
+for(this.trait in traits)
+ run.jags.for.trait("NVS", trait = this.trait, A1 = "MAT", A2 = "MAP")
+
+for(this.trait in traits)
+ run.jags.for.trait("US", trait = this.trait, A1 = "MAT", A2 = "MAP")
+
+for(this.trait in traits)
+ run.jags.for.trait("Spain", trait = this.trait, A1 = "MAT", A2 = "PP")
+
+for(this.trait in traits)
+ run.jags.for.trait("Canada", trait = this.trait, A1 = "MAT", A2 = "MAP")
+
diff --git a/readme.md b/readme.md
index b9e7b4e..3b09209 100644
--- a/readme.md
+++ b/readme.md
@@ -33,6 +33,7 @@ We are making use of the following tools:
- [*markdown*](http://daringfireball.net/projects/markdown/)
- The [*knitr*](http://cran.r-project.org/web/packages/knitr/index.html) package in R, for generating documents with R and markdown. For a n intro see [here](http://nicercode.github.io/guides/reports/)
- version control via git (see next section)
+- JAGS - Just another Gibbs Sampler. Install [here](http://sourceforge.net/projects/mcmc-jags/files/)
### Git
--
GitLab