An error occurred while loading the file. Please try again.
-
fhui28 authoredf85154e5
### MERGE BCI DATA
### Edited by FH
rm(list = ls()); source("./R/format.function.R"); library(reshape)
#########################
## READ DATA
####################
### read individuals tree data
## Requires careful formatting of 7 census datasets
## The raw data is such that, once a tree dies in census X, then it no longer exists in census X+1, X+2 etc...
data.bci1 <- read.table("./data/raw/DataBCI/census1/PlotsDataReport.txt",header=TRUE,stringsAsFactors=FALSE,sep = "\t")
big.bci <- NULL
for(k in 2:7) {
new.directory <- paste("./data/raw/DataBCI/census",k,"/PlotsDataReport.txt",sep="")
data.bci2 <- read.table(new.directory,header=TRUE,stringsAsFactors=FALSE,sep = "\t")
if(!is.null(big.bci)) {
sub.bci <- merge(data.bci1[,c(2:7,11,13)], data.frame(TreeID = data.bci2[["TreeID"]], DBH2 = data.bci2[["DBH"]], Date2 = data.bci2[["Date"]], dead = as.numeric(data.bci2[["Status"]] == "dead")),sort = T, by = "TreeID")
big.bci <- rbind(big.bci, sub.bci) }
if(is.null(big.bci)) {
big.bci <- merge(data.bci1[,c(2:7,11,13)], data.frame(TreeID = data.bci2[["TreeID"]], DBH2 = data.bci2[["DBH"]], Date2 = data.bci2[["Date"]], dead = as.numeric(data.bci2[["Status"]] == "dead")),
sort = T, by = "TreeID") }
data.bci1 <- data.bci2
cat("Census", k, "now included\n") }
rm(data.bci1, data.bci2, data.bci3, sub.bci)
big.bci <- big.bci[order(big.bci$TreeID),]; colnames(big.bci)[c(7:8)] <- c("DBH1","Date1")
data.bci <- big.bci; rm(big.bci)
### read species names
species.clean <- read.table("./data/raw/DataBCI/TaxonomyDataReport.txt",stringsAsFactors=FALSE, header = T, sep = "\t")
## Try to relate SpeciesID in species.clean species names in data.bci
data.bci$sp <- rep(NA,nrow(data.bci))
for(k in 1:nrow(species.clean)) {
latin.name <- apply(species.clean[k,2:3],1,paste,collapse=" ")[1]
data.bci$sp[which(data.bci$Latin == latin.name)] <- species.clean$SpeciesID[k] }
length(unique(data.bci$sp))
######################################
## MASSAGE TRAIT DATA
############################
## THERE ARE HEIGHT VARAIBLES IN THE TRAITS DATA, BUT WE WILL NEED TO DISCUSS WHICH ONES TO USE
##########################################
## FORMAT INDIVIDUAL TREE DATA
#############
data.bci$Date1 <- as.Date(data.bci$Date1); data.bci$Date2 <- as.Date(data.bci$Date2)
#data.bci$yr1 <- format(strptime(data.bci$Date1, format = "%Y-%m-%d"),"%Y")
#data.bci$yr2 <- format(strptime(data.bci$Date2, format = "%Y-%m-%d"),"%Y")
data.bci$year <- as.numeric(difftime(data.bci$Date2, data.bci$Date1, units = "weeks")/52) ## Not rounded
## change unit and names of variables to be the same in all data for the tree
data.bci$G <- 10*(data.bci$DBH1-data.bci$DBH1)/data.bci$year ## diameter growth in mm per year - BASED ON UNROUNDED YEARS
data.bci$D <- data.bci[["DBH1"]];
data.bci$plot <- data.bci[["Quadrat"]] ## plot code?
data.bci$htot <- rep(NA,length(data.bci[["G"]])) ## DOES IT COME FROM THE TRAITS DATA?
data.bci$sp.name <- data.bci$Latin
######################
## ECOREGION
###################
## bci has only 1 eco-region?
######################
## PERCENT DEAD
###################
## variable percent dead/cannot do with since dead variable is missing
## compute numer of dead per plot to remove plot with disturbance
function.perc.dead2 <- function(dead) { out <- sum(dead,na.rm=T)/length(dead[!is.na(dead)]); if(!is.finite(out)) out <- NA; return(out) }
7172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112
perc.dead <- tapply(data.bci[["dead"]],INDEX=data.bci[["plot"]],FUN=function.perc.dead2)
data.bci <- merge(data.bci,data.frame(plot=names(perc.dead),perc.dead=perc.dead), by = "plot", sort=FALSE)
###########################################################
### PLOT SELECTION FOR THE ANALYSIS
###################
## Remove data with dead == 1
table(data.bci$dead)
## Nothing to remove
vec.abio.var.names <- c("MAT","MAP") ## MISSING
vec.basic.var <- c("treeid","sp","sp.name","plot","D","G","dead","year","htot","x","y","perc.dead")
data.tree <- subset(data.bci,select=c(vec.basic.var,vec.abio.var.names))
##############################################
## COMPUTE MATRIX OF COMPETITION INDEX WITH SUM OF BA PER SPECIES IN EACH PLOT in m^2/ha without the target species
###########################
## DON'T KNOW SUBPLOT SIZE!
data.BA.SP <- BA.SP.FUN(id.tree=as.vector(data.bci[["TreeID"]]), diam=as.vector(data.bci[["D"]]),
sp=as.vector(data.bci[["sp"]]), id.plot=as.vector(data.bci[["plot"]]),
weights=1/(pi*(0.5*data.bci$D/100)^2), weight.full.plot=NA)
## change NA and <0 data for 0
data.BA.SP[is.na(data.BA.SP)] <- 0; data.BA.SP[,-1][data.BA.SP[,-1]<0] <- 0
### CHECK IF sp and sp name for column are the same
if(sum(!(names(data.BA.SP)[-1] %in% unique(data.bci[["sp"]]))) >0) stop("competition index sp name not the same as in data.tree")
#### compute BA tot for all competitors
BATOT.COMPET <- apply(data.BA.SP[,-1],1,sum,na.rm=TRUE)
data.BA.SP$BATOT.COMPET <- BATOT.COMPET; rm(BATOT.COMPET)
### create data frame
names(data.BA.SP) <- c("TreeID",names(data.BA.SP)[-1])
data.BA.sp <- merge(data.frame(TreeID=data.bci[["TreeID"]],ecocode=data.bci[["ecocode"]]),data.BA.SP,by="TreeID",sort=FALSE)
## test
if(sum(!data.BA.sp[["TreeID"]] == data.tree[["TreeID"]]) >0) stop("competition index not in the same order than data.tree")
## save everything as a list
list.bci <- list(data.tree=data.tree,data.BA.SP=data.BA.sp,data.traits=data.traits)
save(list.bci,file="./data/process/list.bci.Rdata")