R/ci_ht_regression_ch.R

#' @title Modelisation de Circonference en fonction de la Hauteur
#'
#' @description Cette fonction realise la modele Circonference - Hauteur a partir des tiges echantillons, en priorite selon le modele Coop de donnees, sinon a partir d'une regression lineaire.
#'
#' @param Circ : vecteur des circonferences (cm)
#' @param Haut : vecteur des hauteurs totales (cm). Circ et Haut representent les memes individus, classes dans le meme ordre.
#' @return liste de 7 elements : 
#' - C130.min : la circonference minimale des tiges utilisees pour la modelisation 
#' - C130.max : la circonference maximale des tiges utilisees pour la modelisation,
#' - effectif : l'effectif de tiges utilise pour la modelisation, 
#' - modL : le modele retenu pour la modelisation (coop ou lineaire),
#' - par1 : la valeur du parametre 1 du modele,
#' - par2 : la valeur du parametre 2 du modele,
#' - Observ. : les observations concernant le choix du modele.
#'
#' @author Quentin Girard
#' @references Protocole Coop chene... document d Ingrid et Claudine
#'
#' @seealso  ci_ht_fonction_ch
#' @examples
#' # la regression Coop fonctionne
#' Circ.ex1 = c(59,62,35,47,33,48,33,30,37,37,22,37,42,41,43,36,42,26,27,35,30,30,60,49,33,52,49,52,44)
#' Haut.ex1 = c(1570,1620,1400,1460,1410,1520,1270,1270,1370,1360,1230,1390,1470,1410,1520,1500,1510,1300,1260,1410,1300,NA,1590,1550,1430,1610,1500,1480,1570)
#' regr1 <- ci_ht_regression_ch(Circ.ex1, Haut.ex1)
#' regr1
#' plot(Circ.ex1, Haut.ex1) ; lines(sort(Circ.ex1), ci_ht_fonction_ch(sort(Circ.ex1), regr1$par1, regr1$par2, modL = regr1$modL))
#' # la regression Coop ne fonctionne pas
#' Circ.ex2 = c(1:20)
#' Haut.ex2 = 130 + c(1:20)
#' regr2 <- ci_ht_regression_ch(Circ.ex2, Haut.ex2)
#' regr2
#' plot(Circ.ex2, Haut.ex2) ; lines(sort(Circ.ex2), ci_ht_fonction_ch(sort(Circ.ex2), regr2$par1, regr2$par2, modL = regr2$modL))
#'
#' @keywords function
#'
#' @include coopR-package.R
#' @family coopR
#' @export
ci_ht_regression_ch <- function(Circ, Haut) {
# verifications des donnees
  if (length(Circ) != length(Haut)) {
    stop("Circ et Haut doivent etre de meme longueur")
  }
  data <- data.frame(Circ = Circ, Haut = Haut)
  data <- na.omit(data)
  retour <- list(C130.min = min(data$Circ), 
                 C130.max = max(data$Circ),
                 effectif = nrow(data), 
                 modL     = NA,
                 par1     = NA,
                 par2     = NA,
                 Observ.  = NA)
# en dessous de 10 arbres, on ne fait pas de modele
  if (retour$effectif <= 9){        
    retour[["Observ."]] <- "Modele non calcule : Il faut 10 couples circonference-hauteur au minimum"
  } else {
    regr <- nls(Haut ~ ci_ht_fonction_ch(Circ, a, b, modL = "coop"), data, start=list(a = 1000, b= 0.01), control = nls.control(warnOnly = T))
# la regression 'Coop' fonctionne : on la garde
    if (regr$convInfo$isConv) { #           
      retour[["par1"]] <- summary(regr)$coefficients["a", "Estimate"]
      retour[["par2"]] <- summary(regr)$coefficients["b", "Estimate"]
      retour[["modL"]] <- "coop"
# la regression 'Coop' ne fonctionne pas, on passe en lineaire
    } else {                                                                                          
      regr2 <- lm(Haut - 130 ~ 0 + Circ, data)          
      retour[["par1"]] <- summary.lm(regr2)$coefficients["Circ", "Estimate"]
      retour[["modL"]] <- "lineaire"
      retour[["Observ."]] <- "Regression lineaire : la regression 'Coop' ne convergeait pas"
    }
  } 
  if (!is.na(retour$Observ.)) {
    warning(retour$Observ.)
  }
  return(retour)
}
        
jprenaud-02/coopR documentation built on May 3, 2019, 7:06 p.m.