Analyse BDD copro groupe

Fichier BDD :

data <- read.csv2(
  "coprogroupe(2).csv",
  header = TRUE,
  stringsAsFactors = FALSE,
  na.strings = c("", "NA"))

dim(data)

[1] 168 185

names(data)

  [1] "Mois"       "SGIMultiSP" "SGIPrimiSP" "SGILOTSP"   "SGINulliSP"
  [6] "SGLMultiSP" "SGLPrimiSP" "SGLLOTSP"   "SGLNulliSP" "TRIMultiSP"
 [11] "TRIPrimiSP" "TRILOTSP"   "TRINulliSP" "SKJMultiSP" "SKJPrimiSP"
 [16] "SKJLOTSP"   "SKJNulliSP" "MONMultiSP" "MONPrimiSP" "MONLOTSP"  
 [21] "MONNulliSP" "COCMultiSP" "COCPrimiSP" "COCLOTSP"   "COCNulliSP"
 [26] "GDDMultiSP" "GDDPrimiSP" "GDDLOTSP"   "GDDNulliSP" "PTDMultiSP"
 [31] "PTDPrimiSP" "PTDLOTSP"   "PTDNulliSP" "PARMultiSP" "PARPrimiSP"
 [36] "PARLOTSP"   "PARNulliSP" "MUEMultiSP" "MUEPrimiSP" "MUELOTSP"  
 [41] "MUENulliSP" "EFFMultiSP" "EFFPrimiSP" "EFFTotSP"   "CharSAUSP" 
 [46] "CharSFPSP"  "SGIMultiDP" "SGIPrimiDP" "SGILOTDP"   "SGINulliDP"
 [51] "SGLMultiDP" "SGLPrimiDP" "SGLLOTDP"   "SGLNulliDP" "TRIMultiDP"
 [56] "TRIPrimiDP" "TRILOTDP"   "TRINulliDP" "SKJMultiDP" "SKJPrimiDP"
 [61] "SKJLOTDP"   "SKJNulliDP" "MONMultiDP" "MONPrimiDP" "MONLOTDP"  
 [66] "MONNulliDP" "COCMultiDP" "COCPrimiDP" "COCLOTDP"   "COCNulliDP"
 [71] "GDDMultiDP" "GDDPrimiDP" "GDDLOTDP"   "GDDNulliDP" "PTDMultiDP"
 [76] "PTDPrimiDP" "PTDLOTDP"   "PTDNulliDP" "PARMultiDP" "PARPrimiDP"
 [81] "PARLOTDP"   "PARNulliDP" "MUEMultiDP" "MUEPrimiDP" "MUELOTDP"  
 [86] "MUENulliDP" "EFFMultiDP" "EFFPrimiDP" "EFFTotDP"   "CharSAUDP" 
 [91] "CharSFPDP"  "SGITOTDB"   "SGIMultiDB" "SGIPrimiDB" "SGILOTDB"  
 [96] "SGINulliDB" "SGLTOTDB"   "SGLMultiDB" "SGLPrimiDB" "SGLLOTDB"  
[101] "SGLNulliDB" "TRITOTDB"   "TRIMultiDB" "TRIPrimiDB" "TRILOTDB"  
[106] "TRINulliDB" "SKJTOTDB"   "SKJMultiDB" "SKJPrimiDB" "SKJLOTDB"  
[111] "SKJNulliDB" "MONTOTDB"   "MONMultiDB" "MONPrimiDB" "MONLOTDB"  
[116] "MONNulliDB" "COCTOTDB"   "COCMultiDB" "COCPrimiDB" "COCLOTDB"  
[121] "COCNulliDB" "GDDTOTDB"   "GDDMultiDB" "GDDPrimiDB" "GDDLOTDB"  
[126] "GDDNulliDB" "PTDTOTDB"   "PTDMultiDB" "PTDPrimiDB" "PTDLOTDB"  
[131] "PTDNulliDB" "PARTOTDB"   "PARMultiDB" "PARPrimiDB" "PARLOTDB"  
[136] "PARNulliDB" "MUETOTDB"   "MUEMultiDB" "MUEPrimiDB" "MUELOTDB"  
[141] "MUENulliDB" "EFFMultiDB" "EFFPrimiDB" "EFFTOTDB"   "CharSAUDB" 
[146] "CharSFPDB"  "OesoT"      "TeladoT"    "HaemoT"     "OesoDP"    
[151] "TeladoDP"   "HaemoDP"    "OesoSP"     "TeladoSP"   "HaemoSP"   
[156] "NECSternDB" "NECSternDP" "NECSternSP" "NECLombDB"  "NECLombDP" 
[161] "NECLombSP"  "PoidDB"     "PoidDP"     "PoidSP"     "EffTrSP"   
[166] "ProdCTRLSP" "ProdSP"     "ProdMSP"    "EffTrDP"    "ProdCTRLDP"
[171] "ProdDP"     "ProdMDP"    "EffTrDB"    "ProdCTRLDB" "ProdDB"    
[176] "ProdMDB"    "TempsDP"    "NombreJDP"  "TempsSP"    "NombreJSP" 
[181] "RR_CM"      "TM_MM"      "TN_MM"      "TX_MM"      "RGC_CM"    

# Extraction robuste de l'année à partir de Mois
annee <- 2000 + as.numeric(sub("^[^0-9]+", "", data$Mois))

# Extraction du mois
mois_txt <- tolower(sub("[0-9]+$", "", data$Mois))

# Correction des accents
mois_txt[mois_txt == "février"] <- "fevrier"
mois_txt[mois_txt == "août"] <- "aout"
mois_txt[mois_txt == "décembre"] <- "decembre"

mois_num <- match(
  mois_txt,
  c("janvier", "fevrier", "mars", "avril", "mai", "juin",
    "juillet", "aout", "septembre", "octobre", "novembre", "decembre"))

table(annee, useNA = "ifany")

annee
2013 2014 2015 2016 2017 2018 2019 2020 2021 2022 2023 2024 2025 2026 
  12   12   12   12   12   12   12   12   12   12   12   12   12   12 

table(mois_num, useNA = "ifany")

mois_num
   1    2    3    4    5    6    7    8    9   10   11   12 <NA> 
  14   14   14   14   14   14   14   14   14   14   14   13    1 

Code couleur pour tous les graphiques :

couleur_SP <- "red" 
couleur_DP <- "gold" 
couleur_DB <- "blue"

Evolution annuelle du parasitisme (SP et DP) :

Vérification

length(data$Mois)

[1] 168

length(data$SGILOTSP)

[1] 168

length(data$SGILOTDP)

[1] 168

length(data$SGILOTDB)

[1] 168

Créer l’année à partir du mois

annee_graph <- 2000 + as.numeric(
  sub("^[^0-9]+", "", data$Mois))

Retirer les données manquantes

# SP
ok_SP <- !is.na(annee_graph) &
  !is.na(data$SGILOTSP)

annee_SP <- annee_graph[ok_SP]
sgi_SP <- data$SGILOTSP[ok_SP]

# DP
ok_DP <- !is.na(annee_graph) &
  !is.na(data$SGILOTDP)

annee_DP <- annee_graph[ok_DP]
sgi_DP <- data$SGILOTDP[ok_DP]

# DB
ok_DB <- !is.na(annee_graph) &
  !is.na(data$SGILOTDB)

annee_DB <- annee_graph[ok_DB]
sgi_DB <- data$SGILOTDP[ok_DB]

Moyennes annuelles

moy_SP <- tapply(
  sgi_SP,
  annee_SP,
  mean,
  na.rm = TRUE)

moy_DP <- tapply(
  sgi_DP,
  annee_DP,
  mean,
  na.rm = TRUE)

moy_DB <- tapply(
  sgi_DB,
  annee_DB,
  mean,
  na.rm = TRUE)

Années disponibles

annees <- sort(
  union(
    as.numeric(names(moy_SP)),
    as.numeric(names(moy_DP))))

Enfin le graphique

plot(
  annees,
  moy_SP[as.character(annees)],
  type = "o",
  pch = 16,
  lwd = 2,
  col = "red",
  xlab = "Année",
  ylab = "Nombre moyen d'OPG",
  main = "Evolution de l'infestation par SGI en fonction du temps 
  selon les lots (SP et DP)")

lines(
  annees,
  moy_DP[as.character(annees)],
  type = "o",
  pch = 17,
  lwd = 2,
  col = "gold")

legend(
  "topright",
  legend = c("SP", "DP"),
  col = c("red", "gold"),
  pch = c(16, 17),
  lwd = 2,
  bty = "n")

Diagramme de Tuckey (comparaison SP et DP) :

boxplot(
  data$SGILOTSP,
  data$SGILOTDP,
  names = c("SP", "DP"),
  col = c(couleur_SP, couleur_DP),
  border = c("red", "gold"),
  ylab = "Nombre d'OPG",
  main = "Comparaison du SGI entre SP et DP"
)

Diagramme de Tuckey pour les 3 lots :

boxplot(
  data$SGILOTSP,
  data$SGILOTDP,
  data$SGITOTDB,
  names = c("SP", "DP", "DB"),
  col = c(couleur_SP, couleur_DP, couleur_DB),
  border = c("red", "gold", "blue"),
  ylab = "Nombre d'OPG",
  main = "Comparaison du SGI entre les trois lots"
)

Comparaison SP et DP selon les types de parasites :

Oesophagostomum :

boxplot(
  data$OesoSP,
  data$OesoDP,
  names = c("SP", "DP"),
  col = c("red", "gold"),
  border = c("darkred", "darkgoldenrod"),
  ylab = "Nombre d'OPG d'Oesophagostomum",
  main = "Oesophagostomum : SP vs DP"
)

Teladorsagia :

boxplot(
  data$TeladoSP,
  data$TeladoDP,
  names = c("SP", "DP"),
  col = c("red", "gold"),
  border = c("darkred", "darkgoldenrod"),
  ylab = "Nombre d'OPG de Teladorsagia",
  main = "Teladorsagia : SP vs DP"
)

Haemonchus :

boxplot(
  data$HaemoSP,
  data$HaemoDP,
  names = c("SP", "DP"),
  col = c("red", "gold"),
  border = c("darkred", "darkgoldenrod"),
  ylab = "Nombre d'OPG d'Haemonchus",
  main = "Haemonchus : SP vs DP"
)

Infestation selon la parité et la saison de mise bas :

boxplot(data$SGIMultiSP,
        data$SGIMultiDP,
        data$SGIPrimiSP,
        data$SGIPrimiDP,
        data$SGINulliSP,
        data$SGINulliDP,
        
        names = c("Multi SP", "Multi DP",
                  "Primi SP", "Primi DP",
                  "Nulli SP", "Nulli DP"),
        
        col = c("red", "gold",
                "red", "gold",
                "red", "gold"),
        
        border = c("darkred", "darkgoldenrod",
                   "darkred", "darkgoldenrod",
                   "darkred", "darkgoldenrod"),
        
        las = 2,
        ylab = "Nombre d'OPG"  ,
        main = "Infestation selon la parité et la saison de mise-bas")

Infestation et pâturage :

En fonction du temps :

# SP
ok_SP <- !is.na(data$TempsSP) &
  !is.na(data$SGILOTSP)

plot(data$TempsSP[ok_SP],
     data$SGILOTSP[ok_SP],
     pch = 16,
     col = "red",
     xlab = "Temps au pâturage",
     ylab = "Nombre d'OPG",
     main = "Relation entre temps au pâturage et l'infestation du lot SP")

abline(lm(data$SGILOTSP[ok_SP] ~ data$TempsSP[ok_SP]),
       lwd = 2)

# DP

ok_DP <- !is.na(data$TempsDP) &
  !is.na(data$SGILOTDP)

plot(data$TempsDP[ok_DP],
     data$SGILOTDP[ok_DP],
     pch = 17,
     col = "gold",
     xlab = "Temps au pâturage",
     ylab = "Nombre d'OPG",
     main = "Relation entre temps au pâturage et l'infestation du lot DP")

abline(lm(data$SGILOTDP[ok_DP] ~ data$TempsDP[ok_DP]),
       lwd = 2)

En fonction du nombre de jours au pâturage :

ok_SP <- !is.na(data$NombreJSP) &
  !is.na(data$SGILOTSP)

ok_DP <- !is.na(data$NombreJDP) &
  !is.na(data$SGILOTDP)

plot(data$NombreJSP[ok_SP],
     data$SGILOTSP[ok_SP],
     pch = 16,
     col = "red",
     xlab = "Nombre de jours au pâturage",
     ylab = "Nombre d'OPG",
     main = "Infestation en fonction du nombre de jours au pâturage")

points(data$NombreJDP[ok_DP],
       data$SGILOTDP[ok_DP],
       pch = 17,
       col = "gold")

abline(lm(data$SGILOTSP[ok_SP] ~ data$NombreJSP[ok_SP]),
       lwd = 2,
       col = "red")

abline(lm(data$SGILOTDP[ok_DP] ~ data$NombreJDP[ok_DP]),
       lwd = 2,
       col = "gold")

legend("topleft",
       legend = c("SP", "DP"),
       col = c("red", "gold"),
       pch = c(16, 17),
       bty = "n")

Infestation et température :

ok_SP <- !is.na(data$TM_MM) &
  !is.na(data$SGILOTSP)

ok_DP <- !is.na(data$TM_MM) &
  !is.na(data$SGILOTDP)

plot(data$TM_MM[ok_SP],
     data$SGILOTSP[ok_SP],
     pch = 16,
     col = "red",
     xlab = "Température moyenne (en°C)",
     ylab = "Nombre d'OPG",
     main = "Impact de la température sur l'infestation parasitaire ")

points(data$TM_MM[ok_DP],
       data$SGILOTDP[ok_DP],
       pch = 17,
       col = "gold")

abline(lm(data$SGILOTSP[ok_SP] ~ data$TM_MM[ok_SP]),
       lwd = 2,
       col = "red")

abline(lm(data$SGILOTDP[ok_DP] ~ data$TM_MM[ok_DP]),
       lwd = 2,
       col = "gold")

legend("topleft",
       legend = c("SP", "DP"),
       col = c("red", "gold"),
       pch = c(16, 17),
       bty = "n")

Infestation et précipitation :

ok_SP <- !is.na(data$RR_CM) &
  !is.na(data$SGILOTSP)

ok_DP <- !is.na(data$RR_CM) &
  !is.na(data$SGILOTDP)

plot(data$RR_CM[ok_SP],
     data$SGILOTSP[ok_SP],
     pch = 16,
     col = "red",
     xlab = "Précipitations (en mm)",
     ylab = "Nombre d'OPG",
     main = "Influence des précipitations surl'infestation parasitaire")

points(data$RR_CM[ok_DP],
       data$SGILOTDP[ok_DP],
       pch = 17,
       col = "gold")

abline(lm(data$SGILOTSP[ok_SP] ~ data$RR_CM[ok_SP]),
       lwd = 2,
       col = "red")

abline(lm(data$SGILOTDP[ok_DP] ~ data$RR_CM[ok_DP]),
       lwd = 2,
       col = "gold")

legend("topleft",
       legend = c("SP", "DP"),
       col = c("red", "gold"),
       pch = c(16, 17),
       bty = "n")

Etat corporel et infestation :

NEC lombaire

ok_SP <- !is.na(data$NECLombSP) &
  !is.na(data$SGILOTSP)

ok_DP <- !is.na(data$NECLombDP) &
  !is.na(data$SGILOTDP)

plot(data$NECLombSP[ok_SP],
     data$SGILOTSP[ok_SP],
     pch = 16,
     col = "red",
     xlab = "Note d'état coporelle lombaire",
     ylab = "Nombre d'OPG",
     main = "Nombre d'OPG en fonction de la NEC lombaire")

points(data$NECLombDP[ok_DP],
       data$SGILOTDP[ok_DP],
       pch = 17,
       col = "gold")

abline(lm(data$SGILOTSP[ok_SP] ~ data$NECLombSP[ok_SP]),
       lwd = 2,
       col = "red")

abline(lm(data$SGILOTDP[ok_DP] ~ data$NECLombDP[ok_DP]),
       lwd = 2,
       col = "gold")

legend("topleft",
       legend = c("SP", "DP"),
       col = c("red", "gold"),
       pch = c(16, 17),
       bty = "n")

NEC Sternale

ok_SP <- !is.na(data$NECSternSP) &
  !is.na(data$SGILOTSP)

ok_DP <- !is.na(data$NECSternDP) &
  !is.na(data$SGILOTDP)

plot(data$NECSternSP[ok_SP],
     data$SGILOTSP[ok_SP],
     pch = 16,
     col = "red",
     xlab = "Note d'état coporelle sternale",
     ylab = "Nombre d'OPG",
     main = "Nombre d'OPG en fonction de la NEC sternale")

points(data$NECSternDP[ok_DP],
       data$SGILOTDP[ok_DP],
       pch = 17,
       col = "gold")

abline(lm(data$SGILOTSP[ok_SP] ~ data$NECSternSP[ok_SP]),
       lwd = 2,
       col = "red")

abline(lm(data$SGILOTDP[ok_DP] ~ data$NECSternDP[ok_DP]),
       lwd = 2,
       col = "gold")

legend("topleft",
       legend = c("SP", "DP"),
       col = c("red", "gold"),
       pch = c(16, 17),
       bty = "n")

Poids

ok_SP <- !is.na(data$PoidSP) &
  !is.na(data$SGILOTSP)

ok_DP <- !is.na(data$PoidDP) &
  !is.na(data$SGILOTDP)

plot(data$PoidSP[ok_SP],
     data$SGILOTSP[ok_SP],
     pch = 16,
     col = "red",
     xlab = "Poids",
     ylab = "Nombre d'OPG",
     main = "Nombre d'OPG en fonction du poids")

points(data$PoidDP[ok_DP],
       data$SGILOTDP[ok_DP],
       pch = 17,
       col = "gold")

abline(lm(data$SGILOTSP[ok_SP] ~ data$PoidSP[ok_SP]),
       lwd = 2,
       col = "red")

abline(lm(data$SGILOTDP[ok_DP] ~ data$PoidDP[ok_DP]),
       lwd = 2,
       col = "gold")

legend("topleft",
       legend = c("SP", "DP"),
       col = c("red", "gold"),
       pch = c(16, 17),
       bty = "n")

Différences d’infestation : SP-DP :

ok <- !is.na(mois_num) &
  !is.na(annee) &
  !is.na(data$SGILOTSP) &
  !is.na(data$SGILOTDP)

date_diff <- as.Date(
  paste(annee[ok],
        sprintf("%02d", mois_num[ok]),
        "01",
        sep = "-"))

difference <- data$SGILOTSP[ok] - data$SGILOTDP[ok]

couleurs <- ifelse(difference > 0, "red", "gold")

plot(date_diff,
     difference,
     type = "h",
     lwd = 3,
     col = couleurs,
     xlab = "Date",
     ylab = "OPG SP - OPG DP",
     main = "Différence d'infestation : SP - DP")

abline(h = 0, lty = 2, lwd = 2)

legend("topright",
       legend = c("SP > DP", "SP < DP"),
       col = c("red", "gold"),
       lwd = 3,
       bty = "n")

Ratio SP/DP :

ok <- !is.na(mois_num) &
  !is.na(annee) &
  !is.na(data$SGILOTSP) &
  !is.na(data$SGILOTDP) &
  data$SGILOTDP > 0

date_ratio <- as.Date(
  paste(annee[ok],
        sprintf("%02d", mois_num[ok]),
        "01",
        sep = "-"))

ratio <- data$SGILOTSP[ok] / data$SGILOTDP[ok]

plot(date_ratio,
     ratio,
     type = "o",
     pch = 16,
     col = "red",
     xlab = "Date",
     ylab = "Ratio OPG SP / DP",
     main = "Ratio d'infestation parasitaire SP / DP")

abline(h = 1, lty = 2, lwd = 2)

Infestation moyenne par mois :

mois_nom <- c("Jan", "Fev", "Mar", "Avr", "Mai", "Juin",
              "Juil", "Aout", "Sep", "Oct", "Nov", "Dec")

moy_SP <- numeric(12)
moy_DP <- numeric(12)

for (m in 1:12) {
  valeurs_SP <- data$SGILOTSP[mois_num == m]
  valeurs_DP <- data$SGILOTDP[mois_num == m]
  
  moy_SP[m] <- mean(valeurs_SP, na.rm = TRUE)
  moy_DP[m] <- mean(valeurs_DP, na.rm = TRUE)}

plot(1:12,
     moy_SP,
     type = "o",
     pch = 16,
     lwd = 2,
     col = "red",
     xaxt = "n",
     xlab = "Mois",
     ylab = "Nombre d'OPG moyen",
     main = "Nombre d'OPG moyen selon le mois")

lines(1:12,
      moy_DP,
      type = "o",
      pch = 17,
      lwd = 2,
      col = "gold")

axis(1, at = 1:12, labels = mois_nom)

legend("topright",
       legend = c("SP", "DP"),
       col = c("red", "gold"),
       pch = c(16, 17),
       lty = 1,
       bty = "n")

Heatmap : Infestation SP selon années et mois

annees_uniques <- sort(unique(annee[!is.na(annee)]))

mat_SP <- matrix(NA,
                 nrow = length(annees_uniques),
                 ncol = 12)

rownames(mat_SP) <- annees_uniques
colnames(mat_SP) <- mois_nom

for (i in 1:length(annees_uniques)) {
  
  for (m in 1:12) {
    
    valeurs <- data$SGILOTSP[
      annee == annees_uniques[i] &
        mois_num == m]
    
    if (sum(!is.na(valeurs)) > 0) {
      mat_SP[i, m] <- mean(valeurs, na.rm = TRUE) }}}

heatmap(mat_SP,
        Rowv = NA,
        Colv = NA,
        scale = "none",
        margins = c(6, 6),
        main = "Infestation moyenne du lot SP selon l'année et le mois")

Heatmap : Infestation DP selon années et mois

annees_uniques <- sort(unique(annee[!is.na(annee)]))

mat_DP <- matrix(NA,
                 nrow = length(annees_uniques),
                 ncol = 12)

rownames(mat_DP) <- annees_uniques
colnames(mat_DP) <- mois_nom

for (i in 1:length(annees_uniques)) {
  
  for (m in 1:12) {
    
    valeurs <- data$SGILOTSP[
      annee == annees_uniques[i] &
        mois_num == m]
    
    if (sum(!is.na(valeurs)) > 0) {
      mat_SP[i, m] <- mean(valeurs, na.rm = TRUE)}}}

heatmap(mat_SP,
        Rowv = NA,
        Colv = NA,
        scale = "none",
        margins = c(6, 6),
        main = "Infestation moyenne du lot DP selon l'année et le mois")

Analyses statistiques :

Analyses descriptives

resume <- function(x) {
  x <- x[is.finite(x)]
  c(
    Effectif = length(x),
    Moyenne = mean(x),
    Mediane = median(x),
    Ecart_type = sd(x),
    Minimum = min(x),
    Maximum = max(x))}

SP <- data$SGILOTSP[!is.na(data$SGILOTSP)]
DP <- data$SGILOTDP[!is.na(data$SGILOTDP)]
DB <- data$SGILOTDB[!is.na(data$SGILOTDB)]

cat("Lot SP\n")

Lot SP

print(resume(SP))

  Effectif    Moyenne    Mediane Ecart_type    Minimum    Maximum 
  150.0000   448.0747   350.0000   425.2526     0.0000  2428.0000 

cat("\nLot DP\n")

Lot DP

print(resume(DP))

  Effectif    Moyenne    Mediane Ecart_type    Minimum    Maximum 
  151.0000   180.9073   139.0000   161.9680     0.0000   700.0000 

cat("\nLot DB\n")

Lot DB

print(resume(DB))

  Effectif    Moyenne    Mediane Ecart_type    Minimum    Maximum 
 147.00000   10.84218    0.00000   27.59300    0.00000  122.00000 

Test de Wilcoxon

ok <- !is.na(data$SGILOTSP) &
  !is.na(data$SGILOTDP)

SP <- data$SGILOTSP[ok]
DP <- data$SGILOTDP[ok]

wilcox.test(SP, DP,
            paired = TRUE,
            exact = FALSE)

    Wilcoxon signed rank test with continuity correction

data:  SP and DP
V = 9350, p-value = 7.072e-15
alternative hypothesis: true location shift is not equal to 0