Paracou.R 11 KB
Newer Older
Daniel Falster's avatar
Daniel Falster committed
1
2
#!/usr/bin/env Rscript

Georges Kunstler's avatar
merged    
Georges Kunstler committed
3
### MERGE paracou DATA
Georges Kunstler's avatar
Georges Kunstler committed
4

Georges Kunstler's avatar
Georges Kunstler committed
5
rm(list = ls())
Georges Kunstler's avatar
Georges Kunstler committed
6
source("R/format.data/format-fun.R")
7
dir.create("output/formatted/Paracou", recursive=TRUE,showWarnings=FALSE)
Georges Kunstler's avatar
Georges Kunstler committed
8
library(reshape,quietly=TRUE)
Georges Kunstler's avatar
Georges Kunstler committed
9
 
10
############################ read individuals tree data
11
data.paracou <- read.table("data/raw/Paracou/20130717_paracou_1984_2012.csv",
12
13
                           header=TRUE,stringsAsFactors=FALSE,sep = ";",
                           na.strings = "NULL")
Georges Kunstler's avatar
Georges Kunstler committed
14
# select good columns
Georges Kunstler's avatar
Georges Kunstler committed
15
data.paracou2 <- data.paracou[,c("foret","parcelle","carre","arbre",
16
17
18
19
20
                                "vernaculaire","idtaxon",
                                "x","y","circ_2001","code_2001",
                                "circ_2005","code_2005",
                                "circ_2009","code_2009","campagne_mort",
                                "type_mort")]
Georges Kunstler's avatar
Georges Kunstler committed
21
colnames(data.paracou2) <- c("forest","cluster","plot","tree","vernacular",
22
23
24
                            "taxonid","x","y",
                            "circum2001","code2001","circum2005","code2005",
                            "circum2009",
25
                            "code2009","yeardied","typedeath")
Georges Kunstler's avatar
Georges Kunstler committed
26

Georges Kunstler's avatar
Georges Kunstler committed
27
28
29
30
31
32
33
34
35
36
37
38
39
40
## FOR parcelle 16 assume circum2009 is 2010 and 2001 1998 and include in the year
data.paracou2$circum2009[data.paracou2$cluster==16] <-
    data.paracou$circ_2010[data.paracou$parcelle==16]
data.paracou2$code2009[data.paracou2$cluster==16] <-
    data.paracou$code_2010[data.paracou$parcelle==16]

data.paracou2$circum2001[data.paracou2$cluster==16] <-
    data.paracou$circ_1998[data.paracou$parcelle==16]
data.paracou2$code2001[data.paracou2$cluster==16] <-
    data.paracou$code_1998[data.paracou$parcelle==16]

data.paracou <-  data.paracou2
rm(data.paracou2)

41
### change numeric separator
Georges Kunstler's avatar
Georges Kunstler committed
42
43
numeric.col.name <-   c("x","y","circum2001","code2001","circum2005",
                        "code2005","circum2009","code2009")
Georges Kunstler's avatar
merged    
Georges Kunstler committed
44
for(k in numeric.col.name){ 
45
46
	data.paracou[,k] <- gsub(",",".",data.paracou[,k])
        data.paracou[,k] <- as.numeric(data.paracou[,k])
Georges Kunstler's avatar
merged    
Georges Kunstler committed
47
    } ## Replace all , in decimals with .
48
49
data.paracou$tree.id <- apply(data.paracou[,c("cluster","plot","tree")],
                              1,paste,collapse="_");
Georges Kunstler's avatar
Georges Kunstler committed
50
data.paracou$sp <- paste("sp",data.paracou[["taxonid"]],sep=".")
Georges Kunstler's avatar
merged    
Georges Kunstler committed
51
data.paracou <- data.paracou[,c(ncol(data.paracou),1:(ncol(data.paracou)-1))]	
Georges Kunstler's avatar
Georges Kunstler committed
52

Georges Kunstler's avatar
Georges Kunstler committed
53
############################ SELECT OBSERVATION WITHOUT PROBLEMS
Georges Kunstler's avatar
merged    
Georges Kunstler committed
54
## REMOVE ALL TREES WITH X OR Y >250 m 
Georges Kunstler's avatar
Georges Kunstler committed
55
56
57
58
data.paracou <- subset(data.paracou,subset=(!is.na(data.paracou[["y"]])) &
                        (!is.na(data.paracou[["x"]])) &
                        data.paracou[["x"]]<251 &
                        data.paracou[["y"]]<251)
Georges Kunstler's avatar
merged    
Georges Kunstler committed
59
#### REMOVE PLOTs 16 17 18 ACCORDING TO  GHSILAIN
60
data.paracou <- subset(data.paracou,
Georges Kunstler's avatar
Georges Kunstler committed
61
62
63
64
65
66
67
                       subset=! data.paracou[["cluster"]] %in% 16:18)
## ## keep only tree alive in 2001 and 1998 for cluster 16
## VEC.select <- data.paracou$yeardied
## VEC.select <- NA
## VEC.select[data.paracou$cluster != 16] <-
##     !(data.paracou$yeardied[data.paracou$cluster != 16] <= 2001 &
##       !is.na(data.paracou$yeardied[data.paracou$cluster != 16]))
Georges Kunstler's avatar
Georges Kunstler committed
68

Georges Kunstler's avatar
Georges Kunstler committed
69
70
71
72
73
74
## VEC.select[data.paracou$cluster == 16] <-
##     !(data.paracou$yeardied[data.paracou$cluster == 16] <= 1998 &
##       !is.na(data.paracou$yeardied[data.paracou$cluster == 16]))

## data.paracou <- subset(data.paracou,
##                        subset=VEC.select)
Georges Kunstler's avatar
Georges Kunstler committed
75
76
77

## plot each plot
## pdf("figs/plots.paracou.pdf")
78
79
80
81
## lapply(unique(data.paracou[["cluster"]]),FUN=fun.circles.plot,
## data.paracou[['x']],data.paracou[['y']],
##        data.paracou[["cluster"]],data.paracou[["circum2009"]],
## inches=0.2,fg.l=data.paracou$plot)
Georges Kunstler's avatar
Georges Kunstler committed
82
83
84
85
86
87
88
89
90
91
92
93
94
95
## dev.off()
## #

## create small plot
square.s.t <- 20
make.quad <- do.call("rbind",lapply(unique(data.paracou$cluster),
                           FUN=fun.quadrat.cluster,
                           cluster=data.paracou$cluster,
                           tree.id=data.paracou$tree.id,
                           x=data.paracou$x,
                           y=data.paracou$y,
                           square.s=square.s.t))
data.paracou <- merge(data.paracou, make.quad,by='tree.id')

Georges Kunstler's avatar
Georges Kunstler committed
96
97
Levels.cluster <- levels(as.factor(data.paracou$cluster))
n.cluster <- length(Levels.cluster) # 10 clusters => 10 "plots"
Georges Kunstler's avatar
Georges Kunstler committed
98

Georges Kunstler's avatar
Georges Kunstler committed
99
100
101
102
103
104
for (i in Levels.cluster){
    dat <- subset(data.paracou, subset=data.paracou[["cluster"]]==i)
    x11()
    symbols(dat$x,dat$y,circles=dat$circum2001,inches=0.2,main=i,
            fg=unclass(factor(dat$make.quad)))
}   
Georges Kunstler's avatar
Georges Kunstler committed
105
106
107



108
######################################## MASSAGE TREE DATA
109
### read species names
Georges Kunstler's avatar
Georges Kunstler committed
110

111
112
113
114
species.clean <- read.csv("data/raw/Paracou/20130717_paracou_taxonomie.csv",
                          stringsAsFactors=FALSE, header = T, sep = ";")


115
species.clean$sp <- paste("sp.",species.clean[["idTaxon"]],sep="")
116
117
species.clean$Latin_name <-  paste(species.clean[["Genre"]],
                                   species.clean[["Espece"]],sep=" ")
118
## keep only one row pers idTaxon
119
120
121
species.clean <- subset(species.clean,subset=!duplicated(species.clean[["sp"]]),
                        select=c("sp","Latin_name","Genre","Espece",
                            "Famille","idCIRAD"))
122
## select only species present in data base
123
124
species.clean <-  subset(species.clean,subset=species.clean[["sp"]] %in%
                         data.paracou[["sp"]])
125
## percentage of species with no taxonomic identification 
126
127
length(grep("Indet",species.clean[["Latin_name"]]))/nrow(species.clean)
## 25% agree with Bruno ?
Georges Kunstler's avatar
Georges Kunstler committed
128

Georges Kunstler's avatar
Georges Kunstler committed
129
## GET PALM AND TREE FERN TO REMOVE THEIR GROWTH
130
131
132
133
134
sp.palm.fern <- species.clean$sp[species.clean$Famille %in%
                                 c('Arecaceae','Cyatheaceae','Dicksoniaceae',
                                   'Metaxyaceae', 'Cibotiaceae',
                                   'Loxomataceae', 'Culcitaceae',
                                   'Plagiogyriaceae', 'Thyrsopteridaceae')]
Georges Kunstler's avatar
Georges Kunstler committed
135
136


137
############################################ FORMAT INDIVIDUAL TREE DATA
138
139
140
141
142
143
data.paracou2 <- data.paracou[rep(1:nrow(data.paracou),each=2),
                              c(1:10,(ncol(data.paracou)-2):ncol(data.paracou))]
rownames(data.paracou2) <- 1:nrow(data.paracou2)
data.paracou2 <- as.data.frame(data.paracou2)
data.paracou2$yr1 <- rep(c(2001,2001+4),nrow(data.paracou))
data.paracou2$yr2 <- rep(c(2005,2005+4),nrow(data.paracou))
Georges Kunstler's avatar
Georges Kunstler committed
144
145
146
147
## data.paracou2$yr1[data.paracou2$cluster==16] <- rep(c(1998,1998+7),
##                       nrow(data.paracou[data.paracou2$cluster==16,]))
## data.paracou2$yr2[data.paracou2$cluster==16] <- rep(c(2005,2005+5),
##                       nrow(data.paracou[data.paracou2$cluster==16,])) # For 16
Georges Kunstler's avatar
Georges Kunstler committed
148

Georges Kunstler's avatar
merged    
Georges Kunstler committed
149
data.paracou2$year <- rep(c(4,4),nrow(data.paracou))
Georges Kunstler's avatar
Georges Kunstler committed
150
151
152
data.paracou2$year[data.paracou2$cluster==16] <- rep(c(7,5),
                      nrow(data.paracou[data.paracou2$cluster==16,])) # For 16

153
154
155
156
157
158
159
160
data.paracou2$dbh1 <- c(rbind(data.paracou$circum2001/pi,
                              data.paracou$circum2005/pi))
data.paracou2$dbh2 <- c(rbind(data.paracou$circum2005/pi,
                              data.paracou$circum2009/pi))
data.paracou2$code1 <- c(as.numeric(rbind(data.paracou$code2001,
                                          data.paracou$code2005)))
data.paracou2$code2 <- c(as.numeric(rbind(data.paracou$code2005,
                                          data.paracou$code2009)))
Georges Kunstler's avatar
merged    
Georges Kunstler committed
161
data.paracou2$dead <- rep(0,nrow(data.paracou)*2)
Georges Kunstler's avatar
Georges Kunstler committed
162
163
164
165
data.paracou2$dead[c(as.numeric(data.paracou[["yeardied"]]) %in% 2002:2005 &
                     (!is.na(data.paracou[["yeardied"]])),
                     as.numeric(data.paracou[["yeardied"]]) %in% 2006:2009 &
                     (!is.na(data.paracou[["yeardied"]])))] <- 1
166
## remove tree dead at first census for both date (census 2001-2005 2005-2009)
167
168
169
170
data.paracou <- subset(data.paracou2,
                       subset=!(data.paracou2[['yr1']] ==2005 &
                                (as.numeric(data.paracou[["yeardied"]]) %in%
                                 2002:2005 & (!is.na(data.paracou[["yeardied"]])))))
Georges Kunstler's avatar
merged    
Georges Kunstler committed
171
## change unit and names of variables to be the same in all data for the tree 
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
data.paracou$G <- 10*(data.paracou$dbh2-data.paracou$dbh1)/data.paracou$year
## diameter growth in mm per year
data.paracou$G[data.paracou$code1>0] <- NA
## indivs with code indicating problem in dbh measurment at dbh1
data.paracou$G[data.paracou$code2>0] <- NA
## indivs with code indicating problem in dbh measurment at dbh2
data.paracou$BA.G <- (pi*(data.paracou$dbh2/2)^2-pi*(data.paracou$dbh1/2)^2)/
    data.paracou$year ## BA growth in cm2 per year
data.paracou$BA.G[data.paracou$code1>0] <- NA
## indivs with code indicating problem in dbh measurment at dbh1
data.paracou$BA.G[data.paracou$code2>0] <- NA
## indivs with code indicating problem in dbh measurment at dbh2
data.paracou$G[data.paracou$sp %in% sp.palm.fern] <- NA
## remove palm and fern growth
data.paracou$BA.G[data.paracou$sp %in% sp.palm.fern] <- NA
## remove palm and fern growth
Georges Kunstler's avatar
Georges Kunstler committed
188
data.paracou$census <- data.paracou$yr1
189
190
191
data.paracou$D <- data.paracou[["dbh1"]]
data.paracou$D[data.paracou$D == 0] <- NA ;## diameter in cm
data.paracou$cluster <- paste("p",data.paracou$cluster,sep=".")
192
data.paracou$htot <- rep(NA,length(data.paracou[["G"]])) ## height of tree in m
193
data.paracou$obs.id <- 1:nrow(data.paracou)
Georges Kunstler's avatar
Georges Kunstler committed
194
data.paracou$plot <-  data.paracou$make.quad
Georges Kunstler's avatar
Georges Kunstler committed
195
### delete recruit in 2001 or 2005 for first census
Georges Kunstler's avatar
merged    
Georges Kunstler committed
196
data.paracou <- subset(data.paracou,subset=!is.na(data.paracou$D))
Georges Kunstler's avatar
Georges Kunstler committed
197
## minimum circumfer 30 delete all tree with a dbh <30/pi,
Georges Kunstler's avatar
merged    
Georges Kunstler committed
198
data.paracou <-  subset(data.paracou,subset= data.paracou[["D"]]>(30/pi))
Georges Kunstler's avatar
Georges Kunstler committed
199
## add latin name
200
201
data.paracou <- merge(data.paracou,subset(species.clean,
                                          select=c("sp","Latin_name")),by="sp")
Georges Kunstler's avatar
Georges Kunstler committed
202
203
data.paracou$sp.name <- data.paracou$Latin_name

Georges Kunstler's avatar
Georges Kunstler committed
204
205
206
207
208
209
210
211
212
213
214
215
216
####################
## add lat long
data.paracou[["Lon"]] <- -52.900002
data.paracou[["Lat"]] <- 5.38 
data.paracou[["ecocode"]] <- "tropical"

#### get wc climate 
source("R/utils/climate.R")
clim <- GetClimate(data.paracou$Lat,data.paracou$Lon)
data.paracou$MAT <- clim$MAT
data.paracou$MAP <- clim$MAP

###################### PLOT SELECTION FOR THE ANALYSIS - NEEDS REDOING
217
218
vec.basic.var <- c("obs.id","tree.id", "sp", "sp.name", "cluster","plot",
                   "ecocode", "D", "G","BA.G","year", "dead",
Georges Kunstler's avatar
Georges Kunstler committed
219
    'Lon','Lat',"x", "y", "census",'MAT','MAP')
Georges Kunstler's avatar
Georges Kunstler committed
220
data.tree <- subset(data.paracou, select = c(vec.basic.var))
221
222


223
224
## select tree above 10 cm and last census only
data.tree <- subset(data.tree,subset=data.tree$D>10 & !is.na(data.tree$D))
225
## data.tree <- subset(data.tree,subset=data.tree$census==2005)
226

227
data.tree <- subset(data.tree,subset=!is.na(data.tree$x) & !is.na(data.tree$y))
Georges Kunstler's avatar
Georges Kunstler committed
228
229
## convert var factor in character or numeric
data.tree <- fun.convert.type.B(data.tree)
230

231
write.csv(data.tree,file="output/formatted/Paracou/tree.csv",row.names = FALSE)
Georges Kunstler's avatar
Georges Kunstler committed
232

233
234
235
### write data plot with variables only at the plot level.
### I should delete them from the tree table but so far
### I keep them to not distroy later code
Georges Kunstler's avatar
Georges Kunstler committed
236

Georges Kunstler's avatar
Georges Kunstler committed
237
vec.basic.var.p <- c("plot", "cluster", "Lon","Lat","ecocode",'MAT','MAP')
238
239
data.plot <- subset(data.tree, subset=!duplicated(data.tree$cluster),
                    select = c(vec.basic.var.p))
Georges Kunstler's avatar
Georges Kunstler committed
240
write.csv(data.plot,file="output/formatted/Paracou/plot.csv",row.names = FALSE)
Georges Kunstler's avatar
Georges Kunstler committed
241