library(multcomp) #Complementar ao emmeans, fornece valor de p e letras para os testes de medias.
## Carregando pacotes exigidos: mvtnorm
## Carregando pacotes exigidos: survival
## Carregando pacotes exigidos: TH.data
## Carregando pacotes exigidos: MASS
##
## Attaching package: 'TH.data'
## The following object is masked from 'package:MASS':
##
## geyser
library(emmeans) #Faz anova para modelos mistos.
library(lmerTest) #complementar ao emmeans para funcionar a função lmer
## Carregando pacotes exigidos: lme4
## Carregando pacotes exigidos: Matrix
##
## Attaching package: 'lmerTest'
## The following object is masked from 'package:lme4':
##
## lmer
## The following object is masked from 'package:stats':
##
## step
library(lme4)
library(multcompView)
library(dplyr)
##
## Attaching package: 'dplyr'
## The following object is masked from 'package:MASS':
##
## select
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(tidyr)
##
## Attaching package: 'tidyr'
## The following objects are masked from 'package:Matrix':
##
## expand, pack, unpack
library(readxl)
data = read_excel("C:/Users/Samsung/OneDrive/Documentos/2 - Material de trabalho/7 - Projetos/7.1 Pesquisa/Leguminosas IZ/Plan_PG_Luciana_VO210826.xlsx", sheet = "R table")
head(data)
#fatores
data$Planta=as.factor(data$Planta)
data$Inclusao=as.numeric(data$Inclusao)
data$DMO=as.numeric(data$DMO)
## Warning: NAs introduzidos por coerção
str(data)
## tibble [84 × 11] (S3: tbl_df/tbl/data.frame)
## $ Planta : Factor w/ 8 levels "Calopogonio",..: 8 8 8 5 5 5 3 3 3 1 ...
## $ Inclusao : num [1:84] 100 100 100 0 0 0 100 100 100 100 ...
## $ tempo : chr [1:84] "24h" "24h" "24h" "24h" ...
## $ inoculo : num [1:84] 1 1 1 1 1 1 1 1 1 1 ...
## $ PG_Total : num [1:84] 48.7 55.5 46.5 44.3 92.3 ...
## $ DMO : num [1:84] NA NA NA 242 250 ...
## $ Acetico : num [1:84] 50.3 78.8 22.6 45.7 42.1 ...
## $ Propionico: num [1:84] 5.82 4.82 4.94 24.47 20.42 ...
## $ Butirico : num [1:84] 10.52 11.29 8.67 20.72 17.67 ...
## $ CH4.24h : num [1:84] 2.51 1.89 1.72 5.68 5.26 ...
## $ AP.ratio : num [1:84] 8.65 16.37 4.58 1.87 2.06 ...
View(data)
# Filtrar apenas CM, M e M_CM
data_M = data %>% dplyr::filter(Planta %in% c("Marandu", "Macrotiloma", "Macrotiloma_Marandu"))
table(data_M$Planta)
##
## Calopogonio Calopogonio_Marandu Macrotiloma Macrotiloma_Marandu
## 0 0 6 18
## Marandu Soja_Marandu Soja_perene Tifton_LANA
## 6 0 0 0
View(data_M)
###################################
data_M$Inclusao_sc = (data_M$Inclusao - 50) / 25
variables = c("PG_Total","DMO","Acetico","Propionico","Butirico","CH4.24h","AP.ratio")
for (var in variables) {
cat("\n\n##############################\n")
cat("Variável:", var, "\n")
cat("################################\n")
# Modelo linear
formula_lin = as.formula(paste(var, "~ Inclusao_sc + (1|inoculo)"))
mod_lin = lmer(formula_lin, data = data_M, REML = FALSE)
# Modelo quadrático
formula_quad = as.formula(paste(var, "~ Inclusao_sc + I(Inclusao_sc^2) + (1|inoculo)"))
mod_quad = lmer(formula_quad,data = data_M,REML = FALSE)
cat("\n---------------- MODELO LINEAR ----------------\n")
print(anova(mod_lin))
cat("\n---------------- MODELO QUADRÁTICO ----------------\n")
print(anova(mod_quad))
cat("\n---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------\n")
print(anova(mod_lin, mod_quad))}
##
##
## ##############################
## Variável: PG_Total
## ################################
##
## ---------------- MODELO LINEAR ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 1641.8 1641.8 1 28 7.5012 0.0106 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- MODELO QUADRÁTICO ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 1641.8 1641.8 1 28 12.811 0.0012821 **
## I(Inclusao_sc^2) 2540.0 2540.0 1 28 19.820 0.0001238 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------
## Data: data_M
## Models:
## mod_lin: formula_lin
## mod_quad: formula_quad
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## mod_lin 4 261.24 266.84 -126.62 253.24
## mod_quad 5 248.25 255.26 -119.12 238.25 14.987 1 0.0001083 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## ##############################
## Variável: DMO
## ################################
##
## ---------------- MODELO LINEAR ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 38950 38950 1 19.017 9.0537 0.007212 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- MODELO QUADRÁTICO ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 23100 23100 1 19.016 8.5049 0.008850 **
## I(Inclusao_sc^2) 30600 30600 1 19.009 11.2662 0.003313 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------
## Data: data_M
## Models:
## mod_lin: formula_lin
## mod_quad: formula_quad
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## mod_lin 4 248.24 252.42 -120.12 240.24
## mod_quad 5 241.34 246.56 -115.67 231.34 8.903 1 0.002847 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## ##############################
## Variável: Acetico
## ################################
##
## ---------------- MODELO LINEAR ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 438.97 438.97 1 28 2.726 0.1099
##
## ---------------- MODELO QUADRÁTICO ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 438.97 438.97 1 28 3.2880 0.08052 .
## I(Inclusao_sc^2) 770.67 770.67 1 28 5.7725 0.02315 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------
## Data: data_M
## Models:
## mod_lin: formula_lin
## mod_quad: formula_quad
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## mod_lin 4 252.23 257.83 -122.11 244.23
## mod_quad 5 248.98 255.98 -119.49 238.98 5.2484 1 0.02197 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## ##############################
## Variável: Propionico
## ################################
##
## ---------------- MODELO LINEAR ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 0.1371 0.1371 1 28 0.0021 0.9634
##
## ---------------- MODELO QUADRÁTICO ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 0.14 0.14 1 28 0.0028 0.958437
## I(Inclusao_sc^2) 399.34 399.34 1 28 8.0538 0.008354 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------
## Data: data_M
## Models:
## mod_lin: formula_lin
## mod_quad: formula_quad
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## mod_lin 4 226.03 231.64 -109.02 218.03
## mod_quad 5 220.95 227.96 -105.48 210.95 7.0786 1 0.007801 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## ##############################
## Variável: Butirico
## ################################
##
## ---------------- MODELO LINEAR ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 80.002 80.002 1 28 3.2116 0.08393 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- MODELO QUADRÁTICO ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 80.002 80.002 1 28 3.2304 0.08308 .
## I(Inclusao_sc^2) 4.051 4.051 1 28 0.1636 0.68896
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------
## Data: data_M
## Models:
## mod_lin: formula_lin
## mod_quad: formula_quad
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## mod_lin 4 197.78 203.38 -94.890 189.78
## mod_quad 5 199.62 206.62 -94.808 189.62 0.1631 1 0.6863
##
##
## ##############################
## Variável: CH4.24h
## ################################
## boundary (singular) fit: see help('isSingular')
##
## ---------------- MODELO LINEAR ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 0.39718 0.39718 1 30 0.0973 0.7572
##
## ---------------- MODELO QUADRÁTICO ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 0.397 0.397 1 28 0.1692 0.6839
## I(Inclusao_sc^2) 49.691 49.691 1 28 21.1744 8.238e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------
## Data: data_M
## Models:
## mod_lin: formula_lin
## mod_quad: formula_quad
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## mod_lin 4 135.32 140.93 -63.662 127.32
## mod_quad 5 121.53 128.54 -55.766 111.53 15.792 1 7.071e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## ##############################
## Variável: AP.ratio
## ################################
##
## ---------------- MODELO LINEAR ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 0.3687 0.3687 1 28 12.412 0.001485 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- MODELO QUADRÁTICO ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 0.36870 0.36870 1 28 12.4884 0.001443 **
## I(Inclusao_sc^2) 0.00512 0.00512 1 28 0.1734 0.680249
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------
## Data: data_M
## Models:
## mod_lin: formula_lin
## mod_quad: formula_quad
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## mod_lin 4 -5.7998 -0.1950 6.8999 -13.800
## mod_quad 5 -3.9727 3.0333 6.9863 -13.973 0.1729 1 0.6775
library(lme4)
library(ggplot2)
variables = c("PG_Total","DMO","CH4.24h","AP.ratio")
for (var in variables) {
# Modelo quadrático
formula_quad = as.formula(paste(var, "~ Inclusao_sc + I(Inclusao_sc^2) + (1|inoculo)"))
mod_quad = lmer(formula_quad,data = data_M,REML = FALSE)
# Dados para construir a curva
pred = data.frame(Inclusao = seq(
min(data_M$Inclusao, na.rm = TRUE),
max(data_M$Inclusao, na.rm = TRUE),length.out = 200))
# Mesma transformação utilizada no modelo
pred$Inclusao_sc = (pred$Inclusao - 50) / 25
# Predição apenas dos efeitos fixos
pred$y_pred = predict(mod_quad,newdata = pred,re.form = NA)
# Plot
p = ggplot(data_M,aes(x = Inclusao,y = .data[[var]])) +
geom_point(size = 2.5,alpha = 0.7) +
geom_line(data = pred,aes(x = Inclusao,y = y_pred),linewidth = 1) +
scale_x_continuous(breaks = c(0, 25, 50, 75, 100)) +
labs(x = "Inclusão (%)",y = var) +
theme_classic(base_size = 14)
print(p)}
## Warning: Removed 9 rows containing missing values or values outside the scale range
## (`geom_point()`).
cada variavel teve pesos iguais, as médias foram ajustadas via modelo simples Direções: - DMO: quanto maior melhor; - PG_Total, CH4, AP.ratio: quanto menor melhor (invertido);
EXPLICAÇÃO PARA TEXTOS CIENTIFICOS: Foi construído um índice multicritério para ranquear as inclusões de leguminosas com base na degradabilidade da matéria orgânica (DMO), na produção de gas (PG_Total), produção de metano (CH4.24h) e razão acetato:propionato (AP.ratio). As médias ajustadas foram obtidas por modelo linear considerando inclusão como efeito fixo, seguidas de estimativas de médias marginais (emmeans).Para permitir a comparação entre variáveis em escalas distintas, os valores foram padronizados por transformação min–max (0–1). Para DMO, valores maiores foram considerados mais favoráveis. Para PG_Total, CH4.24h e AP.ratio, a escala foi invertida, de modo que menores produções resultassem em maiores escores padronizados.O índice final foi calculado como a média aritmética simples dos escores padronizados (pesos iguais), e os tratamentos foram ranqueados em ordem decrescente do escore composto.
data_M_ranking = data %>% dplyr::filter(Planta %in% c("Macrotiloma_Marandu"))
View(data_M_ranking)
########################### FUNÇÃO MODELO ###########################
fun_emm = function(data_M_ranking, response){
f <- as.formula(paste0(response, " ~ Inclusao + factor(inoculo)"))
mod <- lm(f, data = data_M_ranking)
emmeans(mod,~ Inclusao,at = list(Inclusao = sort(unique(data_M_ranking$Inclusao)))) %>%
as.data.frame() %>%
mutate(variavel = response)}
########################### MÉDIAS AJUSTADAS (EMMEANS) ###########################
variaveis_modelo = c("PG_Total","DMO","CH4.24h","AP.ratio")
medias1 = bind_rows(lapply(variaveis_modelo, \(v) fun_emm(data_M_ranking, v)))
View(medias1)
########################### ORGANIZAR EM FORMATO LARGO ###########################
emm_wide = medias1 %>% select(Inclusao, variavel, emmean, SE) %>%
pivot_wider(names_from = variavel,values_from = c(emmean, SE))
emm_wide
########################### PADRONIZAÇÃO EM ESCALA ###########################
ranking = emm_wide %>% ungroup() %>% mutate(
# DMO (benefício: maior é melhor)
DMO_score = (emmean_DMO - min(emmean_DMO, na.rm = TRUE)) /
(max(emmean_DMO, na.rm = TRUE) - min(emmean_DMO, na.rm = TRUE)),
# PG_Total (invertido: menor é melhor)
PG_Total_score = (max(emmean_PG_Total, na.rm = TRUE) - emmean_PG_Total) /
(max(emmean_PG_Total, na.rm = TRUE) - min(emmean_PG_Total, na.rm = TRUE)),
# CH4.24h (invertido: menor é melhor)
CH4.24h_score = (max(emmean_CH4.24h, na.rm = TRUE) - emmean_CH4.24h) /
(max(emmean_CH4.24h, na.rm = TRUE) - min(emmean_CH4.24h, na.rm = TRUE)),
# AP.ratio (invertido: menor é melhor)
AP.ratio_score = (max(emmean_AP.ratio, na.rm = TRUE) - emmean_AP.ratio) /
(max(emmean_AP.ratio, na.rm = TRUE) - min(emmean_AP.ratio, na.rm = TRUE)))
ranking
########################### ÍNDICE MULTICRITÉRIO (PESOS IGUAIS) ###########################
ranking = ranking %>%
mutate(Score = (DMO_score + PG_Total_score + CH4.24h_score + AP.ratio_score)/4) %>%
arrange(desc(Score)) %>% mutate(Rank = row_number())
ranking_final = ranking %>% select(Inclusao,
emmean_DMO, DMO_score,
emmean_PG_Total, PG_Total_score,
emmean_CH4.24h, CH4.24h_score,
emmean_AP.ratio, AP.ratio_score,
Score, Rank)
ranking_final
CONCLUSÃO
Para macrotiloma, a melhor inclusão é 25%, com maior DMO e menor PG e CH4.
# Filtrar apenas CM, M e M_CM
data_C = data %>% dplyr::filter(Planta %in% c("Marandu", "Calopogonio", "Calopogonio_Marandu"))
table(data_C$Planta)
##
## Calopogonio Calopogonio_Marandu Macrotiloma Macrotiloma_Marandu
## 6 18 0 0
## Marandu Soja_Marandu Soja_perene Tifton_LANA
## 6 0 0 0
View(data_C)
###################################
data_C$Inclusao_sc = (data_C$Inclusao - 50) / 25
variables = c("PG_Total","DMO","Acetico","Propionico","Butirico","CH4.24h","AP.ratio")
for (var in variables) {
cat("\n\n###########################\n")
cat("Variável:", var, "\n")
cat("##############################\n")
# Modelo linear
formula_lin = as.formula(paste(var, "~ Inclusao_sc + (1|inoculo)"))
mod_lin = lmer(formula_lin, data = data_C, REML = FALSE)
# Modelo quadrático
formula_quad = as.formula(paste(var, "~ Inclusao_sc + I(Inclusao_sc^2) + (1|inoculo)"))
mod_quad = lmer(formula_quad,data = data_C,REML = FALSE)
cat("\n---------------- MODELO LINEAR ----------------\n")
print(anova(mod_lin))
cat("\n---------------- MODELO QUADRÁTICO ----------------\n")
print(anova(mod_quad))
cat("\n---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------\n")
print(anova(mod_lin, mod_quad))}
##
##
## ###########################
## Variável: PG_Total
## ##############################
##
## ---------------- MODELO LINEAR ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 1458.8 1458.8 1 28 5.1057 0.03182 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- MODELO QUADRÁTICO ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 1458.8 1458.8 1 28 5.9789 0.02103 *
## I(Inclusao_sc^2) 1168.6 1168.6 1 28 4.7892 0.03714 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------
## Data: data_C
## Models:
## mod_lin: formula_lin
## mod_quad: formula_quad
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## mod_lin 4 267.99 273.59 -130.00 259.99
## mod_quad 5 265.57 272.57 -127.78 255.57 4.4211 1 0.0355 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## ###########################
## Variável: DMO
## ##############################
##
## ---------------- MODELO LINEAR ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 84207 84207 1 19 12.367 0.002307 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- MODELO QUADRÁTICO ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 76213 76213 1 18.998 11.8834 0.0027 **
## I(Inclusao_sc^2) 7928 7928 1 19.014 1.2361 0.2801
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------
## Data: data_C
## Models:
## mod_lin: formula_lin
## mod_quad: formula_quad
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## mod_lin 4 257.6 261.78 -124.8 249.6
## mod_quad 5 258.4 263.62 -124.2 248.4 1.2004 1 0.2732
##
##
## ###########################
## Variável: Acetico
## ##############################
##
## ---------------- MODELO LINEAR ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 57.753 57.753 1 28 0.4737 0.497
##
## ---------------- MODELO QUADRÁTICO ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 57.753 57.753 1 28 0.4737 0.4969
## I(Inclusao_sc^2) 0.124 0.124 1 28 0.0010 0.9748
##
## ---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------
## Data: data_C
## Models:
## mod_lin: formula_lin
## mod_quad: formula_quad
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## mod_lin 4 244.63 250.24 -118.32 236.63
## mod_quad 5 246.63 253.64 -118.32 236.63 0.001 1 0.9746
##
##
## ###########################
## Variável: Propionico
## ##############################
##
## ---------------- MODELO LINEAR ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 81.571 81.571 1 28 1.0696 0.3099
##
## ---------------- MODELO QUADRÁTICO ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 81.571 81.571 1 28 1.0968 0.3039
## I(Inclusao_sc^2) 52.924 52.924 1 28 0.7116 0.4061
##
## ---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------
## Data: data_C
## Models:
## mod_lin: formula_lin
## mod_quad: formula_quad
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## mod_lin 4 230.85 236.46 -111.42 222.85
## mod_quad 5 232.15 239.15 -111.07 222.15 0.7027 1 0.4019
##
##
## ###########################
## Variável: Butirico
## ##############################
##
## ---------------- MODELO LINEAR ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 29.348 29.348 1 28 1.1839 0.2858
##
## ---------------- MODELO QUADRÁTICO ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 29.3477 29.3477 1 28 1.1854 0.2855
## I(Inclusao_sc^2) 0.8858 0.8858 1 28 0.0358 0.8513
##
## ---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------
## Data: data_C
## Models:
## mod_lin: formula_lin
## mod_quad: formula_quad
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## mod_lin 4 197.96 203.56 -94.980 189.96
## mod_quad 5 199.92 206.93 -94.962 189.92 0.0358 1 0.85
##
##
## ###########################
## Variável: CH4.24h
## ##############################
##
## ---------------- MODELO LINEAR ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 2.3857 2.3857 1 28 0.5707 0.4563
##
## ---------------- MODELO QUADRÁTICO ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 2.3857 2.3857 1 28 0.6815 0.41604
## I(Inclusao_sc^2) 19.0289 19.0289 1 28 5.4360 0.02715 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------
## Data: data_C
## Models:
## mod_lin: formula_lin
## mod_quad: formula_quad
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## mod_lin 4 139.15 144.76 -65.577 131.15
## mod_quad 5 136.19 143.19 -63.093 126.19 4.968 1 0.02582 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## ###########################
## Variável: AP.ratio
## ##############################
##
## ---------------- MODELO LINEAR ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 0.012339 0.012339 1 28 0.5355 0.4704
##
## ---------------- MODELO QUADRÁTICO ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 0.012339 0.012339 1 28 0.6571 0.42442
## I(Inclusao_sc^2) 0.119383 0.119383 1 28 6.3576 0.01766 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------
## Data: data_C
## Models:
## mod_lin: formula_lin
## mod_quad: formula_quad
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## mod_lin 4 -13.843 -8.2381 10.921 -21.843
## mod_quad 5 -17.572 -10.5662 13.786 -27.572 5.7293 1 0.01668 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
variables = c("PG_Total","DMO","CH4.24h","AP.ratio")
for (var in variables) {
# Modelo quadrático
formula_quad = as.formula(paste(var, "~ Inclusao_sc + I(Inclusao_sc^2) + (1|inoculo)"))
mod_quad = lmer(formula_quad,data = data_C,REML = FALSE)
# Dados para construir a curva
pred = data.frame(Inclusao = seq(
min(data_C$Inclusao, na.rm = TRUE),
max(data_C$Inclusao, na.rm = TRUE),length.out = 200))
# Mesma transformação utilizada no modelo
pred$Inclusao_sc = (pred$Inclusao - 50) / 25
# Predição apenas dos efeitos fixos
pred$y_pred = predict(mod_quad,newdata = pred,re.form = NA)
# Plot
p = ggplot(data_C,aes(x = Inclusao,y = .data[[var]])) +
geom_point(size = 2.5,alpha = 0.7) +
geom_line(data = pred,aes(x = Inclusao,y = y_pred),linewidth = 1) +
scale_x_continuous(breaks = c(0, 25, 50, 75, 100)) +
labs(x = "Inclusão (%)",y = var) +
theme_classic(base_size = 14)
print(p)}
## Warning: Removed 9 rows containing missing values or values outside the scale range
## (`geom_point()`).
cada variavel teve pesos iguais, as médias foram ajustadas via modelo simples Direções: - DMO: quanto maior melhor; - PG_Total, CH4, AP.ratio: quanto menor melhor (invertido);
EXPLICAÇÃO PARA TEXTOS CIENTIFICOS: Foi construído um índice multicritério para ranquear as inclusões de leguminosas com base na degradabilidade da matéria orgânica (DMO), na produção de gas (PG_Total), produção de metano (CH4.24h) e razão acetato:propionato (AP.ratio). As médias ajustadas foram obtidas por modelo linear considerando inclusão como efeito fixo, seguidas de estimativas de médias marginais (emmeans).Para permitir a comparação entre variáveis em escalas distintas, os valores foram padronizados por transformação min–max (0–1). Para DMO, valores maiores foram considerados mais favoráveis. Para PG_Total, CH4.24h e AP.ratio, a escala foi invertida, de modo que menores produções resultassem em maiores escores padronizados.O índice final foi calculado como a média aritmética simples dos escores padronizados (pesos iguais), e os tratamentos foram ranqueados em ordem decrescente do escore composto.
data_C_ranking = data %>% dplyr::filter(Planta %in% c("Calopogonio_Marandu"))
View(data_C_ranking)
########################### FUNÇÃO MODELO ###########################
fun_emm = function(data_C_ranking, response){
f <- as.formula(paste0(response, " ~ Inclusao + factor(inoculo)"))
mod <- lm(f, data = data_C_ranking)
emmeans(mod,~ Inclusao,at = list(Inclusao = sort(unique(data_C_ranking$Inclusao)))) %>%
as.data.frame() %>%
mutate(variavel = response)}
########################### MÉDIAS AJUSTADAS (EMMEANS) ###########################
variaveis_modelo = c("PG_Total","DMO","CH4.24h","AP.ratio")
medias1 = bind_rows(lapply(variaveis_modelo, \(v) fun_emm(data_C_ranking, v)))
View(medias1)
########################### ORGANIZAR EM FORMATO LARGO ###########################
emm_wide = medias1 %>% select(Inclusao, variavel, emmean, SE) %>%
pivot_wider(names_from = variavel,values_from = c(emmean, SE))
emm_wide
########################### PADRONIZAÇÃO EM ESCALA ###########################
ranking = emm_wide %>% ungroup() %>% mutate(
# DMO (benefício: maior é melhor)
DMO_score = (emmean_DMO - min(emmean_DMO, na.rm = TRUE)) /
(max(emmean_DMO, na.rm = TRUE) - min(emmean_DMO, na.rm = TRUE)),
# PG_Total (invertido: menor é melhor)
PG_Total_score = (max(emmean_PG_Total, na.rm = TRUE) - emmean_PG_Total) /
(max(emmean_PG_Total, na.rm = TRUE) - min(emmean_PG_Total, na.rm = TRUE)),
# CH4.24h (invertido: menor é melhor)
CH4.24h_score = (max(emmean_CH4.24h, na.rm = TRUE) - emmean_CH4.24h) /
(max(emmean_CH4.24h, na.rm = TRUE) - min(emmean_CH4.24h, na.rm = TRUE)),
# AP.ratio (invertido: menor é melhor)
AP.ratio_score = (max(emmean_AP.ratio, na.rm = TRUE) - emmean_AP.ratio) /
(max(emmean_AP.ratio, na.rm = TRUE) - min(emmean_AP.ratio, na.rm = TRUE)))
ranking
########################### ÍNDICE MULTICRITÉRIO (PESOS IGUAIS) ###########################
ranking = ranking %>%
mutate(Score = (DMO_score + PG_Total_score + CH4.24h_score + AP.ratio_score)/4) %>%
arrange(desc(Score)) %>% mutate(Rank = row_number())
ranking_final = ranking %>% select(Inclusao,
emmean_DMO, DMO_score,
emmean_PG_Total, PG_Total_score,
emmean_CH4.24h, CH4.24h_score,
emmean_AP.ratio, AP.ratio_score,
Score, Rank)
ranking_final
CONCLUSÃO
Para Calopognio, não foi possivel encontrar uma diferença de score, observamos diferenças no score individual de cada variavel, mas o score médio (as 4 variaveis dividido por 4) ficou 0.5 para os tres niveis. Então a melhor inclusão terá que ser escolhida por um fator de maior peso, eu escolheria 25% de inclusão, por ser o nivel com maior DMO e menor A:P ratio.
# Filtrar apenas CM, M e M_CM
data_SP = data %>% dplyr::filter(Planta %in% c("Marandu", "Soja_perene", "Soja_Marandu"))
table(data_SP$Planta)
##
## Calopogonio Calopogonio_Marandu Macrotiloma Macrotiloma_Marandu
## 0 0 0 0
## Marandu Soja_Marandu Soja_perene Tifton_LANA
## 6 18 6 0
View(data_SP)
###################################
data_SP$Inclusao_sc = (data_SP$Inclusao - 50) / 25
variables = c("PG_Total","DMO","Acetico","Propionico","Butirico","CH4.24h","AP.ratio")
for (var in variables) {
cat("\n\n###########################\n")
cat("Variável:", var, "\n")
cat("#########################\n")
# Modelo linear
formula_lin = as.formula(paste(var, "~ Inclusao_sc + (1|inoculo)"))
mod_lin = lmer(formula_lin, data = data_SP, REML = FALSE)
# Modelo quadrático
formula_quad = as.formula(paste(var, "~ Inclusao_sc + I(Inclusao_sc^2) + (1|inoculo)"))
mod_quad = lmer(formula_quad,data = data_SP,REML = FALSE)
cat("\n---------------- MODELO LINEAR ----------------\n")
print(anova(mod_lin))
cat("\n---------------- MODELO QUADRÁTICO ----------------\n")
print(anova(mod_quad))
cat("\n---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------\n")
print(anova(mod_lin, mod_quad))}
##
##
## ###########################
## Variável: PG_Total
## #########################
##
## ---------------- MODELO LINEAR ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 969.48 969.48 1 28 5.1001 0.03191 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- MODELO QUADRÁTICO ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 969.48 969.48 1 28 7.6451 0.0099592 **
## I(Inclusao_sc^2) 1771.86 1771.86 1 28 13.9724 0.0008445 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------
## Data: data_SP
## Models:
## mod_lin: formula_lin
## mod_quad: formula_quad
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## mod_lin 4 257.36 262.96 -124.68 249.36
## mod_quad 5 248.02 255.03 -119.01 238.02 11.335 1 0.0007608 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## ###########################
## Variável: DMO
## #########################
##
## ---------------- MODELO LINEAR ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 80111 80111 1 20 31.754 1.626e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- MODELO QUADRÁTICO ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 67732 67732 1 20 49.764 7.686e-07 ***
## I(Inclusao_sc^2) 23235 23235 1 20 17.071 0.000517 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------
## Data: data_SP
## Models:
## mod_lin: formula_lin
## mod_quad: formula_quad
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## mod_lin 4 247.62 251.98 -119.81 239.62
## mod_quad 5 237.27 242.73 -113.64 227.27 12.342 1 0.0004428 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## ###########################
## Variável: Acetico
## #########################
##
## ---------------- MODELO LINEAR ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 115.66 115.66 1 28 0.3935 0.5355
##
## ---------------- MODELO QUADRÁTICO ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 115.66 115.66 1 28 0.5307 0.472374
## I(Inclusao_sc^2) 2126.42 2126.42 1 28 9.7566 0.004129 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------
## Data: data_SP
## Models:
## mod_lin: formula_lin
## mod_quad: formula_quad
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## mod_lin 4 267.35 272.96 -129.68 259.35
## mod_quad 5 260.98 267.99 -125.49 250.98 8.3707 1 0.003813 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## ###########################
## Variável: Propionico
## #########################
##
## ---------------- MODELO LINEAR ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 0.032648 0.032648 1 28 3e-04 0.9852
##
## ---------------- MODELO QUADRÁTICO ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 0.03 0.03 1 28 0.0006 0.9809934
## I(Inclusao_sc^2) 1034.31 1034.31 1 28 18.3050 0.0001984 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------
## Data: data_SP
## Models:
## mod_lin: formula_lin
## mod_quad: formula_quad
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## mod_lin 4 235.44 241.04 -113.72 227.44
## mod_quad 5 223.35 230.36 -106.67 213.35 14.085 1 0.0001747 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## ###########################
## Variável: Butirico
## #########################
##
## ---------------- MODELO LINEAR ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 99.311 99.311 1 28 2.9097 0.09912 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- MODELO QUADRÁTICO ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 99.311 99.311 1 28 3.3651 0.07724 .
## I(Inclusao_sc^2) 129.324 129.324 1 28 4.3821 0.04550 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------
## Data: data_SP
## Models:
## mod_lin: formula_lin
## mod_quad: formula_quad
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## mod_lin 4 205.82 211.42 -98.908 197.82
## mod_quad 5 203.75 210.75 -96.873 193.75 4.0712 1 0.04362 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## ###########################
## Variável: CH4.24h
## #########################
##
## ---------------- MODELO LINEAR ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 0.000663 0.000663 1 28 1e-04 0.9906
##
## ---------------- MODELO QUADRÁTICO ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 0.0007 0.0007 1 28 0.0002 0.9902
## I(Inclusao_sc^2) 8.8480 8.8480 1 28 2.0379 0.1645
##
## ---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------
## Data: data_SP
## Models:
## mod_lin: formula_lin
## mod_quad: formula_quad
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## mod_lin 4 142.43 148.04 -67.216 134.43
## mod_quad 5 142.47 149.47 -66.233 132.47 1.9671 1 0.1608
##
##
## ###########################
## Variável: AP.ratio
## #########################
##
## ---------------- MODELO LINEAR ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 0.03929 0.03929 1 28 0.4863 0.4913
##
## ---------------- MODELO QUADRÁTICO ----------------
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Inclusao_sc 0.039290 0.039290 1 28 0.4951 0.4875
## I(Inclusao_sc^2) 0.040199 0.040199 1 28 0.5065 0.4825
##
## ---------------- COMPARAÇÃO LINEAR vs QUADRÁTICO ----------------
## Data: data_SP
## Models:
## mod_lin: formula_lin
## mod_quad: formula_quad
## npar AIC BIC logLik deviance Chisq Df Pr(>Chisq)
## mod_lin 4 22.204 27.808 -7.1018 14.204
## mod_quad 5 23.702 30.708 -6.8508 13.702 0.502 1 0.4786
variables = c("PG_Total","DMO","CH4.24h","AP.ratio")
for (var in variables) {
# Modelo quadrático
formula_quad = as.formula(paste(var, "~ Inclusao_sc + I(Inclusao_sc^2) + (1|inoculo)"))
mod_quad = lmer(formula_quad,data = data_SP,REML = FALSE)
# Dados para construir a curva
pred = data.frame(Inclusao = seq(
min(data_SP$Inclusao, na.rm = TRUE),
max(data_SP$Inclusao, na.rm = TRUE),length.out = 200))
# Mesma transformação utilizada no modelo
pred$Inclusao_sc = (pred$Inclusao - 50) / 25
# Predição apenas dos efeitos fixos
pred$y_pred = predict(mod_quad,newdata = pred,re.form = NA)
# Plot
p = ggplot(data_SP,aes(x = Inclusao,y = .data[[var]])) +
geom_point(size = 2.5,alpha = 0.7) +
geom_line(data = pred,aes(x = Inclusao,y = y_pred),linewidth = 1) +
scale_x_continuous(breaks = c(0, 25, 50, 75, 100)) +
labs(x = "Inclusão (%)",y = var) +
theme_classic(base_size = 14)
print(p)}
## Warning: Removed 8 rows containing missing values or values outside the scale range
## (`geom_point()`).
cada variavel teve pesos iguais, as médias foram ajustadas via modelo simples Direções: - DMO: quanto maior melhor; - PG_Total, CH4, AP.ratio: quanto menor melhor (invertido);
EXPLICAÇÃO PARA TEXTOS CIENTIFICOS: Foi construído um índice multicritério para ranquear as inclusões de leguminosas com base na degradabilidade da matéria orgânica (DMO), na produção de gas (PG_Total), produção de metano (CH4.24h) e razão acetato:propionato (AP.ratio). As médias ajustadas foram obtidas por modelo linear considerando inclusão como efeito fixo, seguidas de estimativas de médias marginais (emmeans).Para permitir a comparação entre variáveis em escalas distintas, os valores foram padronizados por transformação min–max (0–1). Para DMO, valores maiores foram considerados mais favoráveis. Para PG_Total, CH4.24h e AP.ratio, a escala foi invertida, de modo que menores produções resultassem em maiores escores padronizados.O índice final foi calculado como a média aritmética simples dos escores padronizados (pesos iguais), e os tratamentos foram ranqueados em ordem decrescente do escore composto.
data_SP_ranking = data %>% dplyr::filter(Planta %in% c("Soja_Marandu"))
View(data_SP_ranking)
########################### FUNÇÃO MODELO ###########################
fun_emm = function(data_SP_ranking, response){
f <- as.formula(paste0(response, " ~ Inclusao + factor(inoculo)"))
mod <- lm(f, data = data_SP_ranking)
emmeans(mod,~ Inclusao,at = list(Inclusao = sort(unique(data_SP_ranking$Inclusao)))) %>%
as.data.frame() %>%
mutate(variavel = response)}
########################### MÉDIAS AJUSTADAS (EMMEANS) ###########################
variaveis_modelo = c("PG_Total","DMO","CH4.24h","AP.ratio")
medias1 = bind_rows(lapply(variaveis_modelo, \(v) fun_emm(data_SP_ranking, v)))
View(medias1)
########################### ORGANIZAR EM FORMATO LARGO ###########################
emm_wide = medias1 %>% select(Inclusao, variavel, emmean, SE) %>%
pivot_wider(names_from = variavel,values_from = c(emmean, SE))
emm_wide
########################### PADRONIZAÇÃO EM ESCALA ###########################
ranking = emm_wide %>% ungroup() %>% mutate(
# DMO (benefício: maior é melhor)
DMO_score = (emmean_DMO - min(emmean_DMO, na.rm = TRUE)) /
(max(emmean_DMO, na.rm = TRUE) - min(emmean_DMO, na.rm = TRUE)),
# PG_Total (invertido: menor é melhor)
PG_Total_score = (max(emmean_PG_Total, na.rm = TRUE) - emmean_PG_Total) /
(max(emmean_PG_Total, na.rm = TRUE) - min(emmean_PG_Total, na.rm = TRUE)),
# CH4.24h (invertido: menor é melhor)
CH4.24h_score = (max(emmean_CH4.24h, na.rm = TRUE) - emmean_CH4.24h) /
(max(emmean_CH4.24h, na.rm = TRUE) - min(emmean_CH4.24h, na.rm = TRUE)),
# AP.ratio (invertido: menor é melhor)
AP.ratio_score = (max(emmean_AP.ratio, na.rm = TRUE) - emmean_AP.ratio) /
(max(emmean_AP.ratio, na.rm = TRUE) - min(emmean_AP.ratio, na.rm = TRUE)))
ranking
########################### ÍNDICE MULTICRITÉRIO (PESOS IGUAIS) ###########################
ranking = ranking %>%
mutate(Score = (DMO_score + PG_Total_score + CH4.24h_score + AP.ratio_score)/4) %>%
arrange(desc(Score)) %>% mutate(Rank = row_number())
ranking_final = ranking %>% select(Inclusao,
emmean_DMO, DMO_score,
emmean_PG_Total, PG_Total_score,
emmean_CH4.24h, CH4.24h_score,
emmean_AP.ratio, AP.ratio_score,
Score, Rank)
ranking_final
CONCLUSÃO
Para sojá perene, o melhor nivel foi 75%, com maior DMO e menor CH4 e A:P ratio.