5  Plant phenotyping

#pkg
library(tidyverse)
library(here)
library(readxl)
library(ggh4x)
library(ggstats)
library(progressr)
library(viridis) 
library(patchwork)
library(ggplot2)
library(progress)


# src
source(here::here("src/function/stat_function/stat_analysis_main.R")) # for make plot 
source(here::here("src/function/fig_export.R")) # This function saves a given plot (plot_x) as both a PDF and a high-resolution PNG file at specified dimensions.
source(here::here("src/function/Evaluate_contrast.R"))

# cosmetics
sulfate_pallet=read_excel(here::here("data/color_palette.xlsm")) %>%
      filter(set == "sulfure_condition") %>%
      dplyr::select(color, treatment) %>%
      pull(color) %>%
      setNames(read_excel(here::here("data/color_palette.xlsm")) %>%
                 filter(set == "sulfure_condition") %>%
                 pull(treatment)
               )

pallet_genotype=read_excel(here::here("data/color_palette.xlsm")) %>%
      filter(set == "genotype") %>%
      dplyr::select(color, treatment) %>%
      pull(color) %>%
      setNames(read_excel(here::here("data/color_palette.xlsm")) %>%
                 filter(set == "genotype") %>%
                 pull(treatment)
               )
pallet_edaphic_condition=read_excel(here::here("data/color_palette.xlsm")) %>%
      filter(set == "edaphic_condition") %>%
      dplyr::select(color, treatment) %>%
      pull(color) %>%
      setNames(read_excel(here::here("data/color_palette.xlsm")) %>%
                 filter(set == "edaphic_condition") %>%
                 pull(treatment)
               )

pallet_compartment <- read_excel(here::here("data/color_palette.xlsm")) %>%
  filter(set == "compartment_plant", treatment %in% c("VL", "R")) %>%
  distinct(treatment, color) %>%        
  deframe()     

pallet_sampling <- read_excel(here::here("data/color_palette.xlsm")) %>%
  filter(set == "sampling") %>%
  distinct(treatment, color) %>%        
  deframe()  

• After phenotyping, put the correct plant name for each photo in the format Side_<plant_id>_20210426-083302_<imgangle>_GEAP207_U4__2L_2684_WW_SS_ECO_Emat_Rep1_1_22_1_2_0.png with this code
• Create a model with ilastic to identify what is plant and what is not
• Segment the image (plant vs. background) and put it in the segmented folder
• Then a script convert the png image with 1 and 2 (corresponding to background and plant) into image compose of RGB (0,0,0) or (255,255,255). For that is use this code.
• Then i tchek the number of pixel acording to the height to see were i can crop with this code.
• Then I created an algorithm that slices the segmented image according to certain pixel heights to clean up pot or stake errors (below and above the plant). For that is use this code.
• Make a convex_hull with my old code. See the code here.
• Then, we create background free image (plant withhout background) (in the folder extracted by using the previous image with this code
• Then, We used this image to extract RGB features only from plant pixels
• After segmenting the plant from the background, we computed color indices on the isolated plant region for all the plant. For that we use the folowing code
• RGB ratios were calculated from plant-only masks to avoid background influence
• We then did the same for each pixel of height of the plant (if, for example, it’s greener for the lower parts than for the upper parts). For that we use the folowing code

Here are the different stages and results in video

Original image

Segmentation post deep learning algorithm

Extraction of colour pixels only for what has been recognised as plant

Measure with the convex shape of the height and width

Measurement of the green, red and blue pixel quantities for each image.

5.1 Measurement to analyze post Python algorithm

• Average plant area (with average for each angle)
• Maximum plant area
• Average convex plant area (with average for each angle)
• Maximum convex plant area
• Average plant height (with average for each angle)
• Average plant width (with average for each angle)
• Max width for each plant

5.1.1 Data importation

• The results are on a hard disk because it takes up too much space. Intermediate results are well taken into account.

• Days after pollination is expressed relative to the date when 80 % of the plants had reached the pollination stage (19 April 2021). The DAS (days after sowing) values, however, differ between cultivars because the Kayanne seeds were disinfected and sown with inoculation on 9 March 2021, whereas Caméor was sown on 16 March 2021

• For rhizotubes they were all sown on 29 March 2021 and i have not found any polinisation dates

• For the analysis below, I only used the values for the 2L pots and the DAP.

# /!\ warning is very very long !!! aproximatly 50 min !!!!!!!!!!!!!!!!!!!!!!!!!!!!
# 3348.91 sec elapsed

## root repertory
phenotyping_dir <- here("data/physio/phenotyping/results/pea_raw_shoot_v2/GEAPS28_0325")
phenotyping_dir <- "E:/Dijon_FILEAS/GEAP_207_C6_U4"

##############test 
phenotyping_dir_test <- "E:/Dijon_FILEAS/GEAP_207_C6_U4/27322/pixelwise_summary.csv"
###################### test 
test = read_csv(file = phenotyping_dir_test, show_col_types = F)

as.numeric(test$MeanRGDiff[20:100])
## petite fonction utilitaire -----------------------------------------------
read_phenotyping <- function(pattern) {
  # Récupère tous les fichier(s) qui correspondent au nom passé en argument,
  # dans n'importe quel sous-dossier (recursive = TRUE)
  list.files(
    phenotyping_dir,
    pattern = pattern,
    recursive = TRUE,
    full.names = TRUE
  ) |>
    map_dfr(~ read_csv(.x, show_col_types = FALSE)
            |> mutate(taskid = basename(dirname(.x))))  # optionnel : identifie la date/lot
}

read_phenotyping <- function(pattern) {
  list.files(
    phenotyping_dir,
    pattern   = pattern,
    recursive = TRUE,
    full.names = TRUE
  ) |>
    purrr::map_dfr(~ {
      df <- readr::read_csv(.x, show_col_types = FALSE)

      # 1. change « , » -> « . » dans toutes les colonnes texte
      df <- dplyr::mutate(
        df,
        across(where(is.character), ~ stringr::str_replace_all(.x, ",", "."))
      )

      # 2. reconvertit automatiquement (nombres, dates, etc.)
      df <- readr::type_convert(df, locale = readr::locale(decimal_mark = "."))

      # 3. ajoute l’identifiant de lot
      df$taskid <- basename(dirname(.x))
      df
    })
}

read_phenotyping <- function(pattern, show_progress) {
  # 1. liste complète des fichiers à traiter
  paths <- list.files(
    phenotyping_dir,
    pattern   = pattern,
    recursive = TRUE,
    full.names = TRUE
  )

  # 2. initialise la barre si demandé
  if (show_progress) {
    pb <- progress::progress_bar$new(
      total   = length(paths),
      format  = " Import [:bar] :current/:total (:percent) | écoulé: :elapsed | restant: :eta",
      clear   = FALSE,   # garde la barre à la fin
      width   = 60
    )
  }

  # 3. boucle de lecture avec mise à jour de la barre
  purrr::map_dfr(paths, ~ {
    if (show_progress) pb$tick()

    df <- readr::read_csv(.x, show_col_types = FALSE)

    df <- dplyr::mutate(
      df,
      across(where(is.character), ~ stringr::str_replace_all(.x, ",", "."))
    )
    df <- readr::type_convert(df, locale = readr::locale(decimal_mark = "."))
    df$taskid <- basename(dirname(.x))
    df
  })
}

## import identique ----------------------------------------------------------
tictoc::tic()
df_black  <- read_phenotyping("result_black_pixel_corrected\\.csv$", show_progress = TRUE) |>
  mutate(Label = str_remove(Label, "_Simple Segmentation_segmented$")) |>
  dplyr::rename(surface = BlackPixels)

df_pixel  <- read_phenotyping("pixelwise_summary\\.csv$", show_progress = TRUE) |>
  mutate(Label = str_remove(Label, "_extracted$"))

df_convex <- read_phenotyping("result_convex_hull\\.csv$", show_progress = TRUE) |>
  mutate(Label = str_remove(Label, "_Simple Segmentation_segmented_cor$")) |>
  dplyr::rename(height = profondeur,
         width  = largeur) |>
  dplyr::select(-num_label)

levels(as.factor(df_convex$taskid))
tictoc::toc()
## import size px for each taskid
#df_px <- read_excel(here::here("data/physio/phenotyping/conversion_px_cm.xlsx"), col_names = TRUE)
df_px <- read.csv(here::here("data/physio/phenotyping/coefficient_mires_side.csv"), sep=";")

## pipeline de fusion --------------------------------------------------------
df_global <- df_black |> 
  left_join(df_pixel,  by = c("Label", "taskid")) |>    # ← ajoute batch
  left_join(df_convex, by = c("Label", "taskid")) |>    # ← ajoute batch
  mutate(
    id        = str_extract(Label, "(?<=Side_)\\d+"),
    angle     = stringr::str_extract(Label,"\\d+(?=_GEAP207_U4_)"),
    timestamp = stringr::str_extract(Label, "\\d{8}-\\d{6}"),
    date      = stringr::str_extract(timestamp, "^\\d{8}"),
    time      = stringr::str_extract(timestamp, "(?<=-)\\d{6}"),
    rt2l  = stringr::str_extract(Label, "(?<=_)(RT|2L)(?=_)"),
    plant_num = ifelse(rt2l == "RT", as.numeric(id)-5000, as.numeric(id)), 
  ) |>
  left_join(
    read_excel(here("data/GEAPSI_EUCLEG_20210426.xlsx"), col_names = TRUE, sheet = "ID_U4") |>
      dplyr::rename(plant_num=pot, 
                    water_condition = TRAITEMENT, 
                    sulfur_condition = Nutrition,
                    genotype = GENO, 
                    line = Ligne,
                    row = Cars, 
                    sampling = PVT
                    ) %>% 
      mutate(position = paste0(line,"_", row),
             edaphic_condition = paste0(water_condition, "_", sulfur_condition), 
             condition = paste0(genotype, "_", edaphic_condition)),
    by = "plant_num"
  ) %>% 
  dplyr::mutate(genotype = case_when(
        genotype %in% "2684" ~ "W78*",
        genotype %in% "4693" ~ "E568K",
        genotype %in% "CAM2684" ~ "WT1",
        genotype %in% "CAM4693" ~ "WT2", 
        genotype %in% "KAY" ~ "KAY"), 
        DAP = ifelse(rt2l == "2L" , as.Date(date, format = "%Y%m%d")-as.Date("20210419", format = "%Y%m%d"), NA), 
        DAS = ifelse(rt2l == "2L" , ifelse(genotype == "KAY", as.Date(date, format = "%Y%m%d")-as.Date("20210309", format = "%Y%m%d"), as.Date(date, format = "%Y%m%d")-as.Date("20210316", format = "%Y%m%d")), as.Date(date, format = "%Y%m%d")-as.Date("20210329", format = "%Y%m%d"))
        ) %>% 
  relocate(
    plant_num, genotype, row, line, condition,sampling,
    id, angle, timestamp, date, time,
    .before = Label
  ) %>% 
#   mutate(genotype = fct_relevel(genotype, "WT1", "W78*", "WT2", "E568K")) |>
#   filter(area > 2e6) %>% 
#   filter(plant_num!= 6) %>% 
#   mutate(DAP = as.numeric(as.Date(date)-as.Date("2025-04-10"))) %>% 
#   mutate(taskid = as.integer(taskid)) %>% 
   left_join(., df_px %>% mutate(taskid = as.character(taskid)), by ="taskid")  %>% 
   mutate(surface_cm_2 = surface * curxscale*curyscale/100, 
         perimeter_cm = perimeter * ((curxscale + curyscale)/2)/10,
         area_cm_2 = area * curxscale*curyscale/100,
         height_cm = height * curyscale/10,
         width_cm = width * curxscale/10)

   # export
write_csv(df_global, here::here("data/physio/phenotyping/df_global.csv"))

5.1.2 Data representation

The average value of all angels was been taken.

5.1.2.1 Examples of interesting ratios to test

• G / (R + B) Green index
• R / G Reddening / senescence tendency
• G / (R + G + B) Percentage of green in the image
• R - G Red-green difference
• (G - R) / (G + R) Normalized Green Red Difference Index (NGRDI) Good indicator of health or stress
• (G - B) / (G + B) Vegetative Index (VI) Often used in simple cases
• -0.5 * [(190*(R-G)) - (120*(R-B))]TGI (Triangular Greenness Index) Estimates chlorophyll without a spectrometer

df_global = read.csv(here::here("data/physio/phenotyping/df_global.csv"),dec = ",",
                     ) %>%
  mutate(genotype=fct_relevel(genotype, "KAY","WT1", "W78*", "WT2", "E568K"),
         across(c(height_cm, width_cm, surface_cm_2, MeanGI, MeanVI, MeanRedden, MeanPctGreen, MeanRGDiff, MeanNGRDI), ~ as.numeric (.x))
  ) %>% 
  filter(rt2l == "2L") %>% 
  filter(genotype != "SP") 

# df_global %>% 
#   filter(taskid == "27489") %>% 
#   ggplot(aes(x = as.factor(plant_num), y = as.numeric(height_cm), fill = edaphic_condition)) +
#   geom_boxplot()+
#   facet_grid(genotype ~ taskid)+
#   scale_color_manual(values = pallet_edaphic_condition) +
#   scale_fill_manual(values = pallet_edaphic_condition) +
#   theme(axis.text.x = element_text(angle = 0, hjust = 1))

# df_global %>% 
#   ggplot(aes(x = as.factor(plant_num), y = MeanGI, fill = edaphic_condition)) +
#   geom_boxplot()+
#   facet_grid(genotype ~ taskid)+
#   scale_color_manual(values = pallet_edaphic_condition) +
#   scale_fill_manual(values = pallet_edaphic_condition)
# 
# df_global %>% 
#   ggplot(aes(x = as.factor(plant_num), y = MeanVI, fill = edaphic_condition)) +
#   geom_boxplot()+
#   facet_grid(genotype ~ taskid)+
#   scale_color_manual(values = pallet_edaphic_condition) +
#   scale_fill_manual(values = pallet_edaphic_condition)

cols_to_average <- c("surface_cm_2", "MeanR", "MeanG", "MeanB", "MeanGI", "MeanRedden", "MeanPctGreen",  "MeanRGDiff", "MeanNGRDI", "MeanVI", "MeanTGI", "perimeter_cm", "area_cm_2", "height_cm", "width_cm")  # remplace par tes vraies colonnes

df_mean <- df_global %>%
  dplyr::group_by(plant_num, genotype, row, line, position, water_condition, sulfur_condition, edaphic_condition, condition, date, sampling, DAP, DAS) %>%
  dplyr::summarise(
    across(
      all_of(cols_to_average),
      list(mean = ~ mean(.x, na.rm = TRUE),
           sd = ~ sd(.x, na.rm = TRUE)),
      .names = "{.fn}_{.col}"
    ),
    .groups = "drop"
  ) %>% mutate(
    genotype = forcats::fct_relevel(genotype, "KAY", "WT1", "W78*", "WT2", "E568K"),
    water_condition = forcats::fct_relevel(water_condition, "WW", "WS"),
    sulfur_condition = forcats::fct_relevel(sulfur_condition, "SS", "SD"),
    sampling = forcats::fct_relevel(sampling, "E0", "E1", "E2"),
    edaphic_condition = forcats::fct_relevel(edaphic_condition, "WW_SS", "WW_SD", "WS_SS", "WS_SD")
  )

# %>%
  # mutate(
  #   DAP_text = factor(
  #     paste0("DAP: ", DAP),
  #     levels = paste0("DAP: ", sort(unique(DAP)))
  #   )
  # ) %>% 
  
# I take the most extreme value to measure the surface area of the plant. 

p_surface = df_mean %>% 
  ggplot(aes(x = edaphic_condition, y = mean_surface_cm_2, fill = edaphic_condition, linetype = genotype)) +
  geom_boxplot(alpha = .6)+
  facet_grid(. ~ paste0("DAP: ",DAP)) + 
  #facet_grid(. ~ paste0("DAS:",DAS), scales = "free_x", space = "free_x") + 
  theme_bw()+
  scale_fill_manual(values = pallet_edaphic_condition)+
  labs(fill = "Edaphic Condition")+
  theme(
    axis.text.x = element_text(angle = 90, vjust = 0.5, size = 12),
    axis.text.y = element_text(size = 12),
    strip.text = element_text(size = 12),
    legend.position = "bottom",
    legend.title = element_text(size = 13),
    legend.text = element_text(size = 14)
  )+
  labs(y = "Surface (cm²)", x = "Edaphic Condition", linetype = "Genotype") ; p_surface

fig_export(here::here(paste0("report/physio/plot/phenotyping/cinetic_surface_DAP")), p_surface, height_i = 12, width_i = 30, res_i = 400)

# with DAS
p_surface = df_mean %>% 
  ggplot(aes(x = edaphic_condition, y = mean_surface_cm_2, fill = edaphic_condition, linetype = genotype)) +
  geom_boxplot(alpha = .6)+
  #facet_grid(. ~ paste0("DAP: ",DAP)) + 
  facet_grid(. ~ paste0("DAS:",DAS), scales = "free_x", space = "free_x") + 
  theme_bw()+
  scale_fill_manual(values = pallet_edaphic_condition)+
  labs(fill = "Edaphic Condition")+
  theme(
    axis.text.x = element_text(angle = 90, vjust = 0.5, size = 12),
    axis.text.y = element_text(size = 12),
    strip.text = element_text(size = 12),
    legend.position = "bottom",
    legend.title = element_text(size = 13),
    legend.text = element_text(size = 14)
  )+
  labs(y = "Surface (cm²)", x = "Edaphic Condition", linetype = "Genotype") ; p_surface

fig_export(here::here(paste0("report/physio/plot/phenotyping/cinetic_surface_DAS")), p_surface, height_i = 12, width_i = 30, res_i = 400)

p_surf_dap <- df_mean %>% 
  ggplot(aes(x = DAP, y = mean_surface_cm_2, color = edaphic_condition, linetype = genotype)) +
  geom_smooth(alpha = .2)+
  theme_bw()+
  scale_color_manual(values = pallet_edaphic_condition) +
  labs(color = "Edaphic Condition", y = "Surface (cm²)", linetype = "Genotype")

p_surf_das <- df_mean %>% 
  ggplot(aes(x = DAS, y = mean_surface_cm_2, color = edaphic_condition, linetype = genotype)) +
  geom_smooth(alpha = .2)+
  theme_bw()+
  scale_color_manual(values = pallet_edaphic_condition) +
  labs(color = "Edaphic Condition", y = "Surface (cm²)", linetype = "Genotype")

p_combine <- p_surf_dap + p_surf_das + plot_layout(guides = "collect")
fig_export(here::here(paste0("report/physio/plot/phenotyping/cinetic_surface_line_DAP_DAS")), p_combine, height_i = 6, width_i = 10, res_i = 500, format = "png")

################################### SURFACE by genotype ######################### witth stats 
#by genotype for surface #############################
p_surf_dap_g <- df_mean %>% 
  ggplot(aes(x = DAP, y = mean_surface_cm_2, color = genotype, fill = genotype)) +
  geom_smooth(alpha = .2)+
  theme_bw()+
  scale_color_manual(values = pallet_genotype) +
  scale_fill_manual(values = pallet_genotype) +
  #facet_wrap(.~edaphic_condition, nrow = 2)+
  ggh4x::facet_wrap2(
    ~ edaphic_condition,
    #scales = "free_y",
    ncol   = 1,                                # ajuste le layout
    strip  = ggh4x::strip_themed(
      background_x = ggh4x::elem_list_rect(fill = pallet_edaphic_condition[c("WW_SS", "WW_SD", "WS_SS", "WS_SD")]),
      text_x       = ggh4x::elem_list_text(colour = "white")
    )
  )+
  labs(y = "Surface (cm²)", fill = "Genotype", color = "Genotype") +
  theme(legend.position = "none"); p_surf_dap_g

p_surf_das_g <- df_mean %>% 
  ggplot(aes(x = DAS, y = mean_surface_cm_2, color = genotype, fill = genotype)) +
  geom_smooth(alpha = .25)+
  theme_bw()+
  scale_color_manual(values = pallet_genotype) +
  scale_fill_manual(values = pallet_genotype) +
  #facet_wrap(.~edaphic_condition, nrow = 2)+
  ggh4x::facet_wrap2(
    ~ edaphic_condition,
    #scales = "free_y",
    ncol   = 1,                                # ajuste le layout
    strip  = ggh4x::strip_themed(
      background_x = ggh4x::elem_list_rect(fill = pallet_edaphic_condition[c("WW_SS", "WW_SD", "WS_SS", "WS_SD")]),
      text_x       = ggh4x::elem_list_text(colour = "white")
    )
  )+
  labs(y = "Surface (cm²)", fill = "Genotype", color = "Genotype") + 
  theme(legend.position = "none"); p_surf_das_g

p_kinetic <- p_surf_das_g + p_surf_dap_g

##### make stats for each edaphic condition at 7 DAP ################################################
v_edaphic_condition <- levels(as.factor(df_mean$edaphic_condition))

df_mean_select = df_mean %>%  
  mutate(
    edaphic_condition_color = unname(pallet_edaphic_condition[edaphic_condition]), 
    couleur_texte_perso = sapply(edaphic_condition_color, evaluate_contrast), 
    couleur_texte_perso = "white"
    ) %>% 
  drop_na(mean_surface_cm_2) %>% 
  filter(DAP == 7) ##################################### !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

plots_list <- lapply(seq_along(v_edaphic_condition), function(i) {
  
  edaphic_condition_select <- v_edaphic_condition[i]
  
  p <- stat_analyse(
    data = df_mean_select %>% 
      as.data.frame() %>% 
      filter(edaphic_condition == edaphic_condition_select) , # !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    
    column_value     = "mean_surface_cm_2",
    category_variables = "genotype",
    grp_var          = "",
    show_plot        = TRUE,
    outlier_show     = FALSE,
    label_outlier    = "plant_num",
    biologist_stats  = TRUE,
    
    Ylab_i = paste0(
      ""
    ),
    
    control_conditions = "",
    strip_normale      = TRUE,
    hex_pallet         = pallet_genotype,
    strip              = "edaphic_condition",
    strip_fill_vec     = "edaphic_condition_color",
    strip_text_vec     = "couleur_texte_perso"
  )[["plot"]]+
    theme_bw()+
    theme(axis.text.x = element_text(angle = 90, hjust = 1))
  
  ## 1-a) On retire les titres d’axe des sous-graphiques
  p <- p + labs(x = NULL, y = NULL)
  
  ## 1-b) Option : on masque aussi les graduations redondantes
  #       (ici : on garde y seulement sur la 1re colonne
  #              et x seulement sur la dernière ligne)
  ncol_layout <- 1
  nrow_layout <- ceiling(length(v_edaphic_condition) / ncol_layout)
  col_i <- ((i - 1) %% ncol_layout) + 1
  row_i <- ceiling(i / ncol_layout)
  
  # if (col_i != 1) {
  #   p
  # }
  if (row_i != nrow_layout) {
    p <- p + theme(axis.text.x  = element_blank())
  }
  
  p +labs(fill = "Genotype", colour = "Genotype")
})

# ------------------------------------------------------------------
# 2) Mosaïque principale -------------------------------------------
# ------------------------------------------------------------------
panel <- wrap_plots(plots_list, ncol = 1, guides = "collect") &
         theme(legend.position = "right", 
               panel.grid.major = element_blank(),
               panel.grid.minor = element_blank()
               )

## ------------
## 2. Labels
## ------------
y_lab <- ggplot() +
         labs(y = "Surface (cm²) at 7 DAP ") +
         theme_void() +
         theme(axis.title.y = element_text(angle = 90, vjust = 0.5,
                                           hjust = 0.5))

x_lab <- ggplot() +
         labs(x = "Condition") +
         theme_void() +
         theme(axis.title.x = element_text(vjust = -0.5, hjust = 0.5))

## ------------
## 3. Assemblage 2 × 2
## ------------
final_plot_1 <- (y_lab  | panel)+
  plot_layout(
    widths  = c(0.02, 0.98),   # 5 % pour la colonne Y, 95 % pour le panel
    guides  = "collect"
  )

final_plot_2 <- (plot_spacer() | x_lab)+
  plot_layout(
    widths  = c(0.02, 0.98),   # 5 % pour la colonne Y, 95 % pour le panel
    guides  = "collect"
  )

p_stats <- (final_plot_1/final_plot_2)+
  plot_layout(
    heights = c(0.999, 0.001),   # 7 % pour le label X
    guides  = "collect"
  )

combine_plot <-  (p_kinetic|p_stats)+
  plot_layout(
    widths = c(0.65, 0.35),   
    guides  = "collect"
  ) 

fig_export(here::here(paste0("report/physio/plot/phenotyping/cinetic_surface_line_DAP_DAS_genotype_stats")), combine_plot, height_i = 9, width_i = 9, res_i = 800, format = "png")




################################### SURFACE by edphyque condition ######################### witth stats 
#by genotype for surface #############################
p_surf_dap_ec <- df_mean %>% 
  ggplot(aes(x = DAP, y = mean_surface_cm_2, color = edaphic_condition, fill = edaphic_condition)) +
  geom_smooth(alpha = .2)+
  theme_bw()+
  scale_color_manual(values = pallet_edaphic_condition) +
  scale_fill_manual(values = pallet_edaphic_condition) +
  #facet_wrap(.~edaphic_condition, nrow = 2)+
  ggh4x::facet_wrap2(
    ~ genotype,
    #scales = "free_y",
    ncol   = 1,                                # ajuste le layout
    strip  = ggh4x::strip_themed(
      background_x = ggh4x::elem_list_rect(fill = pallet_genotype),
      text_x       = ggh4x::elem_list_text(colour = "white")
    )
  )+
  labs(y = "Surface (cm²)", fill = "Edaphic condition", color = "Edaphic condition") +
  theme(legend.position = "none"); p_surf_dap_ec

p_surf_das_ec <- df_mean %>% 
  ggplot(aes(x = DAS, y = mean_surface_cm_2, color = edaphic_condition, fill = edaphic_condition)) +
  geom_smooth(alpha = .25)+
  theme_bw()+
  scale_color_manual(values = pallet_edaphic_condition) +
  scale_fill_manual(values = pallet_edaphic_condition) +
  #facet_wrap(.~edaphic_condition, nrow = 2)+
  ggh4x::facet_wrap2(
    ~ genotype,
    #scales = "free_y",
    ncol   = 1,                                # ajuste le layout
    strip  = ggh4x::strip_themed(
      background_x = ggh4x::elem_list_rect(fill = pallet_genotype),
      text_x       = ggh4x::elem_list_text(colour = "white")
    )
  )+
  labs(y = "Surface (cm²)", fill = "Genotype", color = "Genotype") + 
  theme(legend.position = "none"); p_surf_das_ec

p_kinetic_ec <- p_surf_das_ec + p_surf_dap_ec

##### make stats for each genotype at 7 DAP ################################################
v_genotype <- levels(as.factor(df_mean$genotype))

df_mean_select = df_mean %>%  
  mutate(
    genotype_color = unname(pallet_genotype[genotype]), 
    couleur_texte_perso = sapply(genotype_color, evaluate_contrast), 
    couleur_texte_perso = "white"
    ) %>% 
  drop_na(mean_surface_cm_2) %>% 
  filter(DAP == 7) ##################################### !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

plots_list <- lapply(seq_along(v_genotype), function(i) {
  
  genotype_select <- v_genotype[i]
  
  p <- stat_analyse(
    data = df_mean_select %>% 
      as.data.frame() %>% 
      filter(genotype == genotype_select) , # !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    
    column_value     = "mean_surface_cm_2",
    category_variables = "edaphic_condition",
    grp_var          = "",
    show_plot        = TRUE,
    outlier_show     = FALSE,
    label_outlier    = "plant_num",
    biologist_stats  = TRUE,
    
    Ylab_i = paste0(
      ""
    ),
    
    control_conditions = "",
    strip_normale      = TRUE,
    hex_pallet         = pallet_edaphic_condition,
    strip              = "genotype",
    strip_fill_vec     = "genotype_color",
    strip_text_vec     = "couleur_texte_perso"
  )[["plot"]]+
    theme_bw()+
    theme(axis.text.x = element_text(angle = 90, hjust = 1))
  
  ## 1-a) On retire les titres d’axe des sous-graphiques
  p <- p + labs(x = NULL, y = NULL)
  
  ## 1-b) Option : on masque aussi les graduations redondantes
  #       (ici : on garde y seulement sur la 1re colonne
  #              et x seulement sur la dernière ligne)
  ncol_layout <- 1
  nrow_layout <- ceiling(length(v_edaphic_condition) / ncol_layout)
  col_i <- ((i - 1) %% ncol_layout) + 1
  row_i <- ceiling(i / ncol_layout)
  
  # if (col_i != 1) {
  #   p
  # }
  if (row_i != nrow_layout) {
    p <- p + theme(axis.text.x  = element_blank())
  }
  
  p +labs(fill = "Edaphic condition", colour = "Edaphic condition")
})

# ------------------------------------------------------------------
# 2) Mosaïque principale -------------------------------------------
# ------------------------------------------------------------------
panel <- wrap_plots(plots_list, ncol = 1, guides = "collect") &
         theme(legend.position = "right", 
               panel.grid.major = element_blank(),
               panel.grid.minor = element_blank()
               )

## ------------
## 2. Labels
## ------------
y_lab <- ggplot() +
         labs(y = "Surface (cm²) at 7 DAP ") +
         theme_void() +
         theme(axis.title.y = element_text(angle = 90, vjust = 0.5,
                                           hjust = 0.5))

x_lab <- ggplot() +
         labs(x = "Edaphic condition") +
         theme_void() +
         theme(axis.title.x = element_text(vjust = -0.5, hjust = 0.5))

## ------------
## 3. Assemblage 2 × 2
## ------------
final_plot_1 <- (y_lab  | panel)+
  plot_layout(
    widths  = c(0.02, 0.98),   # 5 % pour la colonne Y, 95 % pour le panel
    guides  = "collect"
  )

final_plot_2 <- (plot_spacer() | x_lab)+
  plot_layout(
    widths  = c(0.02, 0.98),   # 5 % pour la colonne Y, 95 % pour le panel
    guides  = "collect"
  )

p_stats <- (final_plot_1/final_plot_2)+
  plot_layout(
    heights = c(0.999, 0.001),   # 7 % pour le label X
    guides  = "collect"
  )

combine_plot <-  (p_kinetic_ec|p_stats)+
  plot_layout(
    widths = c(0.65, 0.35),   
    guides  = "collect"
  ) 

fig_export(here::here(paste0("report/physio/plot/phenotyping/cinetic_surface_line_DAP_DAS_ec_stats")), combine_plot, height_i = 9, width_i = 9, res_i = 800, format = "png")


#by edaphic condition for mean_MeanVI #############################
p_vi_dap_ec <- df_mean %>% 
  ggplot(aes(x = DAP, y = mean_MeanVI, color = edaphic_condition, fill = edaphic_condition)) +
  geom_smooth(alpha = .2)+
  theme_bw()+
  scale_color_manual(values = pallet_edaphic_condition) +
  scale_fill_manual(values = pallet_edaphic_condition) +
  #facet_wrap(.~edaphic_condition, nrow = 2)+
  ggh4x::facet_wrap2(
    ~ genotype,
    #scales = "free_y",
    nrow   = 1,                                # ajuste le layout
    strip  = ggh4x::strip_themed(
      background_x = ggh4x::elem_list_rect(fill = pallet_genotype),
      text_x       = ggh4x::elem_list_text(colour = "white")
    )
  )+
  theme(legend.position = "bottom")+
  labs(y = "Vegetative index", fill = "Edaphic Condition", color = "Edaphic Condition") ; p_vi_dap_ec

fig_export(here::here(paste0("report/physio/plot/phenotyping/cinetic_VI_DAP_edaphic_condition")), p_vi_dap_ec, height_i = 4, width_i = 8, res_i = 800, format = "png")

#by genotype  for mean_MeanVI #############################
p_vi_dap_g <- df_mean %>% 
  ggplot(aes(x = DAP, y = mean_MeanVI, color = genotype, fill = genotype)) +
  geom_smooth(alpha = .2)+
  theme_bw()+
  scale_color_manual(values = pallet_genotype) +
  scale_fill_manual(values = pallet_genotype) +
  #facet_wrap(.~edaphic_condition, nrow = 2)+
  ggh4x::facet_wrap2(
    ~ edaphic_condition,
    #scales = "free_y",
    nrow   = 1,                                # ajuste le layout
    strip  = ggh4x::strip_themed(
      background_x = ggh4x::elem_list_rect(fill = pallet_edaphic_condition[c("WW_SS", "WW_SD", "WS_SS", "WS_SD")]),
      text_x       = ggh4x::elem_list_text(colour = "white")
    )
  )+
  theme(legend.position = "bottom")+
  labs(y = "Vegetative index", fill = "Genotype", color = "Genotype") ; p_vi_dap_g

fig_export(here::here(paste0("report/physio/plot/phenotyping/cinetic_VI_DAP_genotype")), p_vi_dap_g, height_i = 4, width_i = 8, res_i = 800, format = "png")


#by genotype  for REEDEN #############################
p_red_dap_g <- df_mean %>% 
  ggplot(aes(x = DAP, y = mean_MeanRedden, color = genotype, fill = genotype)) +
  geom_smooth(alpha = .2)+
  theme_bw()+
  scale_color_manual(values = pallet_genotype) +
  scale_fill_manual(values = pallet_genotype) +
  #facet_wrap(.~edaphic_condition, nrow = 2)+
  ggh4x::facet_wrap2(
    ~ edaphic_condition,
    #scales = "free_y",
    nrow   = 1,                                # ajuste le layout
    strip  = ggh4x::strip_themed(
      background_x = ggh4x::elem_list_rect(fill = pallet_edaphic_condition[c("WW_SS", "WW_SD", "WS_SS", "WS_SD")]),
      text_x       = ggh4x::elem_list_text(colour = "white")
    )
  )+
  theme(legend.position = "bottom")+
  labs(y = "Redden index", fill = "Genotype", color = "Genotype") ; p_red_dap_g

fig_export(here::here(paste0("report/physio/plot/phenotyping/cinetic_Redden_DAP_genotype")), p_red_dap_g, height_i = 4, width_i = 8, res_i = 800, format = "png")


p_red_dap_ec <- df_mean %>% 
  ggplot(aes(x = DAP, y = mean_MeanRedden, color = edaphic_condition, fill = edaphic_condition)) +
  geom_smooth(alpha = .2)+
  theme_bw()+
  scale_color_manual(values = pallet_edaphic_condition) +
  scale_fill_manual(values = pallet_edaphic_condition) +
  #facet_wrap(.~edaphic_condition, nrow = 2)+
  ggh4x::facet_wrap2(
    ~ genotype,
    #scales = "free_y",
    nrow   = 1,                                # ajuste le layout
    strip  = ggh4x::strip_themed(
      background_x = ggh4x::elem_list_rect(fill = pallet_genotype),
      text_x       = ggh4x::elem_list_text(colour = "white")
    )
  )+
  theme(legend.position = "bottom")+
  labs(y = "Redden index", fill = "Edaphic Condition", color = "Edaphic Condition") ; p_red_dap_ec

fig_export(here::here(paste0("report/physio/plot/phenotyping/cinetic_Redden_DAP_edaphic_condition")), p_red_dap_ec, height_i = 4, width_i = 8, res_i = 800, format = "png")

Surface of the plant by deep learning

Shows the intensity of green

Redden index a proxy of senescence