Commit 0c2b3aa2 authored by Georges Kunstler's avatar Georges Kunstler
Browse files

new process with dplyr and US new data added with changes to makefile to R CMD...

new process with dplyr and US new data added with changes to makefile to R CMD BATCH instead of Rscript
parent 496bf320
......@@ -43,30 +43,30 @@ $(D2)/TRY/data.TRY.std.rds:
GLOBAL: $(D3)/Done.txt
$(D3)/Done.txt: scripts/process.data/remove.out.R scripts/process.data/summarise.data.R $(D3)/Done.t.txt
R CMD BATCH $<
R CMD BATCH scripts/process.data/summarise.data.R
R CMD BATCH scripts/process.data/plot.data.all.R
R CMD BATCH $< ; rm remove.out.Rout
R CMD BATCH scripts/process.data/summarise.data.R ; rm summarise.data.Rout
R CMD BATCH scripts/process.data/plot.data.all.R ; rm plot.data.all.Rout
$(D3)/Done.t.txt: scripts/process.data/merge.all.processed.data.R $(sites)
R CMD BATCH $<
R CMD BATCH $< ; rm merge.all.processed.data.Rout
#-------------------------------------------------------
TEST.TREE: $(D4)/Done.txt tree.all.sites
$(D4)/Done.txt: scripts/format.data/test.tree.R $(D3tree)
R CMD BATCH $<
R CMD BATCH $< ; rm test.tree.Rout; convert figs/test.format.tree/*.png figs/test.format.tree/all.tree.pdf
#-------------------------------------------------------
TEST.TRAITS: $(D5)/Done.txt #traits.all.sites
$(D5)/Done.txt: scripts/find.trait/test.traits.R $(D3traits)
R CMD BATCH $<
R CMD BATCH $< ; rm test.traits.Rout; convert figs/test.traits/*.pdf figs/test.traits/all.traits.pdf
#-------------------------------------------------------
TEST.CWM: $(D6)/Done.txt
$(D6)/Done.txt: scripts/process.data/test.tree.CWM.R $(D3Done)
R CMD BATCH $<
R CMD BATCH $< ; rm test.tree.CWM.Rout ; convert figs/test.processed/*.p* figs/test.processed/all.CWM.pdf
#-------------------------------------------------------
......
......@@ -80,3 +80,28 @@ names(out) <- lapply(lapply(files,files.details.all),
out <- out[!unlist(lapply(out,FUN=function(x) is.null(x$lmer.summary)))]
saveRDS(out,file='output/list.lmer.out.all.NA.no.log.rds')
## LOOP OVER FILES AND NOT SCAN ALL FILES NO SAMPLE
traits <- c("SLA", "Wood.density", "Max.height", "Leaf.N", "Seed.mass")
type.filling <- 'species'
files <- c()
for (trait in traits){
for(model in c(model.files.lmer.Tf.1,
model.files.lmer.Tf.2)){
source(model, local = TRUE)
model.obj <- load.model()
pathout <- output.dir('lmer', model.obj$name, trait, 'all.no.log',
type.filling=type.filling)
files <- c(files,file.path(pathout,"biomes.results.nolog.all.rds"))
}
}
out <- lapply(files, summarise.lmer.output.all.list)
names(out) <- lapply(lapply(files,files.details.all),
function(x) paste(as.vector(x[names(x) != 'file']),
collapse="_"))
### remove missing
out <- out[!unlist(lapply(out,FUN=function(x) is.null(x$lmer.summary)))]
saveRDS(out,file='output/list.lmer.out.all.biomes.no.log.rds')
......@@ -130,55 +130,6 @@ fun.test.set.B <- function(set, filedir) {
##### function to plot scatter with marginal hist
## from http://www.r-bloggers.com/example-10-3-enhanced-scatterplot-with-marginal-histograms/
scatterhist <- function(x, y, xlab = "", ylab = "", plottitle = "",
xsize = 1, cleanup = TRUE, ...){
# save the old graphics settings-- they may be needed
def.par <- par(no.readonly = TRUE)
zones <- matrix(c(1, 1, 1, 0, 5, 0, 2, 6, 4, 0, 3, 0),
ncol = 3, byrow = TRUE)
layout(zones, widths = c(1, 3, 1), heights = c(3, 3, 6, 2))
#layout(zones)
# tuning to plot histograms nicely
xhist <- hist(x, plot = FALSE)
yhist <- hist(y, plot = FALSE)
top <- max(c(xhist$counts, yhist$counts))
# for all three titles:
# drop the axis titles and omit boxes, set up margins
#par(xaxt = "n", yaxt = "n", bty = "n", mar = c(.3, 2, .3, 0) +.05)
par(xaxt = "n", yaxt = "n", bty = "n")
# fig 1 from the layout
plot(x = 1, y = 1, type = "n", ylim = c(-1, 1), xlim = c(-1, 1))
text(0, 0, paste(plottitle), cex = 2)
# fig 2
plot(x = 1, y = 1, type = "n", ylim = c(-1, 1), xlim = c(-1, 1))
text(0, 0, paste(ylab), cex = 1.5, srt = 90)
# fig 3
plot(x = 1, y = 1, type = "n", ylim = c(-1, 1), xlim = c(-1, 1))
text(0, 0, paste(xlab), cex = 1.5)
# fig 4, the first histogram, needs different margins
# no margin on the left
par(mar = c(2, 0, 1, 1))
barplot(yhist$counts, axes = FALSE, xlim = c(0, top),
space = 0, horiz = TRUE)
# fig 5, other histogram needs no margin on the bottom
par(mar = c(0, 2, 1, 1))
barplot(xhist$counts, axes = FALSE, ylim = c(0, top), space = 0)
# fig 6, finally, the scatterplot-- needs regular axes, different margins
par(mar = c(2, 2, .5, .5), xaxt = "s", yaxt = "s", bty = "n")
# this color allows traparency & overplotting-- useful if a lot of points
#plot(x, y , pch = 19, col = "#00000022", cex = xsize, ...)
plot(x, y , pch = 19, col = "#00000022", cex = xsize)
# reset the graphics, if desired
if(cleanup) {par(def.par)}
}
## Alternative method for scatter plot with marginal histogram; old code from FH
scatterhist2 <- function (x, y, hist.col = gray(0.95), ...) {
......@@ -209,9 +160,7 @@ fun.plot.XY.hist <- function(get.tree, set, var.1, var.2) {
png(file = paste("figs/test.format.tree/tree", var.1,
var.2, set, "png", sep = "."),
width = 880, height = 880)
#scatterhist(get.tree[[var.1]] , get.tree[[var.2]],
# xlab = var.1, ylab = var.2, plottitle=set, xsize=1, cleanup=TRUE)
dfxy = data.frame(get.tree[[var.1]], get.tree[[var.2]])
dfxy = data.frame(get.tree[[var.1]], get.tree[[var.2]])
if(sum(is.na(apply(dfxy, 1, sum))) > 0){
dfxy <- dfxy[-which(is.na(apply(dfxy, 1, sum))), ]
}
......
......@@ -31,12 +31,12 @@ data.max.height <- data.frame(sp = data.max.height$SpCd,
write.csv(data.max.height, file = "output/formatted/US/max.height.csv")
#### READ DATA US
data.us <- read.csv("data/raw/US/FIA51_trees_w_supp.csv",
data.us <- read.csv("data/raw/US/FIA51_trees.csv",
header = TRUE, stringsAsFactors = FALSE)
#MCV remove trees <10 cm dbh
library(dplyr)
data.us <- filter(data.us, InitDbh>=10 & !is.na(InitDbh) )
data.us <- filter(data.us, InitDbh>=10 & !is.na(InitDbh) & IntervalYears > 0 )
## data.us <- filter(data.us, IntervalYears>=2 & !is.na(IntervalYears) )
data.us.plot <- read.csv("data/raw/US/PlotIDInfo3.csv",
......@@ -47,16 +47,25 @@ data.us <- left_join(data.us, data.us.plot, by = 'PlotID')
data.us <- mutate(data.us, PLOT_ID_C = PlotID)
data.us$PlotID <- NULL
## remove plot with harvesting -998
harvest.plots <- group_by(data.us, PLOT_ID_C) %>%
summarise(harvest = sum(FinalDbh == -998)>0) %>%
filter(harvest == TRUE) %>% dplyr::select(PLOT_ID_C)
data.us <- filter(data.us, !(PLOT_ID_C %in% harvest.plots$PLOT_ID_C))
## MASSAGE TRAIT DATA HEIGHT DATA FOR TREE MISSING BRING US DATA FOR HEIGHT OVER
## WHEN WE ANALYZE THAT DATASET LATER ON
## FORMAT INDIVIDUAL TREE DATA
## change unit and names of variables to be the same in all data for the tree
data.us$G <- 10 * (data.us$FinalDbh - data.us$InitDbh)/abs(data.us$IntervalYears)
data.us$G <- 10 * (data.us$FinalDbh - data.us$InitDbh)/(data.us$IntervalYears)
## diameter growth in mm per year
data.us$BA.G <- (pi*(data.us$FinalDbh/2)^2-pi*(data.us$InitDbh/2)^2)/
abs(data.us$IntervalYears) ## ba growth in cm^2/year
data.us$G[which(data.us$InitDbh == 0 | data.us$FinalDbh == -999)] <- NA
data.us$BA.G[which(data.us$InitDbh == 0 | data.us$FinalDbh == -999)] <- NA
(data.us$IntervalYears) ## ba growth in cm^2/year
data.us$G[which(data.us$InitDbh == 0 |
data.us$FinalDbh == -999)] <- NA
data.us$BA.G[which(data.us$InitDbh == 0 |
data.us$FinalDbh == -999)] <- NA
data.us$year <- data.us$IntervalYears ## number of year between measuremen
......
source('R/process.data/process-fun.R')
process_bigplot_dataset('BCI', Rlim=15,std.traits='log')
source('R/process.data/process-fun.R')
process_bigplot_dataset('BCI', Rlim=15,std.traits='no')
source('R/process.data/process-fun.R')
process_inventory_dataset('Canada',std.traits='log')
source('R/process.data/process-fun.R')
process_inventory_dataset('Canada',std.traits='no')
source('R/process.data/process-fun.R')
process_inventory_dataset('France',std.traits='log')
source('R/process.data/process-fun.R')
process_inventory_dataset('France',std.traits='no')
source('R/process.data/process-fun.R')
process_bigplot_dataset('Fushan', Rlim=15,std.traits='log')
source('R/process.data/process-fun.R')
process_bigplot_dataset('Fushan', Rlim=15,std.traits='no')
source('R/process.data/process-fun.R')
process_bigplot_dataset('Japan', Rlim=15,std.traits='log')
source('R/process.data/process-fun.R')
process_bigplot_dataset('Japan', Rlim=15,std.traits='no')
source('R/process.data/process-fun.R')
process_bigplot_dataset('Luquillo', Rlim=15,std.traits='log')
source('R/process.data/process-fun.R')
process_bigplot_dataset('Luquillo', Rlim=15,std.traits='no')
source('R/process.data/process-fun.R')
process_bigplot_dataset('Mbaiki', Rlim=15,std.traits='log')
source('R/process.data/process-fun.R')
process_bigplot_dataset('Mbaiki', Rlim=15,std.traits='no')
source('R/process.data/process-fun.R')
process_inventory_dataset('NSW',std.traits='log')
source('R/process.data/process-fun.R')
process_inventory_dataset('NSW',std.traits='no')
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