210 lines
7.8 KiB
R
210 lines
7.8 KiB
R
###################################################################################################################
|
|
################# 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 #################################################################### |