In VincyaneBadouard/TreeData: Forest Inventories Harmonization & Correction

Packages

library(readr)
library(TreeData)
library(data.table)
library(ggplot2)
library(ggrepel)

Data

PanamaFormated <- read_csv("D:/VSC TmFO/Data/StandardisedData/panama-plots_formated.csv")
setDT(PanamaFormated)

# Mal formaté :
PanamaFormated[, c("IdTree", "IdStem") := list(as.character(IdTree), as.character(IdStem))]
PanamaFormated[, Species := ScientificName]
PanamaFormated[, ScientificName := paste(Genus, Species)]
setnames(PanamaFormated, "SubPlot", "Subplot")
setnames(PanamaFormated, c("Xutm", "Yutm"), c("XTreeUTM", "YTreeUTM"))

Corrections

Bota correction

WFO_Backbone <- file.choose()
load(WFO_Backbone)
Rslt_bota <- BotanicalCorrection(TestData, Source = "WFO", WFOData = WFO_Backbone)

ScfcCor <- unique(Rslt_bota[ScientificNameCor != ScientificName,
                            list(ScientificName, ScientificNameCor,
                                 Family, FamilyCor,
                                 Genus, GenusCor,
                                 Species, SpeciesCor, Subspecies,
                                 BotanicalCorrectionSource, Comment)
])

Life status correction

Rslt_status <- StatusCorrection(PanamaFormated,
                                InvariantColumns = c("Site",
                                                     "Genus",
                                                     "Species",
                                                     "Family",
                                                     "ScientificName"))

pdf("LifeStatusCorrectionPlots.pdf", width = 25, height = 10)
LifeStatusCorrectionPlot(Rslt_status, OnlyCorrected = T, SeveralWindows = F)
dev.off()

Plot - status

Rslt <- Rslt_status
# Order IDs and times in ascending order 
Rslt <- Rslt[order(IdStem, Year)]

IdStemCor <- Rslt[LifeStatus != LifeStatusCor, IdStem] #  corrected stems

RsltCor <- Rslt[IdStem %in% IdStemCor] #  corrected stems

pdf("StatusCorrectionPlots.pdf", width = 25, height = 10)
for(p in 1:(ceiling(length(unique(IdStemCor))/9))){ 
  print(ggplot(Rslt) +
          aes(x = Year, y = LifeStatusCor) +
          geom_point(shape = "circle", size = 3.9, colour = "#46337E") +
          geom_smooth(span = 0.75) +
          theme_minimal() +
          facet_wrap(vars(IdTree), scales = "free") +
          ggforce::facet_wrap_paginate(vars(IdStem), scales = "free", ncol = 3, nrow = 3, page = p) 

  )

}
dev.off()

Taper Correction

Rslt_taper <- TaperCorrection(Data = PanamaFormated)
save(Rslt_taper, file = "PanamaDBHCor_taper.rda")

pdf("TaperCorrectionPlots.pdf", width = 25, height = 10)
DiameterCorrectionPlot(Rslt_taper, OnlyCorrected = TRUE, CorCol = "TaperCorDBH", SeveralWindows = F)
dev.off()

Diameter correction

WhatToCorrect = c("POM change", "punctual", "shift"), CorrectionType = c("quadratic", "linear", "individual", "phylogenetic hierarchical"),

# Rslt_Phylo <- DiameterCorrection(
#   Data = PanamaFormated,
#   WhatToCorrect = c("POM change", "punctual", "shift"),
#   CorrectionType = c("linear", "phylogenetic hierarchical"),
#   MinIndividualNbr = 5)
# 
# save(Rslt_Phylo, file = "PanamaDBHCor_Phylo.rda")
load("D:/VSC TmFO/Data/StandardisedData/PanamaDBHCor_Phylo.rda")


# Rslt_Ind <- DiameterCorrection(
#   PanamaFormated,
#   WhatToCorrect = c("POM change", "punctual", "shift"),
#   CorrectionType = c("linear", "individual"),
#   MinIndividualNbr = 5)
# 
# save(Rslt_Ind, file = "PanamaDBHCor_Ind.rda")
load("D:/VSC TmFO/Data/StandardisedData/PanamaDBHCor_Ind.rda")

Plot - Diameter Correction

# # # Plot phylogenetic hierarchical correction
# pdf("DiameterCorrectionPlots_Phylo.pdf", width = 25, height = 10)
# DiameterCorrectionPlot(Rslt_Phylo, OnlyCorrected = TRUE, SeveralWindows = F)
# dev.off()
# # 
# # # Plot individual correction
# pdf("DiameterCorrectionPlots_Ind.pdf", width = 25, height = 10)
# DiameterCorrectionPlot(Rslt_Ind, OnlyCorrected = TRUE, SeveralWindows = F)
# dev.off()

Machine à debug

# options(warn = 2) # trace warning
# options(warn = 0) # when debug is over


PanamaFormated <- UniqueMeasurement(PanamaFormated, KeepMeas = c("MaxHOM", "MaxDate"), ID = "IdStem")
DataTree <- PanamaFormated[IdStem %in% "443101"]
Data <- PanamaFormated
DefaultHOM = 1.3
MaxDBH = 500
PositiveGrowthThreshold = 5
NegativeGrowthThreshold = -2

Pioneers = NULL
PioneersGrowthThreshold = 7.5

WhatToCorrect = c("POM change", "punctual", "shift")
CorrectionType = c("linear", "phylogenetic hierarchical")

DBHRange = 10
MinIndividualNbr = 5
OtherCrit = NULL
Digits = 1L

coef = 0.9

DetectOnly = FALSE

DiameterCorrectionByTree(
  DataTree, Data,
  WhatToCorrect = c("POM change", "punctual", "shift"),
  CorrectionType = c("linear", "phylogenetic hierarchical")
)

Problèmes rencontrés

# IdStem = "1078509"
DBHCor = c(3.7, 2.8, 2.8, 2.8, NA)
cresc = c(-0.9, 0, 0, NA)
Time = c(2016, 2017, 2018, 2019, 2021)
RegressionInterpolation(Y = cresc, X = Time[-1], CorrectionType = "quadratic") # cresc_Corr
RegressionInterpolation(Y = cresc, X = Time[-1], CorrectionType = "linear") # cresc_Corr
RegressionInterpolation(Y = DBHCor, X = Time, CorrectionType = "quadratic") # cresc_Corr
RegressionInterpolation(Y = DBHCor, X = Time, CorrectionType = "linear") # cresc_Corr


# IdStem = "1078558"
DBHCor = c(1.1, 1.3, 1.5, 1.5, NA)
cresc = c(0.2, 0.2, 0, NA)
Time = c(2016, 2017, 2018, 2019, 2021)
plot(Time[-1], cresc)
RegressionInterpolation(Y = cresc, X = Time[-1], CorrectionType = "quadratic") # cresc_Corr
RegressionInterpolation(Y = cresc, X = Time[-1], CorrectionType = "linear") # cresc_Corr
RegressionInterpolation(Y = DBHCor, X = Time, CorrectionType = "quadratic") # cresc_Corr
RegressionInterpolation(Y = DBHCor, X = Time, CorrectionType = "linear") # cresc_Corr

# IdStem = "443136"
DBHCor = c(17.7, 26.2, NA, 34.6, 34.6, 34.6, 38.0)
cresc = c(0.85, NA, 0.9333333, 0, 0, 1.7)
Time = c(1998, 2008, 2016, 2017, 2018, 2019, 2021)
cresc_abs = ComputeIncrementation(Var = DBHCor, Type = "absolute", Time = Time)

ComputeMeasuredValues(initial_values = DBHCor, cresc_abs = cresc_abs)

cresc_cor <- RegressionInterpolation(Y = cresc, X = Time[-1], CorrectionType = "quadratic") # cresc_Corr
RegressionInterpolation(Y = cresc, X = Time[-1], CorrectionType = "linear") # cresc_Corr
RegressionInterpolation(Y = DBHCor, X = Time, CorrectionType = "quadratic") # cresc_Corr
RegressionInterpolation(Y = DBHCor, X = Time, CorrectionType = "linear") # cresc_Corr

DBHCor[3] <- DBHCor[2] + cresc_cor[2]

# IdStem = "443175"
DBHCor = c(34.5, NA, 34.0, 34.6, 35.0, 34.9, NA)
cresc = c(NA, -0.02777778,  0.6, 0.4, -0.1)
Time = c(1998, 2008, 2016, 2017, 2018, 2019, 2021)
RegressionInterpolation(Y = cresc, X = Time[-1], CorrectionType = "quadratic") # cresc_Corr
RegressionInterpolation(Y = cresc, X = Time[-1], CorrectionType = "linear") # cresc_Corr
RegressionInterpolation(Y = DBHCor, X = Time, CorrectionType = "quadratic") # cresc_Corr
RegressionInterpolation(Y = DBHCor, X = Time, CorrectionType = "linear") # cresc_Corr

# ------------------- Ca marche mieux avec le DBH qu'avec le cresc

# IdStem = "443233"
DBHCor = c(92.9, 77.1, 78.0, 78.7, 78.8, 78.2, 79.2) # ... c'est la 1ere valeur qui semble fausse, mais l'algorithme la considère toujours comme de référence (HOM à 1.3 partout)

# IdStem = "443243" POM change uniquement à la 2e valeur
DBHCor = c(141, 134, 137.5, 138, 139, 139.6, 142.3)
cresc = c(-0.7, 0.4375, 0.5, 1, 0.6, 1.35)
cresc_abs = c(-7.0,  3.5,  0.5,  1.0,  0.6,  2.7)
cresc[-1] = NA
mean(cresc)
DBHCor[2] = NA
Time = c(1998, 2008, 2016, 2017, 2018, 2019, 2021)
plot(Time[-1], cresc)
RegressionInterpolation(Y = cresc, X = Time[-1], CorrectionType = "quadratic") # cresc_Corr là le quadratic est mieux
RegressionInterpolation(Y = cresc, X = Time[-1], CorrectionType = "linear") # cresc_Corr
RegressionInterpolation(Y = DBHCor, X = Time, CorrectionType = "quadratic") # cresc_Corr
RegressionInterpolation(Y = DBHCor, X = Time, CorrectionType = "linear") # cresc_Corr