################################################################################################################### ################# REFACTO COPILOT A VALIDER pour double boucle for sur les vaches ####################### ################################################################################################################### # ========================= # Helpers (normalisations) # ========================= # Normalisation vs stats cheptel: 1 = max norm_chep <- function(x, var) { r <- stats_chep %>% dplyr::filter(var == !!var) if (nrow(r) == 0) return(rep(NA_real_, length(x))) mmin <- r$min[1]; mmax <- r$max[1] 1 - (abs(mmax - x) / abs(mmax - mmin)) } # ========================= # 1) Synthèse veaux par mère # ========================= veaux_summary <- PROD %>% group_by(mere) %>% summarise( n_veaux = n(), # plage de campagnes velages (max - min + 1), robustes aux NA campn_min = {cn <- campn[!is.na(campn)]; if (length(cn)) min(cn) else NA_integer_}, campn_max = {cn <- campn[!is.na(campn)]; if (length(cn)) max(cn) else NA_integer_}, nbcampvel = ifelse(!is.na(campn_min) & !is.na(campn_max), campn_max - campn_min + 1, NA_integer_), # age au 1er vêlage et IVV1 (première occurrence) agevel1 = {x <- agevel[ravelamere == 1]; if (length(x)) x[1] else NA_real_}, ivv1 = {x <- ivv[ravelamere == 2]; if (length(x)) x[1] else NA_real_}, ivv2p = round(mean(ivv[ravelamere > 2], na.rm = TRUE), 1), # date du dernier anais (dernier événement) last_danais = {d <- danais[!is.na(danais)]; if (length(d)) max(d) else as.Date(NA)}, # moyennes nécessaires pour ptgP et précocité mean_devmus = mean(devmus, na.rm = TRUE), mean_devsqe = mean(devsqe, na.rm = TRUE), mean_diff_dev = mean(devsqe - devmus, na.rm = TRUE), # indicateurs produits mort = round(sum(mortsev == "O", na.rm = TRUE) / n_veaux * 100, 1), txrepros = round(sum(repro == "O", na.rm = TRUE) / n_veaux * 100, 1), nbpp = sum(NBPRODIPG, na.rm = TRUE), nbpp_corr = sum(nbpp_corr, na.rm = TRUE), txmales = round(sum(sexbov == "1", na.rm = TRUE) / n_veaux * 100, 1), txvf = round(sum(conais %in% c("1","2"), na.rm = TRUE) / n_veaux * 100, 1), # moyennes par sexe pn_m = round(mean(ponais[sexbov == "1"], na.rm = TRUE), 1), p120_m = round(mean(pat04m[sexbov == "1"], na.rm = TRUE), 1), p210_m = round(mean(pat07m[sexbov == "1"], na.rm = TRUE), 1), pn_f = round(mean(ponais[sexbov == "2"], na.rm = TRUE), 1), p120_f = round(mean(pat04m[sexbov == "2"], na.rm = TRUE), 1), p210_f = round(mean(pat07m[sexbov == "2"], na.rm = TRUE), 1), # versions corrigées pn_corr = round(mean(pn_corr, na.rm = TRUE), 1), p120_corr = round(mean(p120_corr, na.rm = TRUE), 1), p210_corr = round(mean(p210_corr, na.rm = TRUE), 1), .groups = "drop" ) # ========================= # 2) Jointure & calculs vaches # ========================= v_ref <- vaches %>% left_join(veaux_summary, by = c("anim" = "mere")) %>% mutate( # pointage adulte synthétique ptgV = ifelse(!is.na(dmC), round(0.6 * dmC + 0.15 * ds + 0.25 * af, 1), NA_real_), # précocité & alpha precocite = ifelse(!is.na(n_veaux) & n_veaux > 3, round(mean_diff_dev, 2), NA_real_), alpha = ifelse(is.na(precocite), 1.62, 1.62 - 0.01 * precocite), # estimation du poids adulte (pad) pad = dplyr::case_when( !is.na(pat24m) ~ round((pat24m - 50 * exp(-720 * alpha * 10^(-3))) / (1 - exp(-720 * alpha * 10^(-3))), 1), !is.na(pat18m) ~ round((pat18m - 50 * exp(-540 * alpha * 10^(-3))) / (1 - exp(-540 * alpha * 10^(-3))), 1), !is.na(pat12m) ~ round((pat12m - 50 * exp(-360 * alpha * 10^(-3))) / (1 - exp(-360 * alpha * 10^(-3))), 1), TRUE ~ NA_real_ ), # temps improductif (e2, e3, e4) e2 = dplyr::case_when( is.na(ivv1) ~ 0, ivv1 < 390 ~ 0, TRUE ~ ivv1 - 390 ), e3 = dplyr::case_when( is.na(ivv2p) ~ 0, ivv2p < 365 ~ 0, TRUE ~ ivv2p - 365 ), days_since_last = as.numeric(difftime(Sys.Date(), last_danais, units = "days")), e4 = dplyr::case_when( is.na(days_since_last) ~ 0, days_since_last < 365 ~ 0, TRUE ~ days_since_last - 365 ), # temps productif (%) tempsprod = round( (age_days - (agevel1 * 30.4 + e2 + e3 * (nbcampvel - 2) + e4)) / age_days * 100, 1 ), # pointage produits ptgP = round(0.75 * mean_devmus + 0.25 * mean_devsqe, 1) ) %>% # Conversion des NaN en NA sur certaines moyennes mutate(across( c(ptgP, pn_m, pn_f, pn_corr, p120_m, p120_f, p120_corr, p210_m, p210_f, p210_corr), ~ ifelse(is.nan(.), NA_real_, .) )) # ========================= # 3) Normalisations # ========================= v_norm <- v_ref %>% mutate( # age au 1er vêlage normalisé agevel1_n = dplyr::case_when( is.na(agevel1) ~ NA_real_, agevel1 > 48 ~ 0, TRUE ~ round( -2 * (10^(-6)) * (agevel1 * 30.4)^2 + 0.0027 * (agevel1 * 30.4) + 8 * (10^(-15)), 3 ) ), # ivv1 normalisé ivv1_n = dplyr::case_when( is.na(ivv1) ~ NA_real_, ivv1 > 460 ~ 0, ivv1 < 390 ~ 1, TRUE ~ round(1 - abs(390 - ivv1) / abs(390 - 460), 3) ), # ivv2+ normalisé ivv2p_n = dplyr::case_when( is.na(ivv2p) | is.nan(ivv2p) ~ NA_real_, ivv2p > 435 ~ 0, ivv2p < 365 ~ 1, TRUE ~ round(1 - abs(365 - ivv2p) / abs(365 - 435), 3) ), # normalisations "cheptel 1 = max" pad_n = round(norm_chep(pad, "vaches$pad"), 3), ptgv_n = round(norm_chep(ptgV, "vaches$ptgV"), 3), txvf_n = round(norm_chep(txvf, "vaches$txvf"), 3), txm_n = round(norm_chep(txmales, "vaches$txmales"), 3), txrepros_n = round(norm_chep(txrepros, "vaches$txrepros"), 3), nbpp_n = round(norm_chep(nbpp_corr, "vaches$nbpp_corr"), 3), ptgp_n = round(norm_chep(ptgP, "vaches$ptgP"), 3), p120_n = round(norm_chep(p120_corr, "vaches$p120_corr"), 3), p210_n = round(norm_chep(p210_corr, "vaches$p210_corr"), 3), # prolificité prol_n = dplyr::case_when( is.na(prol) ~ NA_real_, prol >= 100 ~ 1, prol < 50 ~ 0, TRUE ~ round(1 - (abs(100 - prol) / abs(100 - 50)), 3) ), # poids naissance corrigé pn_n = dplyr::case_when( is.na(pn_corr) ~ NA_real_, 40 < pn_corr & pn_corr < 50 ~ 1, 22 > pn_corr | pn_corr > 68 ~ 0, 22 < pn_corr & pn_corr < 40 ~ round(0.056 * (pn_corr - 22), 3), TRUE ~ round(1 - 0.056 * (pn_corr - 50), 3) ), # mortalité mort_n = round(1.0 * exp(-0.031 * mort), 3) ) # ========================= # 4) Note carrière (pondérée) # ========================= weights <- as.numeric(Pcar[1, 2:16]) # 15 pondérations norm_cols <- c( "agevel1_n","ivv1_n","ivv2p_n","pad_n","ptgv_n","prol_n","pn_n", "txvf_n","txm_n","mort_n","txrepros_n","nbpp_n","ptgp_n","p120_n","p210_n" ) v_final <- v_norm %>% rowwise() %>% mutate( SOMME_tot = { x <- c_across(all_of(norm_cols)) w <- weights mask <- !is.na(x) & !is.na(w) if (sum(mask) == 0) NA_real_ else (sum(x[mask] * w[mask]) / sum(w[mask])) * 10 }, ecowcarr = if_else( is.na(ptgp_n) & is.na(p120_n) & is.na(p210_n), # règle d'exclusion VA4 NA_real_, round(SOMME_tot * 100, 0) ) ) %>% ungroup() `` ################################################################################################################### ############################## FIN DE REFACTO ####################################################################