Commit dedc93b2 authored by Georges Kunstler's avatar Georges Kunstler
Browse files

re run with A_EV A_D but error in TEST CWM

parent 204db846
#!/usr/bin/env Rscript
rm(list = ls())
source("R/analysis/lmer.output-fun.R")
source("R/utils/plot.R")
library(pander)
## load results all
list.all.results <- readRDS('output/list.lmer.out.nolog.all.rds')
names.param <- unique(unlist(lapply(list.all.results, function(x) names(x$lmer.summary$fixed.coeff.E))))
names.param2 <- names.param
names.param2[names.param2 == "sumBn:MAT"] <- "MAT:sumBn"
names.param2[names.param2 == "sumBn:MAP"] <- "MAP:sumBn"
names.param2[names.param2 == "Tf:MAP"] <- "MAP:Tf"
names.param2[names.param2 == "Tf:MAT"] <- "MAT:Tf"
names(names.param2) <-names.param
DF.global <- do.call( 'rbind', lapply(list.all.results,
fun.merge.results.global,
names.param = names.param2))
if(sum(DF.global$conv != 0) >0) stop( 'bad convergence')
## AIC table
for (i in c('Wood.density', 'SLA', 'Leaf.N', 'Max.height')){
DF.t <- DF.global[DF.global$trait == i, ]
print(i)
print(DF.t$model[DF.t$AIC == min(DF.t$AIC)])
res <- data.frame(model = DF.t$model,delta.aic = DF.t$AIC - min(DF.t$AIC))
print(pander(res))
}
## plot the parameters
pdf('figs/param.global.model.ER.AD.pdf')
par(mfrow = c(2,2))
for (i in c('Wood.density', 'SLA', 'Leaf.N', 'Max.height')){
DF.t <- DF.global[DF.global$trait == i, ]
print(i)
param.vec <- c("Tf", "logD","sumBn", "sumTnBn","sumTfBn", "sumTnTfBn.abs")
param.vec.std <- paste(param.vec, 'Std.Error', sep = '.')
plot(unlist(DF.t[DF.t$model == "lmer.LOGLIN.ER.AD.Tf", param.vec]),
xaxt = 'n', ylab = 'std parameters', xlab = NA,
main = i, ylim = c(-1, 1))
abline(h = 0)
axis(1, 1:length(param.vec), labels = param.vec, las = 3)
fun.plot.error.bar(1:length(param.vec),
DF.t[DF.t$model == "lmer.LOGLIN.ER.AD.Tf", param.vec],
DF.t[DF.t$model == "lmer.LOGLIN.ER.AD.Tf", param.vec.std])
}
dev.off()
for (i in c('Wood.density', 'SLA', 'Leaf.N', 'Max.height')){
DF.t <- DF.global[DF.global$trait == i, ]
print(i)
param.vec <- c("Tf", 'MAT.Tf', 'MAP.Tf', "logD","sumBn", 'MAT.sumBn', 'MAP.sumBn',
"sumTnBn", 'MAT.sumTnBn', 'MAP.sumTnBn', "sumTfBn", 'MAT.sumTfBn',
'MAP.sumTfBn', "sumTnTfBn.abs" , 'MAT.sumTnTfBn.abs',
'MAP.sumTnTfBn.abs')
param.vec.std <- paste(param.vec, 'Std.Error', sep = '.')
plot(unlist(DF.t[DF.t$model == "lmer.LOGLIN.ER.AD.Tf.MAT.MAP", param.vec]),
xaxt = 'n', ylab = 'std parameters', xlab = NA,
main = i, ylim = c(-1, 1))
abline(h = 0)
axis(1, 1:length(param.vec), labels = param.vec, las = 3)
fun.plot.error.bar(1:length(param.vec),
DF.t[DF.t$model == "lmer.LOGLIN.ER.AD.Tf.MAT.MAP", param.vec],
DF.t[DF.t$model == "lmer.LOGLIN.ER.AD.Tf.MAT.MAP", param.vec.std])
}
......@@ -551,36 +551,23 @@ fun.plot.panel.lmer.parameters.boxplot <- function(models, traits, DF.results,
## create a data base from the global data
fun.merge.results.global <- function(list.res,
names.param = c( "(Intercept)", "Tf",
"logD", "sumBn", "sumTnBn",
"sumTfBn", "sumTnTfBn.abs",
"sumTnTfBn.diff")){
df.details <- data.frame(list.res$files.details[1:4],
names.param ){
names.param.vec <- unique(names.param)
df.details <- data.frame(list.res$files.details[1:4],
list.res$lmer.summary[1:6])
dat.t <- data.frame(t(rep(NA, 3*length(names.param))))
names(dat.t) <- c(names.param, paste(names.param, "Std.Error", sep=".")
, paste(names.param, "VAR", sep="."))
dat.t[, match(names(list.res$lmer.summary$fixed.coeff.E), names(dat.t))] <-
list.res$lmer.summary$fixed.coeff.E
dat.t[, length(names.param)+match(names(list.res$lmer.summary$fixed.coeff.E), names(dat.t))] <-
list.res$lmer.summary$fixed.coeff.Std.Error
dat.t[, match(names(list.res$lmer.summary$fixed.var), names(dat.t))] <-
list.res$lmer.summary$fixed.var
dat.t <- data.frame(matrix(rep(NA, 3 * length(names.param.vec)), nrow = 1,
ncol = 3 * length(names.param.vec)))
names(dat.t) <- c(names.param.vec,
paste(names.param.vec, "Std.Error", sep="."),
paste(names.param.vec, "VAR", sep="."))
dat.t[, names.param2[ names(list.res$lmer.summary$fixed.coeff.E)]] <-
list.res$lmer.summary$fixed.coeff.E
dat.t[, paste(names.param2[ names(list.res$lmer.summary$fixed.coeff.E)], "Std.Error", sep=".")] <-
list.res$lmer.summary$fixed.coeff.Std.Error
dat.t[, paste(names.param2[ names(list.res$lmer.summary$fixed.coeff.E)], "VAR", sep=".")] <-
list.res$lmer.summary$fixed.var
res <- data.frame(df.details, dat.t)
## fun.rep.df <- function(x, df = df.details) df
## n_rep <- nrow(list.res$lmer.summary$ecocode.BLUP)
## dat.t <- data.frame(matrix(rep(NA, n_rep * length(names.param)), nrow = n_rep,
## ncol = length(names.param)))
## names(dat.t) <- c(names.param)
## dat.t[, match(names(list.res$lmer.summary$ecocode.BLUP), names(dat.t))] <-
## list.res$lmer.summary$ecocode.BLUP
## df.res <- data.frame(ecocode = rownames(list.res$lmer.summary$ecocode.BLUP),
## do.call( "rbind", lapply(1:n_rep, fun.rep.df)))#,
## ## dat.t)
return(res)
}
......
#### function to analyse lmer output for GLOBAL ANALYSIS
library(lme4)
read.lmer.output <- function(file.name){
tryCatch(readRDS(file.name), error=function(cond)return(NULL)) # Choose a return value in case of error
}
summarise.lmer.all.output <- function(x){
list( nobs = nobs(x),
R2m =Rsquared.glmm.lmer(x)$Marginal,
R2c =Rsquared.glmm.lmer(x)$Conditional,
AIC = AIC(x),
deviance = deviance(x),
conv=x@optinfo$conv,
effect.response.var=variance.fixed.glmm.lmer.effect.and.response(x),
fixed.coeff.E=fixef(x),
fixed.coeff.Std.Error=sqrt(diag(vcov(x))),
fixed.var=variance.fixed.glmm.lmer(x),
set.BLUP=coef(x)$set.id)
}
summarise.lmer.output.all.list <- function(f ){
output.lmer <- read.lmer.output(f)
if(!is.null(output.lmer)){
res <- list(files.details=files.details.all(f),
lmer.summary=summarise.lmer.output( output.lmer))
}else{
res <- NULL
}
return(res)
}
files.details.all <- function(x){
s <- data.frame(t(strsplit(x, "/", fixed=TRUE)[[1]]), stringsAsFactors= FALSE)
names(s) <- c("d1", "d2", "set", "trait", "filling", "model", "file" )
s[-(1:2)]
}
#### R squared functions
# Function rsquared.glmm requires models to be input as a list (can include fixed-
# effects only models,but not a good idea to mix models of class "mer" with models
# of class "lme") FROM http://jslefche.wordpress.com/2013/03/13/r2-for-linear-mixed-effects-models/
Rsquared.glmm.lmer <- function(i){
# Get variance of fixed effects by multiplying coefficients by design matrix
VarF <- var(as.vector(fixef(i) %*% t(i@pp$X)))
# Get variance of random effects by extracting variance components
VarRand <- colSums(do.call(rbind,lapply(VarCorr(i),function(j) j[1])))
# Get residual variance
VarResid <- attr(VarCorr(i),"sc")^2
# Calculate marginal R-squared (fixed effects/total variance)
Rm <- VarF/(VarF+VarRand+VarResid)
# Calculate conditional R-squared (fixed effects+random effects/total variance)
Rc <- (VarF+VarRand)/(VarF+VarRand+VarResid)
# Bind R^2s into a matrix and return with AIC values
Rsquared.mat <- data.frame(Class=class(i),Family="Gaussian",Marginal=Rm,
Conditional=Rc,AIC=AIC(i))
return(Rsquared.mat)
}
variance.fixed.glmm.lmer <- function(i){
# Get variance of for each fixed effects by multiplying coefficients by design matrix
var.vec <- apply(fixef(i) * t(i@pp$X),MARGIN=1,var)
# Get variance of fixed effects by multiplying coefficients by design matrix
VarF <- var(as.vector(fixef(i) %*% t(i@pp$X)))
# Get variance of random effects by extracting variance components
VarRand <- colSums(do.call(rbind,lapply(VarCorr(i),function(j) j[1])))
# Get residual variance
VarResid <- attr(VarCorr(i),"sc")^2
var.vec <- var.vec/(VarF+VarRand+VarResid)
names(var.vec) <- paste(names(var.vec),"VAR",sep=".")
return(var.vec)
}
variance.fixed.glmm.lmer.effect.and.response <- function(i){
if(sum(c("sumTfBn","sumTnBn") %in% names(fixef(i)))==2){
# Get variance of for each fixed effects by multiplying coefficients by design matrix
var.vec <- var(as.vector(fixef(i)[c("sumTfBn","sumTnBn")] %*% t(i@pp$X[,c("sumTfBn","sumTnBn")])))
# Get variance of fixed effects by multiplying coefficients by design matrix
VarF <- var(as.vector(fixef(i) %*% t(i@pp$X)))
# Get variance of random effects by extracting variance components
VarRand <- colSums(do.call(rbind,lapply(VarCorr(i),function(j) j[1])))
# Get residual variance
VarResid <- attr(VarCorr(i),"sc")^2
var.vec <- var.vec/(VarF+VarRand+VarResid)
}else{
var.vec <- NA
}
names(var.vec) <- paste("effect.response","VAR",sep=".")
return(var.vec)
}
## create a data base from the global data
fun.merge.results.global <- function(list.res,
names.param = c( "(Intercept)", "Tf",
"logD", "sumBn", "sumTnBn",
"sumTfBn", "sumTnTfBn.abs",
"sumTnTfBn.diff")){
df.details <- data.frame(list.res$files.details[1:4],
list.res$lmer.summary[1:6])
fun.rep.df <- function(x, df = df.details) df
n_rep <- nrow(list.res$lmer.summary$ecocode.BLUP)
dat.t <- data.frame(matrix(rep(NA,n_rep * length(names.param)), nrow = n_rep,
ncol = length(names.param)))
names(dat.t) <- c(names.param)
dat.t[, match(names(list.res$lmer.summary$ecocode.BLUP), names(dat.t))] <-
list.res$lmer.summary$ecocode.BLUP
df.res <- data.frame(ecocode = rownames(list.res$lmer.summary$ecocode.BLUP),
do.call( "rbind", lapply(1:n_rep, fun.rep.df)),
dat.t)
return(df.res)
}
......@@ -118,9 +118,9 @@ data.tree.tot <- read.csv(file = file.path(base.dir, "data.all.csv"),
fun.data.for.lmer(data.tree.tot, trait, type.filling = type.filling)
}
fun.select.data.for.analysis <- function(data.tree, abs.CWM.tntf, perc.neigh,
BATOT, min.obs,
min.perc.neigh = 0.8, min.BA.G = -60,
fun.select.data.for.analysis <- function(data.tree, abs.CWM.tntf, perc.neigh.sp,
perc.neigh.gs, BATOT, min.obs,
min.perc.neigh.sp = 0.90, min.perc.neigh.gs = 0.98, min.BA.G = -60,
max.BA.G = 500){
## remove tree with NA
data.tree <- subset(data.tree, subset = (!is.na(data.tree[["BA.G"]])) &
......@@ -133,9 +133,12 @@ data.tree <- subset(data.tree, subset = data.tree[["BA.G"]]>min.BA.G &
select.temp <- (!is.na(data.tree[[abs.CWM.tntf]]) &
!is.na(data.tree[[BATOT]]))
data.tree <- data.tree[select.temp, ]
# select obs abov min perc.neigh
data.tree <- subset(data.tree, subset = data.tree[[perc.neigh]] > min.perc.neigh &
!is.na(data.tree[[perc.neigh]]))
# select obs abov min perc.neigh species
data.tree <- subset(data.tree, subset = data.tree[[perc.neigh.sp]] > min.perc.neigh.sp &
!is.na(data.tree[[perc.neigh.sp]]))
# select obs abov min perc.neigh genus
data.tree <- subset(data.tree, subset = data.tree[[perc.neigh.gs]] > min.perc.neigh.gs &
!is.na(data.tree[[perc.neigh.gs]]))
# select species with minimum obs
data.tree <- subset(data.tree, subset = data.tree[["sp"]] %in%
names(table(factor(data.tree[["sp"]])))[
......@@ -202,10 +205,12 @@ CWM.tn <- paste(trait, "CWM", 'fill', sep = ".")
abs.CWM.tntf <- paste(trait, "abs.CWM", 'fill', sep = ".")
tf <- paste(trait, "focal", sep = ".")
BATOT <- "BATOT"
perc.neigh <- paste(trait, "perc", type.filling, sep = ".")
perc.neigh.sp <- paste(trait, "perc", 'species', sep = ".")
perc.neigh.gs <- paste(trait, "perc", 'genus', sep = ".")
data.tree <- fun.select.data.for.analysis(data.tree,
abs.CWM.tntf,
perc.neigh,
perc.neigh,sp,
perc.neigh.gs,
BATOT,
min.obs)
#= DATA LIST FOR LMER
......
load.model <- function () {
list(name="lmer.LOGLIN.HD.Tf",
lmer.formula.tree.id=formula("logG~1+(1|set.id)+(1|species.id)+(1|plot.id)+Tf+logD+sumBn+sumBn:MAT+sumBn:MAP+sumTnTfBn.diff+sumTnTfBn.diff:MAT+sumTnTfBn.diff:MAP+(logD-1|species.id)+(Tf-1|set.id)+(sumBn-1|set.id)+(sumTnTfBn.diff-1|set.id)"))
}
......@@ -105,6 +105,10 @@ data.traits <- merge(data.traits,data.cat.extract[,c("sp","Phylo.group","Pheno.T
data.traits[data.traits$Latin_name %in% c('Vasconcellea cauliflora'),
'Phylo.group'] <- 'Angiosperm'
## create a vector with three categories C A_EV and A_D
## put EV_D in EV
data.traits$cat <- fun.three.cat.C.A_EV.A_D(data.traits$Pheno.T,
data.traits$Phylo.group)
###
......
......@@ -4,56 +4,57 @@
source("R/find.trait/trait-fun.R")
### read species names
data.tree <- read.csv("output/formatted/Canada/tree.csv",
data.tree <- read.csv("output/formatted/Canada/tree.csv",
stringsAsFactors = FALSE)
species.clean <- data.frame(sp=data.tree[!duplicated(data.tree[["sp"]]), "sp"],
Latin_name=data.tree[!duplicated(data.tree[["sp"]])
, "sp.name"],
Latin_name_syn=data.tree[!duplicated(
data.tree[["sp"]]),
"sp.name"],
stringsAsFactors = FALSE)
species.clean <-
data.frame(sp = data.tree[!duplicated(data.tree[["sp"]]), "sp"],
Latin_name = data.tree[!duplicated(data.tree[["sp"]])
, "sp.name"],
Latin_name_syn = data.tree[!duplicated(data.tree[["sp"]]),
"sp.name"],
stringsAsFactors = FALSE)
## read in data
data.TRY.std <- readRDS("output/formatted/TRY/data.TRY.std.rds")
## read us max height
max.height <- read.csv(file="output/formatted/US/max.height.csv",
max.height <- read.csv(file="output/formatted/US/max.height.csv",
stringsAsFactors = FALSE)
#manually assign max height to Alnus spp. (alder),
#manually assign max height to Alnus spp. (alder),
# from: http://plants.usda.gov/java/charProfile?symbol=ALINR
max.height[dim(max.height)[1]+1, ] <- c(dim(max.height)[1]+1, 350, 4.9, NA, NA)
max.height$sp <- paste("sp", max.height$sp, sep=".")
max.height$sp <- paste("sp", max.height$sp, sep = ".")
data.traits <- fun.extract.format.sp.traits.TRY(sp=species.clean[["sp"]],
sp.syno.table=species.clean,
data=data.TRY.std)
data.traits <- fun.extract.format.sp.traits.TRY(sp = species.clean[["sp"]],
sp.syno.table = species.clean,
data = data.TRY.std)
data.traits <- merge(data.traits, subset(max.height,
select=c("sp", "Max.height.mean",
"Max.height.sd")),
by="sp", all.x=TRUE, all.y=FALSE)
data.traits <- merge(data.traits, subset(max.height,
select = c("sp", "Max.height.mean",
"Max.height.sd")),
by = "sp", all.x = TRUE, all.y = FALSE)
data.traits$Max.height.genus <- FALSE
height.genus.DF <- do.call("rbind", lapply(data.traits$Latin_name,
fun.compute.mean.genus, data.traits,
height.genus.DF <- do.call("rbind", lapply(data.traits$Latin_name,
fun.compute.mean.genus, data.traits,
"Max.height.mean"))
data.traits[is.na(data.traits[["Max.height.mean"]]),
data.traits[is.na(data.traits[["Max.height.mean"]]),
c("Max.height.mean", "Max.height.sd", "Max.height.genus")] <-
height.genus.DF[is.na(data.traits[["Max.height.mean"]]), ]
height.genus.DF[is.na(data.traits[["Max.height.mean"]]), ]
#### GET THE ANGIO/CONIF AND EVERGREEN/DECIDUOUS
# read try categrocial data
try.cat <-
read.csv("data/raw/TRY/TRY_Categorical_Traits_Lookup_Table_2012_03_17_TestRelease.csv",
stringsAsFactors=FALSE, na.strings = "")
Pheno.Zanne <- read.csv("data/raw/ZanneNature/GlobalLeafPhenologyDatabase.csv",
stringsAsFactors=FALSE)
read.csv("data/raw/TRY/TRY_Categorical_Traits_Lookup_Table_2012_03_17_TestRelease.csv",
stringsAsFactors = FALSE, na.strings = "")
Pheno.Zanne <- read.csv("data/raw/ZanneNature/GlobalLeafPhenologyDatabase.csv",
stringsAsFactors = FALSE)
# extract
data.cat.extract <- do.call("rbind", lapply(data.traits$sp , fun.get.cat.var.from.try,
data.traits, try.cat, Pheno.Zanne))
data.cat.extract <- do.call("rbind", lapply(data.traits$sp ,
fun.get.cat.var.from.try,
data.traits, try.cat, Pheno.Zanne))
# change category
data.cat.extract <- fun.change.factor.pheno.try(data.cat.extract)
data.cat.extract <- fun.change.factor.angio.try(data.cat.extract)
......@@ -62,15 +63,21 @@ data.cat.extract <- fun.fill.pheno.try.with.zanne(data.cat.extract)
## fix pheno for species with issue
data.cat.extract[data.cat.extract$Latin_name %in%
c('Betula spp.', 'Crataegus spp.', 'Fraxinus spp.',
'Malus spp.', 'Amelanchier spp.', 'Alnus spp.',
'Tilia spp.', 'Ulmus spp.'), 'Pheno.T'] <- 'D'
'Malus spp.', 'Amelanchier spp.', 'Alnus spp.',
'Tilia spp.', 'Ulmus spp.', "Salix spp."), 'Pheno.T'] <- 'D'
data.traits <- merge(data.traits, data.cat.extract[, c("sp", "Phylo.group",
"Pheno.T", 'LeafType.T')],
by="sp")
data.traits <- merge(data.traits,
data.cat.extract[, c("sp", "Phylo.group",
"Pheno.T", 'LeafType.T')],
by = "sp")
###
write.csv(data.traits, file="output/formatted/Canada/traits.csv",
## create a vector with three categories C A_EV and A_D
## put EV_D in EV
data.traits$cat <- fun.three.cat.C.A_EV.A_D(data.traits$Pheno.T,
data.traits$Phylo.group)
###
write.csv(data.traits, file = "output/formatted/Canada/traits.csv",
row.names = FALSE)
......@@ -5,67 +5,91 @@ source("R/find.trait/trait-fun.R")
source("R/format.data/format-fun.R")
### read species names
species.clean <- (read.csv("data/raw/France/species.csv", stringsAsFactors = FALSE))
species.clean$Latin_name <- (gsub("_", " ", species.clean$Latin_name))
species.clean$Latin_name_syn<- (gsub("_", " ", species.clean$Latin_name_syn))
species.clean <- subset(species.clean,subset=!is.na(species.clean$Latin_name_syn))
species.clean$sp <- paste("sp",species.clean$code,sep=".")
species.clean <- read.csv("data/raw/France/species.csv",
stringsAsFactors = FALSE)
species.clean$Latin_name <- gsub("_", " ", species.clean$Latin_name)
species.clean$Latin_name_syn<- gsub("_", " ", species.clean$Latin_name_syn)
species.clean <- subset(species.clean,
subset = !is.na(species.clean$Latin_name_syn))
species.clean$sp <- paste("sp", species.clean$code, sep = ".")
species.clean$code <- NULL
## remove duplicated sp
species.clean <- species.clean[!duplicated(species.clean$sp),]
species.clean <- species.clean[!duplicated(species.clean$sp), ]
## read in data
data.TRY.std <- readRDS("output/formatted/TRY/data.TRY.std.rds")
max.height <- read.csv(file="output/formatted/France/max.height.csv", stringsAsFactors = FALSE)
max.height <- read.csv(file = "output/formatted/France/max.height.csv",
stringsAsFactors = FALSE)
max.height$Max.height.genus <- FALSE
### extract and add height
data.traits <- fun.extract.format.sp.traits.TRY(sp=species.clean[["sp"]],
sp.syno.table=species.clean,
data=data.TRY.std)
data.traits <- fun.extract.format.sp.traits.TRY(sp = species.clean[["sp"]],
sp.syno.table = species.clean,
data = data.TRY.std)
data.traits <- merge(data.traits,
subset(max.height,select=c("sp","Max.height.mean","Max.height.sd","Max.height.genus")),
by="sp",all.x=TRUE,all.y=FALSE)
height.genus.DF <- do.call("rbind",lapply(data.traits$Latin_name,
fun.compute.mean.genus,data.traits,"Max.height.mean"))
subset(max.height,
select = c("sp", "Max.height.mean",
"Max.height.sd", "Max.height.genus")),
by = "sp", all.x = TRUE, all.y = FALSE)
height.genus.DF <- do.call("rbind",
lapply(data.traits$Latin_name,
fun.compute.mean.genus,
data.traits,
"Max.height.mean"))
data.traits[is.na(data.traits[["Max.height.mean"]]),
c("Max.height.mean","Max.height.sd","Max.height.genus")] <-
height.genus.DF[is.na(data.traits[["Max.height.mean"]]),]
c("Max.height.mean", "Max.height.sd", "Max.height.genus")] <-
height.genus.DF[is.na(data.traits[["Max.height.mean"]]), ]
#### GET THE ANGIO/CONIF AND EVERGREEN/DECIDUOUS
# read try categrocial data
try.cat <- read.csv("data/raw/TRY/TRY_Categorical_Traits_Lookup_Table_2012_03_17_TestRelease.csv",
stringsAsFactors=FALSE,na.strings = "")
stringsAsFactors = FALSE, na.strings = "")
Pheno.Zanne <- read.csv("data/raw/ZanneNature/GlobalLeafPhenologyDatabase.csv",
stringsAsFactors=FALSE)
stringsAsFactors = FALSE)
# extract
data.cat.extract <- do.call("rbind",lapply(data.traits$sp ,fun.get.cat.var.from.try,
data.traits,try.cat,Pheno.Zanne))
data.cat.extract <- do.call("rbind", lapply(data.traits$sp ,
fun.get.cat.var.from.try,
data.traits, try.cat, Pheno.Zanne))
# change category
data.cat.extract <- fun.change.factor.pheno.try(data.cat.extract)
data.cat.extract <- fun.change.factor.angio.try(data.cat.extract)
data.cat.extract <- fun.fill.pheno.try.with.zanne(data.cat.extract)
# source Flore forestier Francaise Rameau
## fix pheno for species with issue
data.cat.extract[data.cat.extract$Latin_name %in% c('Populus sp','Tilia europaea',
'Prunus sp','Prunus dulcis',
'Sorbus latifolia','Pyrus spinosa',
'Sorbus semiincisa','Alnus alnobetula',
'Crataegus azarolus','Frangula dodonei',
'Tamarix africana'),'Pheno.T'] <- 'D'
data.cat.extract[data.cat.extract$Latin_name %in% c('Pinus nigra spp. salzmannii')
,'Pheno.T']<- 'EV'
data.traits <- merge(data.traits,data.cat.extract[,c("sp","Phylo.group","Pheno.T",'LeafType.T')],by="sp")
###
write.csv(data.traits,file="output/formatted/France/traits.csv",row.names = FALSE)
data.cat.extract[data.cat.extract$Latin_name %in% c('Populus sp',
'Tilia europaea',
'Prunus sp',
'Prunus dulcis',
'Sorbus latifolia',
'Pyrus spinosa',
'Sorbus semiincisa',
'Alnus alnobetula',
'Crataegus azarolus',
'Frangula dodonei',
'Tamarix africana'),
'Pheno.T'] <- 'D'
data.cat.extract[data.cat.extract$Latin_name %in%
c('Pinus nigra spp. salzmannii',
"Pinus nigra ssp salzmannii",
"Eucalyptus sp")
, 'Pheno.T']<- 'EV'
data.traits <- merge(data.traits, data.cat.extract[, c("sp", "Phylo.group",
"Pheno.T", 'LeafType.T')],
by = "sp")
## create a vector with three categories C A_EV and A_D
## put EV_D in EV
data.traits$cat <- fun.three.cat.C.A_EV.A_D(data.traits$Pheno.T,
data.traits$Phylo.group)
###
write.csv(data.traits, file = "output/formatted/France/traits.csv",
row.names = FALSE)