From 61ae1ff59e6d981781aeab806b2ffc62ba1cf02d Mon Sep 17 00:00:00 2001
From: fhui28 <fhui28@gmail.com>
Date: Sun, 18 Aug 2013 08:01:35 +1000
Subject: [PATCH] modified canada R script for Canada_Data2George_20130816
---
....data.CANADA-fhv1.R => merge.data.CANADA.R | 56 +++++++++----------
1 file changed, 25 insertions(+), 31 deletions(-)
rename merge.data.CANADA-fhv1.R => merge.data.CANADA.R (71%)
diff --git a/merge.data.CANADA-fhv1.R b/merge.data.CANADA.R
similarity index 71%
rename from merge.data.CANADA-fhv1.R
rename to merge.data.CANADA.R
index 1ea8d9d..5f151a9 100644
--- a/merge.data.CANADA-fhv1.R
+++ b/merge.data.CANADA.R
@@ -6,7 +6,8 @@ rm(list = ls()); source("./R/format.function.R"); library(reshape)
## READ DATA
####################
### read individuals tree data
-data.canada <- read.csv("./data/raw/DataCanada/Canada_Data2George_20130815.csv",header=TRUE,stringsAsFactors =FALSE)
+#data.canada <- read.csv("./data/raw/DataCanada/Canada_Data2George_20130815.csv",header=TRUE,stringsAsFactors =FALSE)
+data.canada <- read.csv("./data/raw/DataCanada/Canada_Data2George_20130816.csv",header=TRUE,stringsAsFactors =FALSE)
data.canada <- data.canada[which(!is.na(data.canada$Species)),]
colnames(data.canada)[2] <- "Species"
@@ -31,13 +32,14 @@ species.clean <- read.csv("./data/raw/DataCanada/FIA_REF_SPECIES.csv",stringsAsF
#############
## change unit and names of variables to be the same in all data for the tree
-data.canada$G <- (data.canada[["FinalDBH"]]-data.canada[["InitDBH"]])/data.canada$Interval ## diameter growth in mm per year
-data.canada$year <- data.canada$Interval ## number of year between measurement/missing!
-data.canada$D <- data.canada[["InitDBH"]] ## diameter in mm
-data.canada$dead <- rep(NA,length(data.canada[["Species"]])) ## dummy variable for dead tree 0 alive 1 dead/missing!
+data.canada$G <- 10*(data.canada$FinalDBH-data.canada$InitDBH)/data.canada$Interval ## diameter growth in mm per year
+data.canada$G[which(data.canada$InitDBH == 0 | data.canada$FinalDBH == -999)] <- NA
+data.canada$year <- data.canada$Interval ## number of year between measurement - MISSING
+data.canada$D <- data.canada[["InitDBH"]]; data.canada$D[data.canada$D == 0] <- NA ;## diameter in cm
+data.canada$dead <- as.numeric(data.canada$FinalDBH == -999) ## dummy variable for dead tree 0 alive 1 dead
data.canada$sp <- as.character(data.canada[["Species"]]) ## species code
-data.canada$plot <- (data.canada[["PlotID"]]) ## plot code
-data.canada$htot <- rep(NA,length(data.canada[["Species"]]))## height of tree in m / missing
+data.canada$plot <- (data.canada[["PLOT_ID"]]) ## plot code
+data.canada$htot <- rep(NA,length(data.canada[["Species"]]))## height of tree in m - MISSING
data.canada$tree.id <- data.canada$PLOTTREE ## tree unique id
data.canada$sp.name <- NA;
for(i in 1:length(unique(data.canada$sp))) {
@@ -47,16 +49,10 @@ for(i in 1:length(unique(data.canada$sp))) {
############################
## merge greco to have no ecoregion with low number of observation
-# greco <- read.csv(file = "./data/raw/DataSpain/R_Ecoregion.csv", header = T)
-# greco <- greco[,c("Plot_ID_SFI","BIOME","eco_code")]
-# greco2 <- greco[!duplicated(greco$Plot),];
-# rm(greco)
-#
-# data.canada <- merge(data.canada, greco2, by = "Plot_ID_SFI")
-# rm(greco2)
-#
-# table(data.canada$eco_code)
-# ## There's an eco-region with no code, and one with 55 sites
+greco <- read.csv(file = "./data/raw/DataCanada/EcoregionCodes.csv", header = T, sep = "\t")
+
+table(data.canada$Ecocode)
+## There is only four ecoregions though; do you want to aggregate into divsion still?
#
# library(RColorBrewer); mycols <- brewer.pal(10,"Set3");
# ecoreg <- unclass(data.canada$eco_code);
@@ -67,22 +63,21 @@ for(i in 1:length(unique(data.canada$sp))) {
# data.canada$eco_codemerged <- combine_factor(data.canada$eco_code, c(1:8,6,9))
#######################
-# ## variable percent dead/cannot do with since dead variable is missing
-# ###compute numer of dead per plot to remove plot with disturbance
-# perc.dead <- tapply(data.canada[["dead"]],INDEX=data.canada[["idp"]],FUN=function.perc.dead)
+## variable percent dead/cannot do with since dead variable is missing
+## compute numer of dead per plot to remove plot with disturbance
+perc.dead <- tapply(data.canada[["dead"]],INDEX=data.canada[["plot"]],FUN=function.perc.dead)
# ## VARIABLE TO SELECT PLOT WITH NOT BIG DISTURBANCE KEEP OFTHER VARIABLES IF AVAILABLE (disturbance record)
-# data.canada <- merge(data.canada,data.frame(idp=as.numeric(names(perc.dead)),perc.dead=perc.dead),sort=FALSE)
-data.canada$perc.dead <- NA
+data.canada <- merge(data.canada,data.frame(plot=names(perc.dead),perc.dead=perc.dead), by = "plot", sort=FALSE)
###########################################################
### PLOT SELECTION FOR THE ANALYSIS
###################
-# ## Remove data with mortality == 1 or 2
-# table(data.canada$Mortality_Cut)
-# data.canada <- subset(data.canada,subset= (data.canada[["Mortality_Cut"]] == 0 | data.canada[["Mortality_Cut"]] == ""))
+## Remove data with dead == 1
+table(data.canada$dead)
+## Nothing to remove
-colnames(data.canada)[c(2,5,10,12,14)] <- c("sp","plot","w","ecocode","PP")
-vec.abio.var.names <- c("MAT","PP")
+colnames(data.canada)[c(3,1,11,13)] <- c("sp","plot","w","ecocode")
+vec.abio.var.names <- c("MAT","MAP")
vec.basic.var <- c("tree.id","sp","sp.name","plot","ecocode","D","G","dead","year","htot","Lon","Lat","perc.dead")
data.tree <- subset(data.canada,select=c(vec.basic.var,vec.abio.var.names))
@@ -91,17 +86,16 @@ data.tree <- subset(data.canada,select=c(vec.basic.var,vec.abio.var.names))
###########################
data.BA.SP <- BA.SP.FUN(id.tree=as.vector(data.canada[["tree.id"]]), diam=as.vector(data.canada[["D"]]),
sp=as.vector(data.canada[["sp"]]), id.plot=as.vector(data.canada[["plot"]]),
- weights=1/(10000*data.canada[["SubPlotSize"]]), weight.full.plot=NA)
+ weights=1/(10000*data.canada[["SubPlot_Size"]]), 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
+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.canada[["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],MARGIN=1,FUN=sum,na.rm=TRUE)
+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("tree.id",names(data.BA.SP)[-1])
--
GitLab