Caracterizacion papa
Pedro Lizarazo, Geraldyne Rozo, Jimer Lozano
2025
Librerias usadas
library(readxl)
library(rlang)
library(dplyr)
library(tidyverse)
library(stringi)
library(corrplot)
library(GGally)
library(textshape)
library(FactoMineR)
library(factoextra)
library(rlang)
library(tibble)
library(ggstatsplot)
library(ggsci)
library(agricolae)1. M A N I P U L A C I O N D E L A S B A S E S D E D A T O S
Cargar municipios con aptitud para el cultivo de papa
AP.TOT <- read_excel("MUNICIPIOS CON ALTO DESMPEÑO PRODUCTIVO PAPA 07.10.23.xlsx")
# Quitar las tildes a la columna de municipios
AP.TOT$MUN <- toupper(stri_trans_general(AP.TOT$MUN,"Latin-ASCII"))
AP.CUN <- AP.TOT %>% filter(., DEPARTAMENTO=="CUNDINAMARCA") %>%
mutate(., IDPM= pmin(IDPMi, IDPMii, na.rm = T) )
head(AP.CUN, n = 10)## # A tibble: 10 × 7
## DEPARTAMENTO MUN IDPMi IDPMii ALTITUD TEMPERATURA IDPM
## <chr> <chr> <dbl> <dbl> <chr> <chr> <dbl>
## 1 CUNDINAMARCA BOGOTA 3 3 2625 13.1 3
## 2 CUNDINAMARCA CAJICA 3 3 2558 14 3
## 3 CUNDINAMARCA CARMEN DE CARUPA 1 2 2600 12 1
## 4 CUNDINAMARCA CHIA 2 2 2600 14 2
## 5 CUNDINAMARCA CHIPAQUE 3 3 2400 13 3
## 6 CUNDINAMARCA CHOACHI 3 3 1923 18 3
## 7 CUNDINAMARCA CHOCONTA 1 2 2689 10 1
## 8 CUNDINAMARCA COGUA 1 2 2600 14 1
## 9 CUNDINAMARCA COTA 2 2 2566 14 2
## 10 CUNDINAMARCA CUCUNUBA 2 NA 2590 14 2Cargar base de datos de productividad DATOS ABIERTOS
PRO.TOT <- read_excel("RENDIMIENTO DE CULTIVOS EN COLOMBIA POR AÑO 1.10.23.xlsx", sheet = "RENDIMIENTO DE CULTIVOS EN COLO")
# Quitar las tildes a la columna de municipios
PRO.TOT$MUN <- toupper(stri_trans_general(PRO.TOT$MUN,"Latin-ASCII"))
#
PRO.CUN1 <- PRO.TOT %>%
filter(., Depertamento==c("CUNDINAMARCA", "BOGOTÁ"), Cultivo=="PAPA") %>%
group_by(., MUN) %>%
dplyr::summarise(.,
Asem1=mean(`Area sembrada`, na.rm = TRUE),
Acos1=mean(`Area cosechada`, na.rm = TRUE),
Prod1=mean(`Produccion`, na.rm = TRUE),
Rmax1=max(`Rendimiento`, na.rm = TRUE),
Rmed1=median(`Rendimiento`, na.rm = TRUE),
Rmea1=mean(`Rendimiento`, na.rm = TRUE)) %>%
ungroup() %>%
right_join(., AP.CUN, by="MUN") %>%
dplyr::mutate(.,
ALT=as.numeric(ALTITUD),
TEM=as.numeric(TEMPERATURA),
IDPM=as.numeric(IDPM) ) %>%
select(., !c("IDPMi", "IDPMii", "DEPARTAMENTO", "ALTITUD", "TEMPERATURA") )
head(PRO.CUN1, n=10)## # A tibble: 10 × 10
## MUN Asem1 Acos1 Prod1 Rmax1 Rmed1 Rmea1 IDPM ALT TEM
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 CAJICA 20.1 12.2 251. 30 20 20.7 3 2558 14
## 2 CARMEN DE CARUPA 945. 855. 15280. 32.6 15 16.0 1 2600 12
## 3 CHIA 125. 118. 2035 18 17 17.1 2 2600 14
## 4 CHIPAQUE 118. 116. 1722. 20 14.5 14.0 3 2400 13
## 5 CHOACHI 41.8 41.5 631. 18 15 15.5 3 1923 18
## 6 CHOCONTA 923. 859. 21884. 30.7 20 22.1 1 2689 10
## 7 COGUA 792. 789 15641. 25 19.7 21.1 1 2600 14
## 8 COTA 71.4 71.1 1960. 40 23 23.3 2 2566 14
## 9 CUCUNUBA 158. 158. 2422. 18.1 15.2 15.5 2 2590 14
## 10 EL ROSAL 276. 259. 4601. 18 18 17.4 1 2685 12Cargar base de datos de productividad PAPA TETRAPLOIDE papa de año (AGRONET)
PRO.CUN2 <- read_excel("agronet, 2022. Produccion y rendimiento cundinamarca.xlsx")
# Quitar las tildes a la columna de municipios
PRO.CUN2$MUN <- toupper(stri_trans_general(PRO.CUN2$MUN,"Latin-ASCII"))
#
PRO.CUN2 <- PRO.CUN2 %>%
group_by(., MUN) %>%
dplyr::summarise(.,
Asem2=mean(`Area sembrada`, na.rm = TRUE),
Acos2=mean(`Area cosechada`, na.rm = TRUE),
Prod2=mean(`Produccion`, na.rm = TRUE),
Rmax2=max(`Rendimiento`, na.rm = TRUE),
Rmed2=median(`Rendimiento`, na.rm = TRUE),
Rmea2=mean(`Rendimiento`, na.rm = TRUE)) %>%
ungroup() %>%
right_join(., AP.CUN, by="MUN") %>%
dplyr::mutate(.,
ALT=as.numeric(ALTITUD),
TEM=as.numeric(TEMPERATURA),
IDPM=as.numeric(IDPM) ) %>%
select(., !c("IDPMi", "IDPMii", "DEPARTAMENTO", "ALTITUD", "TEMPERATURA") )
head(PRO.CUN2, n=10)## # A tibble: 10 × 10
## MUN Asem2 Acos2 Prod2 Rmax2 Rmed2 Rmea2 IDPM ALT TEM
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 CAJICA 39.4 26.5 514. 30 19.5 19.5 3 2558 14
## 2 CARMEN DE CARUPA 2276. 2083. 39427 28.7 19.0 19.3 1 2600 12
## 3 CHIA 189. 195. 3465. 25 18 18.2 2 2600 14
## 4 CHIPAQUE 121. 156. 2366. 30 15.7 16.4 3 2400 13
## 5 CHOACHI 70.1 71.8 1080. 18 15.5 15.4 3 1923 18
## 6 CHOCONTA 3448. 3091. 84074. 30.7 28.9 27.2 1 2689 10
## 7 COGUA 1309. 1362. 27124. 25 20.0 20.8 1 2600 14
## 8 COTA 146. 144. 4054. 36.4 26 26.6 2 2566 14
## 9 CUCUNUBA 292. 290. 4856 30 15.7 17.3 2 2590 14
## 10 EL ROSAL 718. 675. 13333. 24 18 19.3 1 2685 12Cargar base de datos de productividad PAPA DIPLOIDE
PRO.CUN3 <- read_excel("agronet, 2022. Produccion y rendimiento cundinamarca criolla.xlsx")
# Quitar las tildes a la columna de municipios
PRO.CUN3$MUN <- toupper(stri_trans_general(PRO.CUN3$MUN,"Latin-ASCII"))
#
PRO.CUN3 <- PRO.CUN3 %>%
group_by(., MUN) %>%
dplyr::summarise(.,
Asem3=mean(`Area sembrada`, na.rm = TRUE),
Acos3=mean(`Area cosechada`, na.rm = TRUE),
Prod3=mean(`Produccion`, na.rm = TRUE),
Rmax3=max(`Rendimiento`, na.rm = TRUE),
Rmed3=median(`Rendimiento`, na.rm = TRUE),
Rmea3=mean(`Rendimiento`, na.rm = TRUE)) %>%
ungroup() %>%
right_join(., AP.CUN, by="MUN") %>%
dplyr::mutate(.,
ALT=as.numeric(ALTITUD),
TEM=as.numeric(TEMPERATURA),
IDPM=as.numeric(IDPM) ) %>%
select(., !c("IDPMi", "IDPMii", "DEPARTAMENTO", "ALTITUD", "TEMPERATURA") )
head(PRO.CUN3, n=10)## # A tibble: 10 × 10
## MUN Asem3 Acos3 Prod3 Rmax3 Rmed3 Rmea3 IDPM ALT TEM
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 CAJICA 10 12 140 11.7 11.7 11.7 3 2558 14
## 2 CARMEN DE CARUPA 77.6 69.4 956. 20 14.0 13.1 1 2600 12
## 3 CHIPAQUE 38.3 46.8 625. 20 11.7 12.5 3 2400 13
## 4 CHOACHI 76.5 70.1 987. 16 14.5 14.4 3 1923 18
## 5 CHOCONTA 150. 138. 2327. 25.6 17.9 18.4 1 2689 10
## 6 COGUA 453. 452. 8403. 19.8 17.9 12.6 1 2600 14
## 7 COTA 22 22 338. 16 15.7 15.6 2 2566 14
## 8 CUCUNUBA 21.2 21.2 353. 22 16.7 17.3 2 2590 14
## 9 EL ROSAL 459 436. 7616 18 18 17.3 1 2685 12
## 10 FACATATIVA 127. 112. 1988. 23.0 17 17.6 1 2586 19Unificamos las bases de datos de productividad
PRO.CUN <- PRO.CUN1 %>%
full_join(., PRO.CUN2) %>%
full_join(., PRO.CUN3) %>%
select(sort(names(.))); PRO.CUN## # A tibble: 47 × 22
## Acos1 Acos2 Acos3 ALT Asem1 Asem2 Asem3 IDPM MUN Prod1 Prod2 Prod3
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <chr> <dbl> <dbl> <dbl>
## 1 12.2 26.5 12 2558 20.1 39.4 10 3 CAJICA 251. 514. 140
## 2 855. 2083. 69.4 2600 945. 2276. 77.6 1 CARMEN… 15280. 39427 956.
## 3 118. 195. NA 2600 125. 189. NA 2 CHIA 2035 3465. NA
## 4 116. 156. 46.8 2400 118. 121. 38.3 3 CHIPAQ… 1722. 2366. 625.
## 5 41.5 71.8 70.1 1923 41.8 70.1 76.5 3 CHOACHI 631. 1080. 987.
## 6 859. 3091. 138. 2689 923. 3448. 150. 1 CHOCON… 21884. 84074. 2327.
## 7 789 1362. 452. 2600 792. 1309. 453. 1 COGUA 15641. 27124. 8403.
## 8 71.1 144. 22 2566 71.4 146. 22 2 COTA 1960. 4054. 338.
## 9 158. 290. 21.2 2590 158. 292. 21.2 2 CUCUNU… 2422. 4856 353.
## 10 259. 675. 436. 2685 276. 718. 459 1 EL ROS… 4601. 13333. 7616
## # ℹ 37 more rows
## # ℹ 10 more variables: Rmax1 <dbl>, Rmax2 <dbl>, Rmax3 <dbl>, Rmea1 <dbl>,
## # Rmea2 <dbl>, Rmea3 <dbl>, Rmed1 <dbl>, Rmed2 <dbl>, Rmed3 <dbl>, TEM <dbl>str(PRO.CUN)## tibble [47 × 22] (S3: tbl_df/tbl/data.frame)
## $ Acos1: num [1:47] 12.2 854.7 118.3 116.3 41.5 ...
## $ Acos2: num [1:47] 26.5 2082.7 195.3 155.5 71.8 ...
## $ Acos3: num [1:47] 12 69.4 NA 46.8 70.1 ...
## $ ALT : num [1:47] 2558 2600 2600 2400 1923 ...
## $ Asem1: num [1:47] 20.1 945.2 124.7 118.3 41.8 ...
## $ Asem2: num [1:47] 39.4 2276.4 189.3 121.3 70.1 ...
## $ Asem3: num [1:47] 10 77.6 NA 38.3 76.5 ...
## $ IDPM : num [1:47] 3 1 2 3 3 1 1 2 2 1 ...
## $ MUN : chr [1:47] "CAJICA" "CARMEN DE CARUPA" "CHIA" "CHIPAQUE" ...
## $ Prod1: num [1:47] 251 15280 2035 1722 631 ...
## $ Prod2: num [1:47] 514 39427 3465 2366 1080 ...
## $ Prod3: num [1:47] 140 956 NA 625 987 ...
## $ Rmax1: num [1:47] 30 32.6 18 20 18 ...
## $ Rmax2: num [1:47] 30 28.7 25 30 18 ...
## $ Rmax3: num [1:47] 11.7 20 NA 20 16 ...
## $ Rmea1: num [1:47] 20.7 16 17.1 14 15.5 ...
## $ Rmea2: num [1:47] 19.5 19.3 18.2 16.4 15.4 ...
## $ Rmea3: num [1:47] 11.7 13.1 NA 12.5 14.4 ...
## $ Rmed1: num [1:47] 20 15 17 14.5 15 ...
## $ Rmed2: num [1:47] 19.5 19 18 15.7 15.5 ...
## $ Rmed3: num [1:47] 11.7 13.9 NA 11.7 14.5 ...
## $ TEM : num [1:47] 14 12 14 13 18 10 14 14 14 12 ...Cargar base de datos de analisis de suelos
RAS.TOT <- read_excel("RESULTADOS ANALISIS DE SUELOS 1.10.23.xlsx", sheet = "RAS")
# Quitar las tildes a la columna de municipios
RAS.TOT$MUN <- toupper(stri_trans_general(RAS.TOT$MUN,"Latin-ASCII"))
head(RAS.TOT, n=10)## # A tibble: 10 × 31
## DEP MUN CUL EST TIE TOP DRE RIE FER FEC PH MAT FOS
## <chr> <chr> <chr> <chr> <chr> <chr> <chr> <chr> <chr> <lgl> <chr> <dbl> <chr>
## 1 CUND… FUNZA Uchu… POR … NO A… Ondu… Bueno No i… No i… NA 5.66 9.71 5.62…
## 2 CUND… BITU… Citr… POR … NO A… Ondu… Bueno No c… No i… NA 8.08 3.42 7.57…
## 3 CUND… VILL… Past… ESTA… NO I… Ondu… Bueno No c… No i… NA 5.87 2.34 16.9…
## 4 CUND… VILL… Past… ESTA… NO I… Ondu… Bueno No c… No i… NA 5.56 6.85 32.0…
## 5 CUND… BOGO… Papa… POR … NO A… Pend… Bueno No c… 15-1… NA 4.87 16.0 64.1…
## 6 HUILA GIGA… Agua… POR … NO A… Pend… Bueno No c… No h… NA 5.82 2.34 7.09…
## 7 META ACAC… Citr… ESTA… mas … Plano Bueno No c… No i… NA 4.68 2.12 9.16…
## 8 META ACAC… Agua… ESTA… mas … Plano Regu… No c… Cal … NA 4.51… 1.51 36.1…
## 9 META ACAC… Cacao ESTA… mas … Plano Bueno No c… 18-1… NA 5.15 2.10 12.0…
## 10 META ACAC… Cacao ESTA… de 5… Plano Bueno No c… Calf… NA 4.79 1.22 7.81…
## # ℹ 18 more variables: AZU <dbl>, ACI <chr>, ALU <chr>, CAL <chr>, MAG <chr>,
## # POT <chr>, SOD <chr>, CAP <dbl>, CON <dbl>, HIE <chr>, COB <chr>,
## # MAN <chr>, ZIN <chr>, BOR <dbl>, HIE2 <chr>, COB2 <chr>, MAN2 <chr>,
## # ZIN2 <chr>Extraemos los analisis de suelos para municipios con potencial de producción de papa
cultivos.calido <- c("Plátano", "Cacao", "Aguacate", "Pastos-Estrella", "Yuca", "Caña panelera/azucar", "Maracuyá", "Citricos","Arroz", "Pastos-brachiaria", "Ñame", "Piña", "Algodón", "Sacha Inchi", "Citricos-Limón", "Citricos-Naranjo", "Guayaba", "Mango", "Melón", "Palma de aceite", "Caucho", "Chirimoya", "Guadua", "Guanabana", "Sábila", "Papaya", "Ají", "Balu", "Pimentón", "Pastos-King grass","Pastos-Angleton"
)
RAS.CUN <- RAS.TOT %>%
mutate_all(~ ifelse(. == "ND", 0.0, .)) %>%
mutate_all(~ ifelse(. == "<0,01", 0.005, .)) %>%
mutate_all(~ ifelse(. == "<0,06", 0.030, .)) %>%
mutate_all(~ ifelse(. == "<0,09", 0.045, .)) %>%
mutate_all(~ ifelse(. == "<0,14", 0.070, .)) %>%
mutate_all(~ ifelse(. == "<0,20", 0.100, .)) %>%
mutate_all(~ ifelse(. == "<0,50", 0.250, .)) %>%
mutate_all(~ ifelse(. == "<0,55", 0.275, .)) %>%
mutate_all(~ ifelse(. == "<0,59", 0.295, .)) %>%
mutate_all(~ ifelse(. == "<1,00", 0.500, .)) %>%
mutate_all(~ ifelse(. == "<3,87", 1.935, .)) %>%
mutate_all(~ ifelse(. == "<3,80", 1.900, .)) %>%
mutate_all(~ ifelse(. == "<5,00", 2.500, .)) %>%
mutate_all(~ ifelse(. == "<10,00", 5.000, .)) %>%
mutate_all(~ ifelse(. == ">12,88", 12.88, .)) %>%
select(., !c(DEP, EST, TIE, TOP, DRE, RIE, FER, FEC, HIE2, COB2, MAN2, ZIN2) ) %>%
dplyr::right_join(., PRO.CUN, by = "MUN") %>%
mutate(.,
MUN=as.factor(MUN),
CUL=as.factor(CUL),
IDPM=as.numeric(IDPM),
ALT=as.numeric(ALT),
TEM=as.numeric(TEM),
Asem1=as.numeric(Asem1),
Asem2=as.numeric(Asem2),
Asem3=as.numeric(Asem3),
Acos1=as.numeric(Acos1),
Acos2=as.numeric(Acos2),
Acos3=as.numeric(Acos3),
Prod1=as.numeric(Prod1),
Prod2=as.numeric(Prod2),
Prod3=as.numeric(Prod3),
Rmax1=as.numeric(Rmax1),
Rmax2=as.numeric(Rmax2),
Rmax3=as.numeric(Rmax3),
Rmea1=as.numeric(Rmea1),
Rmea2=as.numeric(Rmea2),
Rmea3=as.numeric(Rmea3),
Rmed1=as.numeric(Rmed1),
Rmed2=as.numeric(Rmed2),
Rmed3=as.numeric(Rmed3),
pH=as.numeric(PH),
CEC=as.numeric(CAP),
EC=as.numeric(CON),
OM=as.numeric(MAT),
P=as.numeric(FOS),
K=as.numeric(POT),
Ca=as.numeric(CAL),
S=as.numeric(AZU),
Mg=as.numeric(MAG),
Na=as.numeric(SOD),
Fe=as.numeric(HIE),
Cu=as.numeric(COB),
Mn=as.numeric(MAN),
Zn=as.numeric(ZIN),
B=as.numeric(BOR),
EA=as.numeric(ACI),
Al=as.numeric(ALU) ) %>%
mutate(.,
Pav=(P*2.29), # FOSFORO DISPONIBLE EN PPM P2O5
N.tot=(OM*0.05), # NITROGENO EN PORCENTAJE
N.dis=((OM*0.05)*0.015)*10000, # NITROGENO DISPONIBLE CLIMA FRÍO EN ppm
BS=((Ca+Mg+K+Na)/CEC)*100, # SATURACION DE BASES
Alsat=(Al/CEC)*100, # SATURACION DE ALUMINIO
S.Ca=(Ca/CEC)*100, # SATURACION DE CALCIO
S.Mg=(Mg/CEC)*100, # SATURACION DE MAGNESIO
S.K=(K/CEC)*100, # SATURACION DE POTASIO
S.Na=(Na/CEC)*100, # SATURACION DE SODIO
Ca_Mg=(Ca/Mg), # RELACION CALCIO MAGNESIO
Mg_K=(Mg/K), # RELACION MAGNESIO POTASIO
Ca_K=(Ca/K), # RELACION CALCIO POTASIO
Ca.Mg_K=((Ca+Mg)/K), # RELACION CALCIO MAGNESIO CON POTASIO
Ca_B=(Ca/B), # RELACION CALCIO BORO
Fe_Mn=(Fe/Mn), # RELACION HIERRO MANGANESO
P_Zn=(P/Zn), # RELACION FOSFORO ZINC
Fe_Zn=(Fe/Zn) ) %>% # RELACION HIERRO ZINC
select(., !c(PH, MAT, FOS, AZU, ACI, ALU, CAL, MAG, POT, SOD, CAP, CON, HIE, COB, MAN, ZIN, BOR ) ) %>%
filter(., !(CUL %in% cultivos.calido)
)
head(RAS.CUN, n=10)## # A tibble: 10 × 57
## MUN CUL Acos1 Acos2 Acos3 ALT Asem1 Asem2 Asem3 IDPM Prod1 Prod2
## <fct> <fct> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 FUNZA Uchu… 408. 586. NA 2548 492. 651. NA 1 8161. 11779.
## 2 BOGOTA Papa… NA NA NA 2625 NA NA NA 3 NA NA
## 3 CHIA Frut… 118. 195. NA 2600 125. 189. NA 2 2035 3465.
## 4 CARMEN D… Papa… 855. 2083. 69.4 2600 945. 2276. 77.6 1 15280. 39427
## 5 CARMEN D… Papa… 855. 2083. 69.4 2600 945. 2276. 77.6 1 15280. 39427
## 6 CARMEN D… Papa… 855. 2083. 69.4 2600 945. 2276. 77.6 1 15280. 39427
## 7 CARMEN D… Papa… 855. 2083. 69.4 2600 945. 2276. 77.6 1 15280. 39427
## 8 CARMEN D… Arve… 855. 2083. 69.4 2600 945. 2276. 77.6 1 15280. 39427
## 9 CARMEN D… Past… 855. 2083. 69.4 2600 945. 2276. 77.6 1 15280. 39427
## 10 CARMEN D… Past… 855. 2083. 69.4 2600 945. 2276. 77.6 1 15280. 39427
## # ℹ 45 more variables: Prod3 <dbl>, Rmax1 <dbl>, Rmax2 <dbl>, Rmax3 <dbl>,
## # Rmea1 <dbl>, Rmea2 <dbl>, Rmea3 <dbl>, Rmed1 <dbl>, Rmed2 <dbl>,
## # Rmed3 <dbl>, TEM <dbl>, pH <dbl>, CEC <dbl>, EC <dbl>, OM <dbl>, P <dbl>,
## # K <dbl>, Ca <dbl>, S <dbl>, Mg <dbl>, Na <dbl>, Fe <dbl>, Cu <dbl>,
## # Mn <dbl>, Zn <dbl>, B <dbl>, EA <dbl>, Al <dbl>, Pav <dbl>, N.tot <dbl>,
## # N.dis <dbl>, BS <dbl>, Alsat <dbl>, S.Ca <dbl>, S.Mg <dbl>, S.K <dbl>,
## # S.Na <dbl>, Ca_Mg <dbl>, Mg_K <dbl>, Ca_K <dbl>, Ca.Mg_K <dbl>, …names(RAS.CUN)## [1] "MUN" "CUL" "Acos1" "Acos2" "Acos3" "ALT" "Asem1"
## [8] "Asem2" "Asem3" "IDPM" "Prod1" "Prod2" "Prod3" "Rmax1"
## [15] "Rmax2" "Rmax3" "Rmea1" "Rmea2" "Rmea3" "Rmed1" "Rmed2"
## [22] "Rmed3" "TEM" "pH" "CEC" "EC" "OM" "P"
## [29] "K" "Ca" "S" "Mg" "Na" "Fe" "Cu"
## [36] "Mn" "Zn" "B" "EA" "Al" "Pav" "N.tot"
## [43] "N.dis" "BS" "Alsat" "S.Ca" "S.Mg" "S.K" "S.Na"
## [50] "Ca_Mg" "Mg_K" "Ca_K" "Ca.Mg_K" "Ca_B" "Fe_Mn" "P_Zn"
## [57] "Fe_Zn"Conteo de observaciones por cultivo
head(sort(table(RAS.CUN$CUL), decreasing = T), n = 80)##
## Pastos Papa de año
## 999 427
## Hortalizas Frijol
## 153 116
## Papa Criolla Arveja
## 105 99
## Uchuva Café
## 96 95
## Maíz Pastos-raigrás
## 93 91
## Cebolla de bulbo Mora
## 72 65
## Tomate Fresa
## 63 59
## Quinua Pastos-kikuyo
## 51 44
## Zanahoria Forestal
## 39 34
## No indica Tomate de arbol
## 32 30
## Avena Cebolla de rama
## 29 24
## Gulupa Ajo
## 20 19
## Lechuga Cebolla
## 14 12
## Granadilla Frutales
## 12 11
## Lulo Arandanos
## 11 10
## Aromaticas Aromaticas-Hierbabuena
## 10 10
## Cannabis Feijoa
## 10 6
## Frambuesa Eucalipto
## 6 5
## Pastos-Kikuyo Soya
## 5 5
## Alfalfa Aromaticas-Calendula
## 4 4
## Caducifolios Calas
## 4 4
## Cebada Cebolla puerro
## 4 3
## flores Ornamentales Ahuyama
## 3 2
## Alcachofa Aromaticas-Laurel
## 2 2
## Arracacha Calabaza
## 2 2
## Citricos-Tangelo Coliflor
## 2 2
## Espinaca Follaje Ornamental
## 2 2
## Guisantes agrosilvopastoril
## 2 1
## Aromaticas-Estragón Aromaticas-Tomillo
## 1 1
## Banco de proteina forrajero Calabacín
## 1 1
## Cañamo cebolla
## 1 1
## Cebollin Cilantro
## 1 1
## Cubio Epifitas
## 1 1
## flores industrial-Rosa flores Ornamentales-Crisantemo
## 1 1
## flores Ornamentales-Hortencias flores Ornamentales-Snapdragon
## 1 1
## Follaje Ornamental-Brillantina follaje Ornamental-Ruscus
## 1 1
## Forestal-Pino Haba
## 1 1
## Helecho Huerta
## 1 1
## Solidago Stevia
## 1 1
## Aguacate Ají
## 0 0### NOTA ELIMINAR ANALISIS DE SUELOS RELACIONADOS CON CULTIVOS DE CLIMA CALIDOConteo de observaciones por municipio
head(sort(table(RAS.CUN$MUN), decreasing = T), n = 50)##
## BOGOTA CHOCONTA GUASCA CARMEN DE CARUPA
## 322 156 144 128
## SIBATE GUATAVITA CUCUNUBA MANTA
## 109 106 102 102
## GRANADA TAUSA SIMIJACA GUTIERREZ
## 97 95 92 90
## LENGUAZAQUE SUESCA LA CALERA PASCA
## 90 89 84 83
## SAN CAYETANO CHOACHI SAN BERNARDO FACATATIVA
## 82 73 73 70
## SUSA GUACHETA UNE CHIPAQUE
## 65 63 63 59
## SUBACHOQUE FUQUENE FUNZA SESQUILE
## 55 54 49 49
## PACHO COGUA VILLAPINZON FOSCA
## 42 35 34 32
## JUNIN UBAQUE TENJO TABIO
## 32 28 26 22
## ZIPAQUIRA CHIA SUTATAUSA EL ROSAL
## 22 21 21 20
## ZIPACON MADRID SOACHA SOPO
## 20 17 14 6
## CAJICA COTA MACHETA
## 5 5 5Agrupar variables por municipio
p.truncado <- 0.05
GR.CUN <- RAS.CUN %>% group_by(., MUN) %>%
dplyr::summarise(.,
pH=mean(pH, na.rm = TRUE, trim = p.truncado),
OM=mean(OM, na.rm = TRUE, trim = p.truncado),
N.tot=mean(N.tot, na.rm = TRUE, trim = p.truncado),
N.dis=mean(N.dis, na.rm = TRUE, trim = p.truncado),
P=mean(P, na.rm = TRUE, trim = p.truncado),
Pav=mean(Pav, na.rm = TRUE, trim = p.truncado),
S=mean(S, na.rm = TRUE, trim = p.truncado),
EA=mean(EA, na.rm = TRUE, trim = p.truncado),
Al=mean(Al, na.rm = TRUE, trim = p.truncado),
Ca=mean(Ca, na.rm = TRUE, trim = p.truncado),
Mg=mean(Mg, na.rm = TRUE, trim = p.truncado),
K=mean(K, na.rm = TRUE, trim = p.truncado),
Na=mean(Na, na.rm = TRUE, trim = p.truncado),
CEC=mean(CEC, na.rm = TRUE, trim = p.truncado),
EC=mean(EC, na.rm = TRUE, trim = p.truncado),
Fe=mean(Fe, na.rm = TRUE, trim = p.truncado),
Cu=mean(Cu, na.rm = TRUE, trim = p.truncado),
Mn=mean(Mn, na.rm = TRUE, trim = p.truncado),
Zn=mean(Zn, na.rm = TRUE, trim = p.truncado),
B=mean(B, na.rm = TRUE, trim = p.truncado),
BS=mean(BS, na.rm = TRUE, trim = p.truncado),
Alsat=mean(Alsat, na.rm = TRUE, trim = p.truncado),
S.Ca=mean(S.Ca, na.rm = TRUE, trim = p.truncado),
S.Mg=mean(S.Mg, na.rm = TRUE, trim = p.truncado),
S.K=mean(S.K, na.rm = TRUE, trim = p.truncado),
S.Na=mean(S.Na, na.rm = TRUE, trim = p.truncado),
`Ca/Mg`=mean(Ca_Mg, na.rm = TRUE, trim = p.truncado),
`Mg/K`=mean(Mg_K, na.rm = TRUE, trim = p.truncado),
`Ca/K`=mean(Ca_K, na.rm = TRUE, trim = p.truncado),
`Ca+Mg/K`=mean(Ca.Mg_K, na.rm = TRUE, trim = p.truncado),
`Ca/B`=mean(Ca_B, na.rm = TRUE, trim = p.truncado),
`Fe/Mn`=mean(Fe_Mn, na.rm = TRUE, trim = p.truncado),
`Fe/Zn`=mean(Fe_Zn, na.rm = TRUE, trim = p.truncado)) %>%
right_join(., PRO.CUN, by="MUN") %>%
textshape::column_to_rownames(., 1) %>%
select(., c(Rmed1, Rmed2, Rmed3, Rmea1, Rmea2, Rmea3, Rmax1, Rmax2, Rmax3,
Asem1, Asem2, Asem3, Acos1, Acos2, Acos3, Prod1, Prod2, Prod3,
IDPM, ALT, TEM,
pH, EA, Al, EC, CEC, OM,
N.tot, N.dis, P, Pav, K, Ca, Mg, S, Na, Fe, Cu, Mn, Zn, B,
BS, Alsat, S.Ca, S.Mg, S.K, S.Na,
`Ca/Mg`, `Mg/K`, `Ca/K`, `Ca+Mg/K`,
`Ca/B`, `Fe/Mn`, `Fe/Zn`) )
head(GR.CUN, n=10)## Rmed1 Rmed2 Rmed3 Rmea1 Rmea2 Rmea3 Rmax1 Rmax2
## BOGOTA NA NA NA NA NA NA NA NA
## CAJICA 20.00 19.510 11.670 20.72727 19.52375 11.67000 30.00 30.00
## CARMEN DE CARUPA 15.00 19.040 13.950 15.99538 19.32875 13.07273 32.57 28.74
## CHIA 17.00 18.000 NA 17.10222 18.21562 NA 18.00 25.00
## CHIPAQUE 14.50 15.670 11.690 14.04067 16.39133 12.45067 20.00 30.00
## CHOACHI 15.00 15.460 14.470 15.53412 15.41750 14.40231 18.00 18.00
## CHOCONTA 20.00 28.945 17.885 22.10111 27.16437 18.35500 30.69 30.71
## COGUA 19.67 20.020 17.890 21.05467 20.84250 12.57333 25.00 25.00
## COTA 23.00 26.000 15.725 23.27875 26.62875 15.61250 40.00 36.37
## CUCUNUBA 15.24 15.685 16.660 15.49250 17.27750 17.33000 18.13 30.00
## Rmax3 Asem1 Asem2 Asem3 Acos1 Acos2
## BOGOTA NA NA NA NA NA NA
## CAJICA 11.67 20.09091 39.35625 10.00000 12.18182 26.49375
## CARMEN DE CARUPA 20.00 945.23077 2276.43750 77.63636 854.69231 2082.68750
## CHIA NA 124.66667 189.31250 NA 118.33333 195.31250
## CHIPAQUE 20.00 118.33333 121.30000 38.33333 116.33333 155.50000
## CHOACHI 16.00 41.82353 70.06250 76.53846 41.47059 71.75000
## CHOCONTA 25.64 922.88889 3448.18750 149.56250 858.70370 3091.25000
## COGUA 19.83 792.26667 1308.68750 453.33333 789.00000 1361.68750
## COTA 16.00 71.37500 145.75438 22.00000 71.06250 144.26687
## CUCUNUBA 22.00 157.50000 291.50000 21.25000 157.50000 289.62500
## Acos3 Prod1 Prod2 Prod3 IDPM ALT TEM
## BOGOTA NA NA NA NA 3 2625 13.1
## CAJICA 12.00000 250.5455 513.6344 140.0000 3 2558 14.0
## CARMEN DE CARUPA 69.36364 15279.7692 39427.0000 956.2727 1 2600 12.0
## CHIA NA 2035.0000 3465.1875 NA 2 2600 14.0
## CHIPAQUE 46.80000 1722.1333 2366.3500 625.2013 3 2400 13.0
## CHOACHI 70.07692 630.7059 1079.5000 987.4615 3 1923 18.0
## CHOCONTA 138.43750 21884.0741 84073.7125 2326.6875 1 2689 10.0
## COGUA 452.33333 15640.6667 27124.3125 8403.3333 1 2600 14.0
## COTA 22.00000 1959.8750 4053.8387 338.5000 2 2566 14.0
## CUCUNUBA 21.25000 2422.5000 4856.0000 352.7500 2 2590 14.0
## pH EA Al EC CEC OM
## BOGOTA 5.559103 1.48842672 1.14813493 0.4744674 12.685202 9.778200
## CAJICA 5.316000 1.10836635 0.87816355 0.4607279 10.384645 8.094967
## CARMEN DE CARUPA 5.062155 2.06117455 1.63019473 0.3139291 6.417927 12.586212
## CHIA 5.918421 0.08288328 0.05627741 0.8380108 14.984885 7.087336
## CHIPAQUE 5.771273 0.76391403 0.61289557 0.5017144 12.069595 6.234734
## CHOACHI 5.151642 2.73579100 2.20641876 0.2094300 7.449557 7.019542
## CHOCONTA 5.114296 1.57266228 1.31436399 0.3881219 7.185019 9.893829
## COGUA 5.452424 1.13583515 0.94955086 0.3825780 8.582600 10.649070
## COTA 5.014000 1.68571680 1.16897724 1.2191779 11.399998 8.561286
## CUCUNUBA 5.555217 0.50128910 0.36834854 0.3394613 7.801094 6.126479
## N.tot N.dis P Pav K Ca
## BOGOTA 0.4889100 73.33650 55.97627 128.18565 0.8372324 7.609837
## CAJICA 0.4047484 60.71225 29.65147 67.90187 0.6229017 6.828402
## CARMEN DE CARUPA 0.6293106 94.39659 50.34010 115.27883 0.6813839 2.703054
## CHIA 0.3543668 53.15502 59.40440 136.03607 1.3264158 10.384599
## CHIPAQUE 0.3117367 46.76050 106.19090 243.17715 0.6996663 8.785649
## CHOACHI 0.3509771 52.64657 20.21929 46.30217 0.4339747 2.978155
## CHOCONTA 0.4946914 74.20372 51.97780 119.02915 0.7592882 3.426372
## COGUA 0.5324535 79.86803 48.10996 110.17181 0.5279060 5.137754
## COTA 0.4280643 64.20965 37.84977 86.67598 1.3810860 5.196928
## CUCUNUBA 0.3063240 45.94859 21.26958 48.70733 0.8593057 4.451287
## Mg S Na Fe Cu Mn
## BOGOTA 1.7951393 10.507277 0.12319214 569.1982 2.330431 7.651463
## CAJICA 1.5094093 13.102876 0.30956664 668.5700 2.600000 7.291400
## CARMEN DE CARUPA 0.7181425 8.190727 0.08674512 532.3260 1.506336 5.883974
## CHIA 2.5320848 31.323270 0.27498420 488.5985 3.551789 6.985933
## CHIPAQUE 1.4130209 12.127141 0.10290401 669.0530 4.207126 7.183408
## CHOACHI 0.8782449 6.669945 0.07743290 793.1288 1.318992 6.973152
## CHOCONTA 0.9684186 8.993108 0.17511828 520.5033 2.564419 9.027619
## COGUA 1.4431673 15.964870 0.13715575 763.5126 3.872457 6.770941
## COTA 2.7790936 20.586415 0.35717320 795.6744 4.328800 10.535400
## CUCUNUBA 1.5760503 11.428284 0.12308086 624.6912 2.737838 5.479770
## Zn B BS Alsat S.Ca S.Mg
## BOGOTA 5.463535 0.3275570 80.05188 14.7452755 55.16836 13.56820
## CAJICA 6.398600 0.3558244 80.67935 15.4288655 57.30971 13.82068
## CARMEN DE CARUPA 1.920793 0.2525678 66.48783 26.3535645 41.91207 11.13792
## CHIA 11.472134 0.4177221 99.31008 0.4654782 68.81959 16.91423
## CHIPAQUE 5.334073 0.3241120 89.93421 7.9391507 69.18755 12.30076
## CHOACHI 3.039539 0.1272001 59.04389 31.9970562 39.42074 11.50512
## CHOCONTA 2.880188 0.2669020 76.23524 19.5441738 47.32400 13.26507
## COGUA 5.599856 0.2983707 85.37701 12.1882056 59.12271 17.29694
## COTA 28.979400 0.4398446 72.98395 19.4199499 39.26746 19.91903
## CUCUNUBA 2.904986 0.2615827 92.33921 5.3468453 56.74571 20.73359
## S.K S.Na Ca/Mg Mg/K Ca/K Ca+Mg/K
## BOGOTA 7.243108 1.320546 4.207091 2.613154 10.679175 13.431593
## CAJICA 6.534641 3.014313 4.650541 4.041579 23.226207 27.267786
## CARMEN DE CARUPA 10.871190 1.392300 3.952252 1.244362 4.686444 5.946376
## CHIA 8.629590 1.929716 4.150889 2.791340 13.325539 16.116879
## CHIPAQUE 6.428829 0.959038 6.157040 2.107519 13.378390 15.714203
## CHOACHI 5.908447 1.127415 3.411689 2.389101 8.001158 10.422192
## CHOCONTA 10.551657 2.568885 4.155828 1.600995 5.741427 7.352359
## COGUA 6.084893 1.773415 4.019314 3.544209 11.898310 15.463167
## COTA 11.360041 2.437411 2.313750 1.995038 4.270404 6.265442
## CUCUNUBA 11.023580 1.583335 3.124588 2.223052 5.858726 8.184432
## Ca/B Fe/Mn Fe/Zn
## BOGOTA 22.91193 76.11250 159.48352
## CAJICA 22.10321 99.79597 224.82995
## CARMEN DE CARUPA 13.35808 101.88795 348.51124
## CHIA 28.28309 82.06782 44.59346
## CHIPAQUE 31.35667 106.15040 192.95755
## CHOACHI 44.79601 127.84451 469.15464
## CHOCONTA 15.30352 64.50734 279.89002
## COGUA 21.59959 127.99761 231.20179
## COTA 10.59602 137.53683 58.32030
## CUCUNUBA 23.36961 124.07392 283.95300Tabla resumen
tab.resumen <- GR.CUN %>%
select(., c(IDPM,TEM,ALT,Acos2,Acos3,Prod2,Prod3,Rmea2,Rmea3) ) %>%
mutate(., IDPM=as.numeric(IDPM) )
head(round(tab.resumen, digits = 1) )## IDPM TEM ALT Acos2 Acos3 Prod2 Prod3 Rmea2 Rmea3
## BOGOTA 3 13.1 2625 NA NA NA NA NA NA
## CAJICA 3 14.0 2558 26.5 12.0 513.6 140.0 19.5 11.7
## CARMEN DE CARUPA 1 12.0 2600 2082.7 69.4 39427.0 956.3 19.3 13.1
## CHIA 2 14.0 2600 195.3 NA 3465.2 NA 18.2 NA
## CHIPAQUE 3 13.0 2400 155.5 46.8 2366.3 625.2 16.4 12.5
## CHOACHI 3 18.0 1923 71.8 70.1 1079.5 987.5 15.4 14.42. A N A L I S I S E X P L O R A T O R I O
Boxplot exploratorios
# Obtener los nombres de las variables
variables <- colnames(GR.CUN)
par(mfrow=c(1,3))
# Iterar a través de las columnas y crear los gráficos de caja
for (i in variables) {
# Crear un gráfico de caja para la columna i
boxplot(x=GR.CUN[, i], main = i, ylab = paste("Valores de", i), horizontal = F)
}
Correlaciones con variables numericas y agrupadas por municipio
# Mediana
Rmed = cor( data.frame(GR.CUN[-c(1:3,20:54)]), method = "spearman", use = "complete.obs")
round(Rmed, digits = 2)## Rmea1 Rmea2 Rmea3 Rmax1 Rmax2 Rmax3 Asem1 Asem2 Asem3 Acos1 Acos2 Acos3
## Rmea1 1.00 0.80 0.64 0.61 0.47 0.38 0.37 0.40 0.14 0.39 0.39 0.13
## Rmea2 0.80 1.00 0.53 0.75 0.71 0.42 0.42 0.53 0.29 0.44 0.51 0.27
## Rmea3 0.64 0.53 1.00 0.27 0.27 0.72 0.25 0.32 0.23 0.27 0.32 0.25
## Rmax1 0.61 0.75 0.27 1.00 0.71 0.23 0.35 0.36 0.18 0.35 0.36 0.17
## Rmax2 0.47 0.71 0.27 0.71 1.00 0.23 0.23 0.30 -0.02 0.24 0.30 -0.03
## Rmax3 0.38 0.42 0.72 0.23 0.23 1.00 0.19 0.22 0.21 0.20 0.23 0.23
## Asem1 0.37 0.42 0.25 0.35 0.23 0.19 1.00 0.96 0.42 1.00 0.97 0.39
## Asem2 0.40 0.53 0.32 0.36 0.30 0.22 0.96 1.00 0.42 0.96 1.00 0.39
## Asem3 0.14 0.29 0.23 0.18 -0.02 0.21 0.42 0.42 1.00 0.43 0.43 0.99
## Acos1 0.39 0.44 0.27 0.35 0.24 0.20 1.00 0.96 0.43 1.00 0.97 0.40
## Acos2 0.39 0.51 0.32 0.36 0.30 0.23 0.97 1.00 0.43 0.97 1.00 0.40
## Acos3 0.13 0.27 0.25 0.17 -0.03 0.23 0.39 0.39 0.99 0.40 0.40 1.00
## Prod1 0.48 0.54 0.35 0.40 0.33 0.25 0.98 0.97 0.43 0.99 0.98 0.41
## Prod2 0.46 0.58 0.38 0.42 0.34 0.27 0.95 0.99 0.45 0.95 0.99 0.43
## Prod3 0.17 0.28 0.35 0.18 -0.01 0.32 0.40 0.41 0.98 0.41 0.42 0.98
## IDPM -0.44 -0.53 -0.38 -0.29 -0.30 -0.11 -0.79 -0.85 -0.42 -0.80 -0.84 -0.40
## Prod1 Prod2 Prod3 IDPM
## Rmea1 0.48 0.46 0.17 -0.44
## Rmea2 0.54 0.58 0.28 -0.53
## Rmea3 0.35 0.38 0.35 -0.38
## Rmax1 0.40 0.42 0.18 -0.29
## Rmax2 0.33 0.34 -0.01 -0.30
## Rmax3 0.25 0.27 0.32 -0.11
## Asem1 0.98 0.95 0.40 -0.79
## Asem2 0.97 0.99 0.41 -0.85
## Asem3 0.43 0.45 0.98 -0.42
## Acos1 0.99 0.95 0.41 -0.80
## Acos2 0.98 0.99 0.42 -0.84
## Acos3 0.41 0.43 0.98 -0.40
## Prod1 1.00 0.98 0.43 -0.83
## Prod2 0.98 1.00 0.45 -0.85
## Prod3 0.43 0.45 1.00 -0.40
## IDPM -0.83 -0.85 -0.40 1.00resumen.Rmed=colMeans(abs(Rmed))
# Media
Rmea = cor( data.frame(GR.CUN[4:6], GR.CUN[20:54]), method = "spearman", use = "complete.obs")
round(Rmea, digits = 2)## Rmea1 Rmea2 Rmea3 ALT TEM pH EA Al EC CEC OM N.tot
## Rmea1 1.00 0.80 0.64 0.24 -0.12 -0.06 -0.01 -0.03 0.18 0.02 0.24 0.24
## Rmea2 0.80 1.00 0.53 0.29 -0.11 0.02 -0.07 -0.09 0.11 0.00 0.28 0.28
## Rmea3 0.64 0.53 1.00 0.33 -0.17 0.11 -0.19 -0.22 0.14 0.01 0.22 0.22
## ALT 0.24 0.29 0.33 1.00 -0.60 -0.31 0.19 0.18 -0.09 -0.24 0.19 0.19
## TEM -0.12 -0.11 -0.17 -0.60 1.00 0.24 -0.08 -0.08 -0.18 -0.03 -0.24 -0.24
## pH -0.06 0.02 0.11 -0.31 0.24 1.00 -0.84 -0.81 0.40 0.54 -0.02 -0.02
## EA -0.01 -0.07 -0.19 0.19 -0.08 -0.84 1.00 0.99 -0.47 -0.48 -0.06 -0.06
## Al -0.03 -0.09 -0.22 0.18 -0.08 -0.81 0.99 1.00 -0.48 -0.47 -0.08 -0.08
## EC 0.18 0.11 0.14 -0.09 -0.18 0.40 -0.47 -0.48 1.00 0.87 0.35 0.35
## CEC 0.02 0.00 0.01 -0.24 -0.03 0.54 -0.48 -0.47 0.87 1.00 0.14 0.14
## OM 0.24 0.28 0.22 0.19 -0.24 -0.02 -0.06 -0.08 0.35 0.14 1.00 1.00
## N.tot 0.24 0.28 0.22 0.19 -0.24 -0.02 -0.06 -0.08 0.35 0.14 1.00 1.00
## N.dis 0.24 0.28 0.22 0.19 -0.24 -0.02 -0.06 -0.08 0.35 0.14 1.00 1.00
## P 0.33 0.17 0.25 0.04 -0.10 0.18 -0.14 -0.14 0.45 0.22 0.24 0.24
## Pav 0.33 0.17 0.25 0.04 -0.10 0.18 -0.14 -0.14 0.45 0.22 0.24 0.24
## K 0.30 0.14 0.38 -0.10 -0.10 0.31 -0.37 -0.36 0.64 0.52 0.20 0.20
## Ca 0.06 0.11 0.11 -0.19 -0.02 0.74 -0.73 -0.72 0.82 0.90 0.15 0.15
## Mg -0.02 -0.04 0.04 -0.20 -0.06 0.50 -0.62 -0.62 0.63 0.78 -0.08 -0.08
## S 0.07 -0.07 0.06 0.08 -0.22 0.14 -0.34 -0.34 0.75 0.64 0.33 0.33
## Na 0.12 0.02 0.02 0.07 -0.29 0.19 -0.28 -0.27 0.76 0.75 0.08 0.08
## Fe 0.17 -0.03 -0.05 0.19 -0.09 -0.47 0.50 0.48 -0.28 -0.13 -0.37 -0.37
## Cu 0.33 0.25 0.25 -0.09 -0.04 0.37 -0.49 -0.52 0.52 0.42 0.01 0.01
## Mn 0.43 0.40 0.26 0.28 -0.14 -0.31 0.16 0.11 -0.10 -0.12 -0.12 -0.12
## Zn 0.12 0.10 0.08 -0.10 -0.07 0.38 -0.57 -0.58 0.85 0.79 0.26 0.26
## B 0.20 0.05 0.17 0.04 -0.37 0.21 -0.35 -0.35 0.87 0.69 0.40 0.40
## BS 0.04 0.07 0.26 -0.17 0.02 0.80 -0.95 -0.95 0.48 0.52 0.10 0.10
## Alsat -0.07 -0.09 -0.30 0.14 -0.04 -0.78 0.95 0.96 -0.47 -0.49 -0.10 -0.10
## S.Ca 0.03 0.10 0.13 -0.20 0.06 0.84 -0.87 -0.86 0.52 0.56 0.13 0.13
## S.Mg -0.03 -0.07 0.08 -0.05 -0.09 0.24 -0.46 -0.47 0.09 0.18 -0.22 -0.22
## S.K 0.38 0.17 0.47 0.12 -0.10 -0.11 -0.05 -0.06 0.10 -0.13 0.17 0.17
## S.Na 0.10 0.02 -0.05 0.22 -0.34 -0.15 -0.03 -0.02 0.40 0.30 0.09 0.09
## Ca.Mg 0.17 0.24 0.13 -0.02 0.04 0.45 -0.30 -0.29 0.44 0.35 0.33 0.33
## Mg.K -0.31 -0.17 -0.31 -0.16 0.04 0.30 -0.32 -0.30 0.09 0.33 -0.26 -0.26
## Ca.K -0.24 -0.04 -0.28 -0.06 0.01 0.55 -0.48 -0.45 0.35 0.48 0.00 0.00
## Ca.Mg.K -0.30 -0.10 -0.29 -0.11 0.05 0.56 -0.49 -0.46 0.29 0.48 -0.07 -0.07
## Ca.B -0.21 -0.07 -0.11 -0.39 0.50 0.65 -0.39 -0.39 0.05 0.43 -0.30 -0.30
## Fe.Mn -0.11 -0.15 -0.29 0.14 -0.04 -0.41 0.54 0.55 -0.34 -0.15 -0.44 -0.44
## Fe.Zn -0.08 -0.15 -0.08 0.22 -0.04 -0.51 0.62 0.61 -0.76 -0.58 -0.43 -0.43
## N.dis P Pav K Ca Mg S Na Fe Cu Mn Zn
## Rmea1 0.24 0.33 0.33 0.30 0.06 -0.02 0.07 0.12 0.17 0.33 0.43 0.12
## Rmea2 0.28 0.17 0.17 0.14 0.11 -0.04 -0.07 0.02 -0.03 0.25 0.40 0.10
## Rmea3 0.22 0.25 0.25 0.38 0.11 0.04 0.06 0.02 -0.05 0.25 0.26 0.08
## ALT 0.19 0.04 0.04 -0.10 -0.19 -0.20 0.08 0.07 0.19 -0.09 0.28 -0.10
## TEM -0.24 -0.10 -0.10 -0.10 -0.02 -0.06 -0.22 -0.29 -0.09 -0.04 -0.14 -0.07
## pH -0.02 0.18 0.18 0.31 0.74 0.50 0.14 0.19 -0.47 0.37 -0.31 0.38
## EA -0.06 -0.14 -0.14 -0.37 -0.73 -0.62 -0.34 -0.28 0.50 -0.49 0.16 -0.57
## Al -0.08 -0.14 -0.14 -0.36 -0.72 -0.62 -0.34 -0.27 0.48 -0.52 0.11 -0.58
## EC 0.35 0.45 0.45 0.64 0.82 0.63 0.75 0.76 -0.28 0.52 -0.10 0.85
## CEC 0.14 0.22 0.22 0.52 0.90 0.78 0.64 0.75 -0.13 0.42 -0.12 0.79
## OM 1.00 0.24 0.24 0.20 0.15 -0.08 0.33 0.08 -0.37 0.01 -0.12 0.26
## N.tot 1.00 0.24 0.24 0.20 0.15 -0.08 0.33 0.08 -0.37 0.01 -0.12 0.26
## N.dis 1.00 0.24 0.24 0.20 0.15 -0.08 0.33 0.08 -0.37 0.01 -0.12 0.26
## P 0.24 1.00 1.00 0.44 0.30 -0.18 0.35 0.07 -0.17 0.50 -0.10 0.26
## Pav 0.24 1.00 1.00 0.44 0.30 -0.18 0.35 0.07 -0.17 0.50 -0.10 0.26
## K 0.20 0.44 0.44 1.00 0.47 0.40 0.49 0.50 -0.28 0.32 -0.05 0.46
## Ca 0.15 0.30 0.30 0.47 1.00 0.76 0.57 0.61 -0.28 0.57 -0.11 0.80
## Mg -0.08 -0.18 -0.18 0.40 0.76 1.00 0.58 0.72 -0.03 0.39 0.05 0.70
## S 0.33 0.35 0.35 0.49 0.57 0.58 1.00 0.64 -0.05 0.48 -0.06 0.74
## Na 0.08 0.07 0.07 0.50 0.61 0.72 0.64 1.00 -0.04 0.19 -0.09 0.62
## Fe -0.37 -0.17 -0.17 -0.28 -0.28 -0.03 -0.05 -0.04 1.00 0.01 0.48 -0.27
## Cu 0.01 0.50 0.50 0.32 0.57 0.39 0.48 0.19 0.01 1.00 0.19 0.63
## Mn -0.12 -0.10 -0.10 -0.05 -0.11 0.05 -0.06 -0.09 0.48 0.19 1.00 0.01
## Zn 0.26 0.26 0.26 0.46 0.80 0.70 0.74 0.62 -0.27 0.63 0.01 1.00
## B 0.40 0.54 0.54 0.63 0.64 0.50 0.84 0.70 -0.19 0.48 -0.13 0.72
## BS 0.10 0.16 0.16 0.40 0.73 0.64 0.36 0.30 -0.39 0.53 -0.09 0.56
## Alsat -0.10 -0.14 -0.14 -0.39 -0.71 -0.63 -0.35 -0.28 0.38 -0.54 0.04 -0.57
## S.Ca 0.13 0.35 0.35 0.25 0.82 0.49 0.29 0.24 -0.40 0.59 -0.17 0.59
## S.Mg -0.22 -0.51 -0.51 0.09 0.21 0.72 0.24 0.32 0.07 0.20 0.22 0.28
## S.K 0.17 0.33 0.33 0.73 -0.13 -0.07 0.16 -0.02 -0.09 0.12 0.14 -0.02
## S.Na 0.09 -0.13 -0.13 0.23 0.18 0.41 0.41 0.81 -0.03 -0.08 -0.05 0.31
## Ca.Mg 0.33 0.76 0.76 0.20 0.52 -0.13 0.16 0.03 -0.34 0.39 -0.19 0.30
## Mg.K -0.26 -0.57 -0.57 -0.38 0.36 0.60 0.11 0.35 0.01 0.02 -0.09 0.27
## Ca.K 0.00 0.00 0.00 -0.25 0.65 0.38 0.20 0.30 -0.23 0.22 -0.24 0.45
## Ca.Mg.K -0.07 -0.13 -0.13 -0.29 0.63 0.47 0.17 0.31 -0.19 0.20 -0.24 0.42
## Ca.B -0.30 -0.20 -0.20 -0.08 0.46 0.38 -0.19 0.04 -0.08 0.04 -0.05 0.13
## Fe.Mn -0.44 -0.23 -0.23 -0.39 -0.30 -0.12 -0.15 -0.02 0.79 -0.18 0.00 -0.32
## Fe.Zn -0.43 -0.36 -0.36 -0.49 -0.70 -0.42 -0.51 -0.39 0.66 -0.53 0.19 -0.79
## B BS Alsat S.Ca S.Mg S.K S.Na Ca.Mg Mg.K Ca.K Ca.Mg.K
## Rmea1 0.20 0.04 -0.07 0.03 -0.03 0.38 0.10 0.17 -0.31 -0.24 -0.30
## Rmea2 0.05 0.07 -0.09 0.10 -0.07 0.17 0.02 0.24 -0.17 -0.04 -0.10
## Rmea3 0.17 0.26 -0.30 0.13 0.08 0.47 -0.05 0.13 -0.31 -0.28 -0.29
## ALT 0.04 -0.17 0.14 -0.20 -0.05 0.12 0.22 -0.02 -0.16 -0.06 -0.11
## TEM -0.37 0.02 -0.04 0.06 -0.09 -0.10 -0.34 0.04 0.04 0.01 0.05
## pH 0.21 0.80 -0.78 0.84 0.24 -0.11 -0.15 0.45 0.30 0.55 0.56
## EA -0.35 -0.95 0.95 -0.87 -0.46 -0.05 -0.03 -0.30 -0.32 -0.48 -0.49
## Al -0.35 -0.95 0.96 -0.86 -0.47 -0.06 -0.02 -0.29 -0.30 -0.45 -0.46
## EC 0.87 0.48 -0.47 0.52 0.09 0.10 0.40 0.44 0.09 0.35 0.29
## CEC 0.69 0.52 -0.49 0.56 0.18 -0.13 0.30 0.35 0.33 0.48 0.48
## OM 0.40 0.10 -0.10 0.13 -0.22 0.17 0.09 0.33 -0.26 0.00 -0.07
## N.tot 0.40 0.10 -0.10 0.13 -0.22 0.17 0.09 0.33 -0.26 0.00 -0.07
## N.dis 0.40 0.10 -0.10 0.13 -0.22 0.17 0.09 0.33 -0.26 0.00 -0.07
## P 0.54 0.16 -0.14 0.35 -0.51 0.33 -0.13 0.76 -0.57 0.00 -0.13
## Pav 0.54 0.16 -0.14 0.35 -0.51 0.33 -0.13 0.76 -0.57 0.00 -0.13
## K 0.63 0.40 -0.39 0.25 0.09 0.73 0.23 0.20 -0.38 -0.25 -0.29
## Ca 0.64 0.73 -0.71 0.82 0.21 -0.13 0.18 0.52 0.36 0.65 0.63
## Mg 0.50 0.64 -0.63 0.49 0.72 -0.07 0.41 -0.13 0.60 0.38 0.47
## S 0.84 0.36 -0.35 0.29 0.24 0.16 0.41 0.16 0.11 0.20 0.17
## Na 0.70 0.30 -0.28 0.24 0.32 -0.02 0.81 0.03 0.35 0.30 0.31
## Fe -0.19 -0.39 0.38 -0.40 0.07 -0.09 -0.03 -0.34 0.01 -0.23 -0.19
## Cu 0.48 0.53 -0.54 0.59 0.20 0.12 -0.08 0.39 0.02 0.22 0.20
## Mn -0.13 -0.09 0.04 -0.17 0.22 0.14 -0.05 -0.19 -0.09 -0.24 -0.24
## Zn 0.72 0.56 -0.57 0.59 0.28 -0.02 0.31 0.30 0.27 0.45 0.42
## B 1.00 0.36 -0.34 0.35 0.07 0.22 0.46 0.36 -0.02 0.19 0.15
## BS 0.36 1.00 -0.99 0.89 0.49 0.06 0.00 0.27 0.26 0.40 0.42
## Alsat -0.34 -0.99 1.00 -0.87 -0.49 -0.08 0.01 -0.26 -0.25 -0.38 -0.40
## S.Ca 0.35 0.89 -0.87 1.00 0.18 -0.16 -0.12 0.60 0.25 0.64 0.61
## S.Mg 0.07 0.49 -0.49 0.18 1.00 0.02 0.33 -0.60 0.56 0.03 0.18
## S.K 0.22 0.06 -0.08 -0.16 0.02 1.00 0.00 -0.07 -0.70 -0.70 -0.74
## S.Na 0.46 0.00 0.01 -0.12 0.33 0.00 1.00 -0.24 0.34 0.11 0.14
## Ca.Mg 0.36 0.27 -0.26 0.60 -0.60 -0.07 -0.24 1.00 -0.27 0.47 0.33
## Mg.K -0.02 0.26 -0.25 0.25 0.56 -0.70 0.34 -0.27 1.00 0.67 0.79
## Ca.K 0.19 0.40 -0.38 0.64 0.03 -0.70 0.11 0.47 0.67 1.00 0.98
## Ca.Mg.K 0.15 0.42 -0.40 0.61 0.18 -0.74 0.14 0.33 0.79 0.98 1.00
## Ca.B -0.25 0.39 -0.38 0.46 0.11 -0.42 -0.27 0.19 0.44 0.48 0.54
## Fe.Mn -0.23 -0.46 0.48 -0.41 -0.07 -0.29 0.02 -0.31 0.15 -0.05 0.00
## Fe.Zn -0.63 -0.56 0.56 -0.62 -0.10 -0.07 -0.18 -0.45 -0.10 -0.38 -0.33
## Ca.B Fe.Mn Fe.Zn
## Rmea1 -0.21 -0.11 -0.08
## Rmea2 -0.07 -0.15 -0.15
## Rmea3 -0.11 -0.29 -0.08
## ALT -0.39 0.14 0.22
## TEM 0.50 -0.04 -0.04
## pH 0.65 -0.41 -0.51
## EA -0.39 0.54 0.62
## Al -0.39 0.55 0.61
## EC 0.05 -0.34 -0.76
## CEC 0.43 -0.15 -0.58
## OM -0.30 -0.44 -0.43
## N.tot -0.30 -0.44 -0.43
## N.dis -0.30 -0.44 -0.43
## P -0.20 -0.23 -0.36
## Pav -0.20 -0.23 -0.36
## K -0.08 -0.39 -0.49
## Ca 0.46 -0.30 -0.70
## Mg 0.38 -0.12 -0.42
## S -0.19 -0.15 -0.51
## Na 0.04 -0.02 -0.39
## Fe -0.08 0.79 0.66
## Cu 0.04 -0.18 -0.53
## Mn -0.05 0.00 0.19
## Zn 0.13 -0.32 -0.79
## B -0.25 -0.23 -0.63
## BS 0.39 -0.46 -0.56
## Alsat -0.38 0.48 0.56
## S.Ca 0.46 -0.41 -0.62
## S.Mg 0.11 -0.07 -0.10
## S.K -0.42 -0.29 -0.07
## S.Na -0.27 0.02 -0.18
## Ca.Mg 0.19 -0.31 -0.45
## Mg.K 0.44 0.15 -0.10
## Ca.K 0.48 -0.05 -0.38
## Ca.Mg.K 0.54 0.00 -0.33
## Ca.B 1.00 -0.10 -0.06
## Fe.Mn -0.10 1.00 0.61
## Fe.Zn -0.06 0.61 1.00resumen.Rmea=colMeans(abs(Rmea))
# Maximo
Rmax = cor( data.frame(GR.CUN[7:9], GR.CUN[20:54]), method = "spearman", use = "complete.obs")
round(Rmax, digits = 2)## Rmax1 Rmax2 Rmax3 ALT TEM pH EA Al EC CEC OM N.tot
## Rmax1 1.00 0.71 0.23 0.24 -0.06 -0.17 0.10 0.07 -0.12 -0.25 0.24 0.24
## Rmax2 0.71 1.00 0.23 0.27 -0.13 0.03 -0.09 -0.11 -0.10 -0.14 0.08 0.08
## Rmax3 0.23 0.23 1.00 0.25 -0.10 0.23 -0.29 -0.32 0.00 -0.13 0.09 0.09
## ALT 0.24 0.27 0.25 1.00 -0.60 -0.31 0.19 0.18 -0.09 -0.24 0.19 0.19
## TEM -0.06 -0.13 -0.10 -0.60 1.00 0.24 -0.08 -0.08 -0.18 -0.03 -0.24 -0.24
## pH -0.17 0.03 0.23 -0.31 0.24 1.00 -0.84 -0.81 0.40 0.54 -0.02 -0.02
## EA 0.10 -0.09 -0.29 0.19 -0.08 -0.84 1.00 0.99 -0.47 -0.48 -0.06 -0.06
## Al 0.07 -0.11 -0.32 0.18 -0.08 -0.81 0.99 1.00 -0.48 -0.47 -0.08 -0.08
## EC -0.12 -0.10 0.00 -0.09 -0.18 0.40 -0.47 -0.48 1.00 0.87 0.35 0.35
## CEC -0.25 -0.14 -0.13 -0.24 -0.03 0.54 -0.48 -0.47 0.87 1.00 0.14 0.14
## OM 0.24 0.08 0.09 0.19 -0.24 -0.02 -0.06 -0.08 0.35 0.14 1.00 1.00
## N.tot 0.24 0.08 0.09 0.19 -0.24 -0.02 -0.06 -0.08 0.35 0.14 1.00 1.00
## N.dis 0.24 0.08 0.09 0.19 -0.24 -0.02 -0.06 -0.08 0.35 0.14 1.00 1.00
## P 0.07 -0.09 0.32 0.04 -0.10 0.18 -0.14 -0.14 0.45 0.22 0.24 0.24
## Pav 0.07 -0.09 0.32 0.04 -0.10 0.18 -0.14 -0.14 0.45 0.22 0.24 0.24
## K -0.06 0.00 0.19 -0.10 -0.10 0.31 -0.37 -0.36 0.64 0.52 0.20 0.20
## Ca -0.17 -0.04 0.11 -0.19 -0.02 0.74 -0.73 -0.72 0.82 0.90 0.15 0.15
## Mg -0.22 -0.07 -0.04 -0.20 -0.06 0.50 -0.62 -0.62 0.63 0.78 -0.08 -0.08
## S -0.23 -0.30 -0.03 0.08 -0.22 0.14 -0.34 -0.34 0.75 0.64 0.33 0.33
## Na -0.20 -0.07 -0.18 0.07 -0.29 0.19 -0.28 -0.27 0.76 0.75 0.08 0.08
## Fe 0.11 0.01 -0.20 0.19 -0.09 -0.47 0.50 0.48 -0.28 -0.13 -0.37 -0.37
## Cu 0.05 -0.05 0.33 -0.09 -0.04 0.37 -0.49 -0.52 0.52 0.42 0.01 0.01
## Mn 0.47 0.31 0.29 0.28 -0.14 -0.31 0.16 0.11 -0.10 -0.12 -0.12 -0.12
## Zn -0.15 -0.11 -0.02 -0.10 -0.07 0.38 -0.57 -0.58 0.85 0.79 0.26 0.26
## B -0.11 -0.14 0.01 0.04 -0.37 0.21 -0.35 -0.35 0.87 0.69 0.40 0.40
## BS -0.12 0.08 0.34 -0.17 0.02 0.80 -0.95 -0.95 0.48 0.52 0.10 0.10
## Alsat 0.11 -0.10 -0.37 0.14 -0.04 -0.78 0.95 0.96 -0.47 -0.49 -0.10 -0.10
## S.Ca -0.13 0.01 0.28 -0.20 0.06 0.84 -0.87 -0.86 0.52 0.56 0.13 0.13
## S.Mg -0.07 0.05 0.08 -0.05 -0.09 0.24 -0.46 -0.47 0.09 0.18 -0.22 -0.22
## S.K 0.18 0.11 0.31 0.12 -0.10 -0.11 -0.05 -0.06 0.10 -0.13 0.17 0.17
## S.Na -0.07 0.01 -0.21 0.22 -0.34 -0.15 -0.03 -0.02 0.40 0.30 0.09 0.09
## Ca.Mg 0.07 0.02 0.27 -0.02 0.04 0.45 -0.30 -0.29 0.44 0.35 0.33 0.33
## Mg.K -0.24 -0.04 -0.23 -0.16 0.04 0.30 -0.32 -0.30 0.09 0.33 -0.26 -0.26
## Ca.K -0.21 -0.01 -0.07 -0.06 0.01 0.55 -0.48 -0.45 0.35 0.48 0.00 0.00
## Ca.Mg.K -0.25 -0.04 -0.09 -0.11 0.05 0.56 -0.49 -0.46 0.29 0.48 -0.07 -0.07
## Ca.B -0.21 -0.02 -0.04 -0.39 0.50 0.65 -0.39 -0.39 0.05 0.43 -0.30 -0.30
## Fe.Mn -0.06 0.01 -0.40 0.14 -0.04 -0.41 0.54 0.55 -0.34 -0.15 -0.44 -0.44
## Fe.Zn 0.07 0.06 -0.07 0.22 -0.04 -0.51 0.62 0.61 -0.76 -0.58 -0.43 -0.43
## N.dis P Pav K Ca Mg S Na Fe Cu Mn Zn
## Rmax1 0.24 0.07 0.07 -0.06 -0.17 -0.22 -0.23 -0.20 0.11 0.05 0.47 -0.15
## Rmax2 0.08 -0.09 -0.09 0.00 -0.04 -0.07 -0.30 -0.07 0.01 -0.05 0.31 -0.11
## Rmax3 0.09 0.32 0.32 0.19 0.11 -0.04 -0.03 -0.18 -0.20 0.33 0.29 -0.02
## ALT 0.19 0.04 0.04 -0.10 -0.19 -0.20 0.08 0.07 0.19 -0.09 0.28 -0.10
## TEM -0.24 -0.10 -0.10 -0.10 -0.02 -0.06 -0.22 -0.29 -0.09 -0.04 -0.14 -0.07
## pH -0.02 0.18 0.18 0.31 0.74 0.50 0.14 0.19 -0.47 0.37 -0.31 0.38
## EA -0.06 -0.14 -0.14 -0.37 -0.73 -0.62 -0.34 -0.28 0.50 -0.49 0.16 -0.57
## Al -0.08 -0.14 -0.14 -0.36 -0.72 -0.62 -0.34 -0.27 0.48 -0.52 0.11 -0.58
## EC 0.35 0.45 0.45 0.64 0.82 0.63 0.75 0.76 -0.28 0.52 -0.10 0.85
## CEC 0.14 0.22 0.22 0.52 0.90 0.78 0.64 0.75 -0.13 0.42 -0.12 0.79
## OM 1.00 0.24 0.24 0.20 0.15 -0.08 0.33 0.08 -0.37 0.01 -0.12 0.26
## N.tot 1.00 0.24 0.24 0.20 0.15 -0.08 0.33 0.08 -0.37 0.01 -0.12 0.26
## N.dis 1.00 0.24 0.24 0.20 0.15 -0.08 0.33 0.08 -0.37 0.01 -0.12 0.26
## P 0.24 1.00 1.00 0.44 0.30 -0.18 0.35 0.07 -0.17 0.50 -0.10 0.26
## Pav 0.24 1.00 1.00 0.44 0.30 -0.18 0.35 0.07 -0.17 0.50 -0.10 0.26
## K 0.20 0.44 0.44 1.00 0.47 0.40 0.49 0.50 -0.28 0.32 -0.05 0.46
## Ca 0.15 0.30 0.30 0.47 1.00 0.76 0.57 0.61 -0.28 0.57 -0.11 0.80
## Mg -0.08 -0.18 -0.18 0.40 0.76 1.00 0.58 0.72 -0.03 0.39 0.05 0.70
## S 0.33 0.35 0.35 0.49 0.57 0.58 1.00 0.64 -0.05 0.48 -0.06 0.74
## Na 0.08 0.07 0.07 0.50 0.61 0.72 0.64 1.00 -0.04 0.19 -0.09 0.62
## Fe -0.37 -0.17 -0.17 -0.28 -0.28 -0.03 -0.05 -0.04 1.00 0.01 0.48 -0.27
## Cu 0.01 0.50 0.50 0.32 0.57 0.39 0.48 0.19 0.01 1.00 0.19 0.63
## Mn -0.12 -0.10 -0.10 -0.05 -0.11 0.05 -0.06 -0.09 0.48 0.19 1.00 0.01
## Zn 0.26 0.26 0.26 0.46 0.80 0.70 0.74 0.62 -0.27 0.63 0.01 1.00
## B 0.40 0.54 0.54 0.63 0.64 0.50 0.84 0.70 -0.19 0.48 -0.13 0.72
## BS 0.10 0.16 0.16 0.40 0.73 0.64 0.36 0.30 -0.39 0.53 -0.09 0.56
## Alsat -0.10 -0.14 -0.14 -0.39 -0.71 -0.63 -0.35 -0.28 0.38 -0.54 0.04 -0.57
## S.Ca 0.13 0.35 0.35 0.25 0.82 0.49 0.29 0.24 -0.40 0.59 -0.17 0.59
## S.Mg -0.22 -0.51 -0.51 0.09 0.21 0.72 0.24 0.32 0.07 0.20 0.22 0.28
## S.K 0.17 0.33 0.33 0.73 -0.13 -0.07 0.16 -0.02 -0.09 0.12 0.14 -0.02
## S.Na 0.09 -0.13 -0.13 0.23 0.18 0.41 0.41 0.81 -0.03 -0.08 -0.05 0.31
## Ca.Mg 0.33 0.76 0.76 0.20 0.52 -0.13 0.16 0.03 -0.34 0.39 -0.19 0.30
## Mg.K -0.26 -0.57 -0.57 -0.38 0.36 0.60 0.11 0.35 0.01 0.02 -0.09 0.27
## Ca.K 0.00 0.00 0.00 -0.25 0.65 0.38 0.20 0.30 -0.23 0.22 -0.24 0.45
## Ca.Mg.K -0.07 -0.13 -0.13 -0.29 0.63 0.47 0.17 0.31 -0.19 0.20 -0.24 0.42
## Ca.B -0.30 -0.20 -0.20 -0.08 0.46 0.38 -0.19 0.04 -0.08 0.04 -0.05 0.13
## Fe.Mn -0.44 -0.23 -0.23 -0.39 -0.30 -0.12 -0.15 -0.02 0.79 -0.18 0.00 -0.32
## Fe.Zn -0.43 -0.36 -0.36 -0.49 -0.70 -0.42 -0.51 -0.39 0.66 -0.53 0.19 -0.79
## B BS Alsat S.Ca S.Mg S.K S.Na Ca.Mg Mg.K Ca.K Ca.Mg.K
## Rmax1 -0.11 -0.12 0.11 -0.13 -0.07 0.18 -0.07 0.07 -0.24 -0.21 -0.25
## Rmax2 -0.14 0.08 -0.10 0.01 0.05 0.11 0.01 0.02 -0.04 -0.01 -0.04
## Rmax3 0.01 0.34 -0.37 0.28 0.08 0.31 -0.21 0.27 -0.23 -0.07 -0.09
## ALT 0.04 -0.17 0.14 -0.20 -0.05 0.12 0.22 -0.02 -0.16 -0.06 -0.11
## TEM -0.37 0.02 -0.04 0.06 -0.09 -0.10 -0.34 0.04 0.04 0.01 0.05
## pH 0.21 0.80 -0.78 0.84 0.24 -0.11 -0.15 0.45 0.30 0.55 0.56
## EA -0.35 -0.95 0.95 -0.87 -0.46 -0.05 -0.03 -0.30 -0.32 -0.48 -0.49
## Al -0.35 -0.95 0.96 -0.86 -0.47 -0.06 -0.02 -0.29 -0.30 -0.45 -0.46
## EC 0.87 0.48 -0.47 0.52 0.09 0.10 0.40 0.44 0.09 0.35 0.29
## CEC 0.69 0.52 -0.49 0.56 0.18 -0.13 0.30 0.35 0.33 0.48 0.48
## OM 0.40 0.10 -0.10 0.13 -0.22 0.17 0.09 0.33 -0.26 0.00 -0.07
## N.tot 0.40 0.10 -0.10 0.13 -0.22 0.17 0.09 0.33 -0.26 0.00 -0.07
## N.dis 0.40 0.10 -0.10 0.13 -0.22 0.17 0.09 0.33 -0.26 0.00 -0.07
## P 0.54 0.16 -0.14 0.35 -0.51 0.33 -0.13 0.76 -0.57 0.00 -0.13
## Pav 0.54 0.16 -0.14 0.35 -0.51 0.33 -0.13 0.76 -0.57 0.00 -0.13
## K 0.63 0.40 -0.39 0.25 0.09 0.73 0.23 0.20 -0.38 -0.25 -0.29
## Ca 0.64 0.73 -0.71 0.82 0.21 -0.13 0.18 0.52 0.36 0.65 0.63
## Mg 0.50 0.64 -0.63 0.49 0.72 -0.07 0.41 -0.13 0.60 0.38 0.47
## S 0.84 0.36 -0.35 0.29 0.24 0.16 0.41 0.16 0.11 0.20 0.17
## Na 0.70 0.30 -0.28 0.24 0.32 -0.02 0.81 0.03 0.35 0.30 0.31
## Fe -0.19 -0.39 0.38 -0.40 0.07 -0.09 -0.03 -0.34 0.01 -0.23 -0.19
## Cu 0.48 0.53 -0.54 0.59 0.20 0.12 -0.08 0.39 0.02 0.22 0.20
## Mn -0.13 -0.09 0.04 -0.17 0.22 0.14 -0.05 -0.19 -0.09 -0.24 -0.24
## Zn 0.72 0.56 -0.57 0.59 0.28 -0.02 0.31 0.30 0.27 0.45 0.42
## B 1.00 0.36 -0.34 0.35 0.07 0.22 0.46 0.36 -0.02 0.19 0.15
## BS 0.36 1.00 -0.99 0.89 0.49 0.06 0.00 0.27 0.26 0.40 0.42
## Alsat -0.34 -0.99 1.00 -0.87 -0.49 -0.08 0.01 -0.26 -0.25 -0.38 -0.40
## S.Ca 0.35 0.89 -0.87 1.00 0.18 -0.16 -0.12 0.60 0.25 0.64 0.61
## S.Mg 0.07 0.49 -0.49 0.18 1.00 0.02 0.33 -0.60 0.56 0.03 0.18
## S.K 0.22 0.06 -0.08 -0.16 0.02 1.00 0.00 -0.07 -0.70 -0.70 -0.74
## S.Na 0.46 0.00 0.01 -0.12 0.33 0.00 1.00 -0.24 0.34 0.11 0.14
## Ca.Mg 0.36 0.27 -0.26 0.60 -0.60 -0.07 -0.24 1.00 -0.27 0.47 0.33
## Mg.K -0.02 0.26 -0.25 0.25 0.56 -0.70 0.34 -0.27 1.00 0.67 0.79
## Ca.K 0.19 0.40 -0.38 0.64 0.03 -0.70 0.11 0.47 0.67 1.00 0.98
## Ca.Mg.K 0.15 0.42 -0.40 0.61 0.18 -0.74 0.14 0.33 0.79 0.98 1.00
## Ca.B -0.25 0.39 -0.38 0.46 0.11 -0.42 -0.27 0.19 0.44 0.48 0.54
## Fe.Mn -0.23 -0.46 0.48 -0.41 -0.07 -0.29 0.02 -0.31 0.15 -0.05 0.00
## Fe.Zn -0.63 -0.56 0.56 -0.62 -0.10 -0.07 -0.18 -0.45 -0.10 -0.38 -0.33
## Ca.B Fe.Mn Fe.Zn
## Rmax1 -0.21 -0.06 0.07
## Rmax2 -0.02 0.01 0.06
## Rmax3 -0.04 -0.40 -0.07
## ALT -0.39 0.14 0.22
## TEM 0.50 -0.04 -0.04
## pH 0.65 -0.41 -0.51
## EA -0.39 0.54 0.62
## Al -0.39 0.55 0.61
## EC 0.05 -0.34 -0.76
## CEC 0.43 -0.15 -0.58
## OM -0.30 -0.44 -0.43
## N.tot -0.30 -0.44 -0.43
## N.dis -0.30 -0.44 -0.43
## P -0.20 -0.23 -0.36
## Pav -0.20 -0.23 -0.36
## K -0.08 -0.39 -0.49
## Ca 0.46 -0.30 -0.70
## Mg 0.38 -0.12 -0.42
## S -0.19 -0.15 -0.51
## Na 0.04 -0.02 -0.39
## Fe -0.08 0.79 0.66
## Cu 0.04 -0.18 -0.53
## Mn -0.05 0.00 0.19
## Zn 0.13 -0.32 -0.79
## B -0.25 -0.23 -0.63
## BS 0.39 -0.46 -0.56
## Alsat -0.38 0.48 0.56
## S.Ca 0.46 -0.41 -0.62
## S.Mg 0.11 -0.07 -0.10
## S.K -0.42 -0.29 -0.07
## S.Na -0.27 0.02 -0.18
## Ca.Mg 0.19 -0.31 -0.45
## Mg.K 0.44 0.15 -0.10
## Ca.K 0.48 -0.05 -0.38
## Ca.Mg.K 0.54 0.00 -0.33
## Ca.B 1.00 -0.10 -0.06
## Fe.Mn -0.10 1.00 0.61
## Fe.Zn -0.06 0.61 1.00resumen.Rmax=colMeans(abs(Rmax))
#
head(round(rbind(resumen.Rmed, resumen.Rmea, resumen.Rmax), digits = 2) )## Rmea1 Rmea2 Rmea3 ALT TEM pH EA Al EC CEC OM N.tot
## resumen.Rmed 0.45 0.54 0.41 0.41 0.34 0.33 0.60 0.63 0.44 0.61 0.63 0.43
## resumen.Rmea 0.22 0.18 0.23 0.20 0.17 0.39 0.41 0.41 0.43 0.40 0.26 0.26
## resumen.Rmax 0.19 0.13 0.21 0.19 0.16 0.39 0.42 0.41 0.43 0.41 0.25 0.25
## N.dis P Pav K Ca Mg S Na Fe Cu Mn Zn B
## resumen.Rmed 0.65 0.66 0.45 0.58 0.45 0.54 0.41 0.41 0.34 0.33 0.60 0.63 0.44
## resumen.Rmea 0.26 0.30 0.30 0.36 0.47 0.39 0.35 0.31 0.27 0.32 0.18 0.42 0.40
## resumen.Rmax 0.25 0.30 0.30 0.34 0.47 0.39 0.36 0.32 0.27 0.31 0.18 0.42 0.39
## BS Alsat S.Ca S.Mg S.K S.Na Ca.Mg Mg.K Ca.K Ca.Mg.K Ca.B Fe.Mn
## resumen.Rmed 0.61 0.63 0.43 0.65 0.66 0.45 0.58 0.45 0.54 0.41 0.41 0.34
## resumen.Rmea 0.41 0.40 0.43 0.26 0.23 0.21 0.33 0.32 0.32 0.35 0.29 0.29
## resumen.Rmax 0.41 0.40 0.43 0.26 0.22 0.21 0.33 0.31 0.32 0.34 0.29 0.29
## Fe.Zn
## resumen.Rmed 0.33
## resumen.Rmea 0.42
## resumen.Rmax 0.42Gráficos dispersión exploratorios - BD1 rendimiento
A <- data.frame(GR.CUN[1], GR.CUN[20:54])
# Obtener los nombres de las variables
variables <- colnames(A)
par(mfrow=c(1,3))
# Iterar a través de las columnas y crear los gráficos de dispersion
for (i in variables) {
# Crear un gráfico de caja para la columna i
plot(y = A$Rmed1, x = A[, i], main = i, ylab = "Rendimiento", xlab = i)
}
Correlograma - Todas las variables
A.cor <- cor( A, method = "spearman", use = "complete.obs")
# Variables que mas se correlacionan con otras
sort(colMeans(abs(A.cor)))## Mn Rmed1 TEM ALT S.K OM N.tot N.dis
## 0.1209878 0.1595547 0.1750694 0.1928242 0.2353007 0.2383485 0.2383485 0.2383485
## S.Na Ca.Mg Fe Cu S.Mg P Pav Ca.B
## 0.2662123 0.2916740 0.2969598 0.2987915 0.3022473 0.3034871 0.3034871 0.3140793
## Fe.Mn Mg.K Ca.K Ca.Mg.K K Na S Mg
## 0.3153580 0.3273464 0.3629852 0.3668755 0.3778980 0.3810856 0.3923306 0.4205345
## pH B Zn CEC S.Ca Al Alsat EA
## 0.4244952 0.4325956 0.4486457 0.4500061 0.4504486 0.4518154 0.4547798 0.4577166
## BS Fe.Zn EC Ca
## 0.4603993 0.4629116 0.4629762 0.5002974#
corrplot(A.cor, method = 'number', tl.cex = 0.5, number.cex = 0.4)
corrplot.mixed(A.cor, lower.col = "black", tl.cex = 0.5, number.cex = 0.4)
Correlograma - Variables más importantes
# Retiro variables que poco aportan al PCA.
names(A)## [1] "Rmed1" "ALT" "TEM" "pH" "EA" "Al" "EC"
## [8] "CEC" "OM" "N.tot" "N.dis" "P" "Pav" "K"
## [15] "Ca" "Mg" "S" "Na" "Fe" "Cu" "Mn"
## [22] "Zn" "B" "BS" "Alsat" "S.Ca" "S.Mg" "S.K"
## [29] "S.Na" "Ca.Mg" "Mg.K" "Ca.K" "Ca.Mg.K" "Ca.B" "Fe.Mn"
## [36] "Fe.Zn"B <- A[,-c(1,2,3,6,10,11,12,19,21,22,26,27,28,29,30,31,32,33,34)];B## pH EA EC CEC OM Pav
## BOGOTA 5.559103 1.48842672 0.4744674 12.685202 9.778200 128.18565
## CAJICA 5.316000 1.10836635 0.4607279 10.384645 8.094967 67.90187
## CARMEN DE CARUPA 5.062155 2.06117455 0.3139291 6.417927 12.586212 115.27883
## CHIA 5.918421 0.08288328 0.8380108 14.984885 7.087336 136.03607
## CHIPAQUE 5.771273 0.76391403 0.5017144 12.069595 6.234734 243.17715
## CHOACHI 5.151642 2.73579100 0.2094300 7.449557 7.019542 46.30217
## CHOCONTA 5.114296 1.57266228 0.3881219 7.185019 9.893829 119.02915
## COGUA 5.452424 1.13583515 0.3825780 8.582600 10.649070 110.17181
## COTA 5.014000 1.68571680 1.2191779 11.399998 8.561286 86.67598
## CUCUNUBA 5.555217 0.50128910 0.3394613 7.801094 6.126479 48.70733
## EL ROSAL 5.579444 0.47765807 0.6773613 12.920994 20.790809 67.66688
## FACATATIVA 5.761250 0.32234394 1.2174820 16.797049 13.113757 121.00601
## FOSCA 5.394667 0.86581904 0.3239092 6.237568 5.311936 256.24507
## FUNZA 5.596444 0.29146909 1.1630190 13.420920 11.586402 162.55062
## FUQUENE 5.097600 1.15825959 0.3985827 8.620903 8.077060 46.89706
## GRANADA 5.338798 0.60767044 1.0046892 9.646722 13.933129 109.24778
## GUACHETA 5.192105 0.97107427 0.4200620 8.585921 5.944584 44.54421
## GUASCA 5.509538 0.47896412 0.2686508 5.720285 10.018231 69.44832
## GUATAVITA 5.104896 2.05703561 0.2192190 5.516139 12.252528 62.43887
## GUTIERREZ 5.304146 1.36450159 0.3135042 7.407471 7.487618 196.51523
## JUNIN 5.245667 2.73902422 0.1621760 6.131229 8.237882 25.11463
## LA CALERA 5.470921 1.46741721 0.4327048 9.397545 11.516739 83.48464
## LENGUAZAQUE 4.933537 2.64988291 0.3815322 6.748871 7.177118 160.74707
## MACHETA 5.224000 3.36183015 0.1341186 7.422145 4.311126 32.29613
## MADRID 5.772353 0.34159585 1.4138369 17.466278 13.073099 229.74511
## MANTA 5.731957 0.54565071 0.2626550 10.327537 4.042129 96.06190
## PACHO 5.222105 1.68835696 0.2557147 4.981400 10.007553 21.33420
## PASCA 5.687200 0.44322089 0.6390085 10.783257 10.289681 242.61607
## SAN BERNARDO 5.237612 1.79245649 0.3380505 8.861645 6.340610 180.82794
## SAN CAYETANO 5.159459 2.01875624 0.2280057 6.290645 7.793149 32.82863
## SESQUILE 5.184222 2.03168523 0.3121670 6.733806 12.484487 98.62468
## SIBATE 5.218485 1.35931642 0.6461802 8.908929 17.425073 138.12108
## SIMIJACA 5.044405 1.79354525 0.9808050 11.942048 7.753307 64.74293
## SOACHA 5.767143 0.60816555 2.1790074 14.701346 8.185902 98.50659
## SOPO 5.533333 0.74140222 0.4027289 9.028509 5.764675 43.50582
## SUBACHOQUE 5.296275 1.11352096 0.4407986 7.670090 20.534207 82.05546
## SUESCA 5.359383 0.72004084 0.2968667 6.033355 6.763315 66.38734
## SUSA 5.174068 1.93916819 0.2559407 6.598838 11.852425 94.47381
## SUTATAUSA 5.585789 0.52922592 1.1215687 11.402752 5.880491 181.30399
## TABIO 5.374000 1.22141038 0.4345830 9.137376 10.278675 49.97889
## TAUSA 5.010575 3.25329783 0.4933287 7.664277 19.107391 231.03623
## TENJO 5.531667 0.65449398 0.8992575 11.832478 13.849393 230.37548
## UBAQUE 5.477308 0.67162119 0.4094572 7.157293 12.631683 253.55491
## UNE 5.901754 0.56291347 0.4648777 11.402627 6.630656 390.24844
## VILLAPINZON 5.264062 1.22915100 0.3099199 6.455969 11.234601 167.86534
## ZIPACON 5.543889 0.25472939 0.3463251 7.790373 13.774041 22.88857
## ZIPAQUIRA 5.420000 1.29654018 0.4560806 9.809536 13.956820 104.95754
## K Ca Mg S Na Cu
## BOGOTA 0.8372324 7.609837 1.7951393 10.507277 0.12319214 2.3304314
## CAJICA 0.6229017 6.828402 1.5094093 13.102876 0.30956664 2.6000000
## CARMEN DE CARUPA 0.6813839 2.703054 0.7181425 8.190727 0.08674512 1.5063362
## CHIA 1.3264158 10.384599 2.5320848 31.323270 0.27498420 3.5517887
## CHIPAQUE 0.6996663 8.785649 1.4130209 12.127141 0.10290401 4.2071258
## CHOACHI 0.4339747 2.978155 0.8782449 6.669945 0.07743290 1.3189923
## CHOCONTA 0.7592882 3.426372 0.9684186 8.993108 0.17511828 2.5644189
## COGUA 0.5279060 5.137754 1.4431673 15.964870 0.13715575 3.8724567
## COTA 1.3810860 5.196928 2.7790936 20.586415 0.35717320 4.3288000
## CUCUNUBA 0.8593057 4.451287 1.5760503 11.428284 0.12308086 2.7378376
## EL ROSAL 1.2002920 8.727308 2.1647076 17.861430 0.17910686 2.7410294
## FACATATIVA 1.1939526 11.861185 2.4470635 39.566737 0.25014833 2.4855832
## FOSCA 0.5307785 3.943650 0.7771380 8.819503 0.05958466 7.1078333
## FUNZA 0.6941090 8.937923 2.6222886 49.984800 0.40958406 3.7488916
## FUQUENE 0.4861748 4.639825 2.0289622 37.360956 0.13548251 3.6976594
## GRANADA 0.5508807 6.193626 1.8223795 46.500499 0.11860993 5.2264668
## GUACHETA 0.6005536 4.809975 1.9723353 16.462326 0.15544412 3.3075683
## GUASCA 0.4989278 3.367402 1.1208034 8.403813 0.09926487 2.6483480
## GUATAVITA 0.3936100 1.962767 0.7815846 7.676379 0.05321787 2.6560143
## GUTIERREZ 0.4522899 4.620345 0.7560211 8.590721 0.09168774 2.5071049
## JUNIN 0.2749652 2.092673 0.6138866 5.284479 0.07137625 1.5084053
## LA CALERA 0.5700326 5.960935 0.9350799 8.644076 0.09284027 1.9973696
## LENGUAZAQUE 0.5608173 2.382194 0.7358587 14.960065 0.10778068 1.7988788
## MACHETA 0.4614810 2.317661 1.1483899 2.926136 0.08478355 0.8968649
## MADRID 1.9763222 10.806383 3.7989980 80.936942 0.52650822 2.4339097
## MANTA 0.6681469 7.370387 1.3783723 6.815903 0.13806267 4.6780870
## PACHO 0.2527059 1.932011 0.6736951 7.971077 0.06530890 1.7057403
## PASCA 0.7231711 7.924545 1.1939198 14.293321 0.12028547 8.7542173
## SAN BERNARDO 0.3785195 5.036355 1.0463044 9.607265 0.09865306 3.2564606
## SAN CAYETANO 0.1826725 2.847023 0.7954409 5.257374 0.10541479 0.8345135
## SESQUILE 0.5565116 2.827452 0.6986989 9.848105 0.05951783 2.7688667
## SIBATE 0.6264456 5.325352 1.0166715 14.500787 0.18124861 4.2744133
## SIMIJACA 0.7921099 6.725803 2.1156070 76.475056 0.23171513 1.1451890
## SOACHA 0.8523538 9.432303 3.0191177 73.127804 0.77869144 3.2276655
## SOPO 0.5194166 5.295229 2.2712253 11.543862 0.19123587 1.9216120
## SUBACHOQUE 0.8484872 4.075166 1.1870670 13.485496 0.10165910 2.2254359
## SUESCA 0.7465279 3.135905 1.1041059 8.328884 0.08185894 2.9441468
## SUSA 0.6221803 2.962146 0.8112063 17.523490 0.08375177 2.0400633
## SUTATAUSA 0.6200010 8.159469 1.3669953 34.351701 0.50158465 2.3745789
## TABIO 0.8974949 4.789854 1.6674089 13.560286 0.15346145 2.8903118
## TAUSA 0.6459231 2.759601 0.6489134 14.953695 0.16554767 1.3604383
## TENJO 1.2233185 7.490951 1.8092387 25.496444 0.23576234 3.0664202
## UBAQUE 0.7351836 4.676258 0.8740120 11.937510 0.06705591 3.3484615
## UNE 0.7137648 8.240352 1.2919375 12.952596 0.10805196 7.4151957
## VILLAPINZON 0.6554738 3.531951 0.7974836 14.355854 0.07051152 2.9893624
## ZIPACON 0.4305246 5.514894 1.4249624 9.759495 0.08402422 2.0710813
## ZIPAQUIRA 0.5026438 6.082949 1.4028510 13.693114 0.14889878 3.4091941
## B BS Alsat Fe.Mn Fe.Zn
## BOGOTA 0.3275570 80.05188 14.7452755 76.11250 159.48352
## CAJICA 0.3558244 80.67935 15.4288655 99.79597 224.82995
## CARMEN DE CARUPA 0.2525678 66.48783 26.3535645 101.88795 348.51124
## CHIA 0.4177221 99.31008 0.4654782 82.06782 44.59346
## CHIPAQUE 0.3241120 89.93421 7.9391507 106.15040 192.95755
## CHOACHI 0.1272001 59.04389 31.9970562 127.84451 469.15464
## CHOCONTA 0.2669020 76.23524 19.5441738 64.50734 279.89002
## COGUA 0.2983707 85.37701 12.1882056 127.99761 231.20179
## COTA 0.4398446 72.98395 19.4199499 137.53683 58.32030
## CUCUNUBA 0.2615827 92.33921 5.3468453 124.07392 283.95300
## EL ROSAL 0.3866722 92.95790 4.5533721 48.45113 31.38065
## FACATATIVA 0.4082495 95.29772 2.8526592 45.69982 39.99380
## FOSCA 0.2029506 83.60708 12.1007762 73.69159 121.38584
## FUNZA 0.5725170 96.31206 1.9473541 77.67510 37.83787
## FUQUENE 0.2827392 85.86268 10.1669683 119.74120 291.30422
## GRANADA 0.5594030 89.49959 7.1625878 57.03536 58.06217
## GUACHETA 0.2952467 87.07609 9.5278237 124.07656 235.46899
## GUASCA 0.2513500 90.86439 5.9007369 82.78848 271.42576
## GUATAVITA 0.1513918 62.78067 23.6880716 99.78590 406.21095
## GUTIERREZ 0.3029657 78.26410 17.2662794 62.19801 181.54670
## JUNIN 0.1574059 48.76613 41.3027097 120.11062 543.41485
## LA CALERA 0.2263380 75.56990 19.0462365 106.75522 252.41292
## LENGUAZAQUE 0.3167189 58.54142 34.4287355 194.46502 604.75564
## MACHETA 0.1479174 53.18053 40.2219427 105.21371 770.00062
## MADRID 0.5474817 94.83883 3.9872994 54.92473 30.35978
## MANTA 0.3111185 91.89935 5.9956611 93.91041 269.76120
## PACHO 0.1343537 63.66178 27.1144823 89.83600 155.79598
## PASCA 0.3379247 91.86695 5.9067966 81.70473 53.13072
## SAN BERNARDO 0.2715315 76.61698 18.8071387 90.54170 174.79516
## SAN CAYETANO 0.1695861 61.76864 32.0957150 151.18552 266.55311
## SESQUILE 0.1569658 65.19560 22.4494314 91.64701 328.82178
## SIBATE 0.3672236 82.59789 12.1732038 98.58501 104.90479
## SIMIJACA 0.3752028 81.78400 14.8497602 146.26137 291.09793
## SOACHA 0.4015325 83.76125 12.2199911 59.29421 60.45141
## SOPO 0.2216387 86.61411 9.8155378 130.78752 387.17109
## SUBACHOQUE 0.3047218 81.68475 14.0101016 74.09232 119.03846
## SUESCA 0.2140021 88.01743 8.2137716 66.46467 266.70547
## SUSA 0.2305126 69.28075 23.8886218 124.09790 355.08541
## SUTATAUSA 0.3312990 93.87815 4.6078133 65.38077 113.22886
## TABIO 0.2634338 85.63516 11.4733695 79.87481 175.23262
## TAUSA 0.3908919 55.09451 35.6446294 81.64860 164.86202
## TENJO 0.4902641 92.25186 5.8883004 71.50533 112.37959
## UBAQUE 0.2824154 88.07622 8.9062406 62.11595 70.69824
## UNE 0.3007872 92.26753 6.1445993 92.23563 85.98830
## VILLAPINZON 0.3706109 77.13117 18.5653449 100.98522 248.48556
## ZIPACON 0.1585402 93.37620 4.1878722 48.52643 74.34853
## ZIPAQUIRA 0.3252102 83.87497 12.9462299 124.80835 229.29622names(B) <- c("pH","EA","EC","CEC","OM","Pav","K","Ca","Mg","S","Na","Cu","B","BS","Alsat","Fe/Mn","Fe/Zn")
# [,-c(1,2,4,6,10,11,12,13,20,22,29,30,32,35,36,37)] -> 40.80+15.09 = 55.89 (Pedro)
# [,-c(1,2,3,4,6,7,11,12,13,20,21,22,28,29,30,31,32,33)] -> 47.14+15.57 = 62.71 (Mariam)
# [,-c(1,2,4,6,7,10,11,13,20,21,22,28,30,31,32,33,36,37)] -> 43.09+14.27 = 57.36 (Geraldine y Jimer)
# [,-c(1,2,3,6,10,11,12,19,21,22,26,27,28,29,30,31,32,33,34)] -> 52.50+14.80 = 67.3 (Pedro 2)
B.cor <- cor( B, method = "spearman", use = "complete.obs")
round(B.cor, digits = 2)## pH EA EC CEC OM Pav K Ca Mg S Na Cu
## pH 1.00 -0.83 0.43 0.61 -0.07 0.32 0.40 0.77 0.52 0.22 0.33 0.40
## EA -0.83 1.00 -0.52 -0.55 -0.05 -0.24 -0.45 -0.74 -0.64 -0.43 -0.43 -0.51
## EC 0.43 -0.52 1.00 0.86 0.29 0.44 0.68 0.82 0.71 0.83 0.80 0.36
## CEC 0.61 -0.55 0.86 1.00 0.05 0.31 0.61 0.94 0.82 0.68 0.80 0.32
## OM -0.07 -0.05 0.29 0.05 1.00 0.13 0.20 0.05 -0.02 0.30 0.04 -0.03
## Pav 0.32 -0.24 0.44 0.31 0.13 1.00 0.36 0.39 -0.04 0.34 0.16 0.47
## K 0.40 -0.45 0.68 0.61 0.20 0.36 1.00 0.57 0.56 0.55 0.57 0.27
## Ca 0.77 -0.74 0.82 0.94 0.05 0.39 0.57 1.00 0.79 0.63 0.71 0.46
## Mg 0.52 -0.64 0.71 0.82 -0.02 -0.04 0.56 0.79 1.00 0.67 0.77 0.34
## S 0.22 -0.43 0.83 0.68 0.30 0.34 0.55 0.63 0.67 1.00 0.71 0.32
## Na 0.33 -0.43 0.80 0.80 0.04 0.16 0.57 0.71 0.77 0.71 1.00 0.18
## Cu 0.40 -0.51 0.36 0.32 -0.03 0.47 0.27 0.46 0.34 0.32 0.18 1.00
## B 0.34 -0.45 0.86 0.76 0.30 0.54 0.65 0.72 0.60 0.82 0.76 0.40
## BS 0.80 -0.96 0.56 0.62 0.03 0.28 0.49 0.78 0.68 0.50 0.49 0.52
## Alsat -0.79 0.96 -0.55 -0.59 -0.05 -0.27 -0.48 -0.76 -0.67 -0.48 -0.46 -0.52
## Fe/Mn -0.48 0.55 -0.36 -0.24 -0.36 -0.36 -0.35 -0.37 -0.17 -0.19 -0.12 -0.19
## Fe/Zn -0.58 0.67 -0.75 -0.62 -0.38 -0.47 -0.51 -0.71 -0.49 -0.55 -0.47 -0.48
## B BS Alsat Fe/Mn Fe/Zn
## pH 0.34 0.80 -0.79 -0.48 -0.58
## EA -0.45 -0.96 0.96 0.55 0.67
## EC 0.86 0.56 -0.55 -0.36 -0.75
## CEC 0.76 0.62 -0.59 -0.24 -0.62
## OM 0.30 0.03 -0.05 -0.36 -0.38
## Pav 0.54 0.28 -0.27 -0.36 -0.47
## K 0.65 0.49 -0.48 -0.35 -0.51
## Ca 0.72 0.78 -0.76 -0.37 -0.71
## Mg 0.60 0.68 -0.67 -0.17 -0.49
## S 0.82 0.50 -0.48 -0.19 -0.55
## Na 0.76 0.49 -0.46 -0.12 -0.47
## Cu 0.40 0.52 -0.52 -0.19 -0.48
## B 1.00 0.48 -0.47 -0.30 -0.67
## BS 0.48 1.00 -0.99 -0.48 -0.62
## Alsat -0.47 -0.99 1.00 0.51 0.63
## Fe/Mn -0.30 -0.48 0.51 1.00 0.67
## Fe/Zn -0.67 -0.62 0.63 0.67 1.00# Variables que mas se correlacionan con otras
sort(colMeans(abs(B.cor)))## OM Pav Fe/Mn Cu K Na pH S
## 0.1955910 0.3599132 0.3942088 0.3969500 0.5107539 0.5164608 0.5220452 0.5414853
## Mg EA B Alsat Fe/Zn BS CEC EC
## 0.5584494 0.5874327 0.5950305 0.5999483 0.6028459 0.6045805 0.6102057 0.6356791
## Ca
## 0.6583637# Correlograma
corrplot(B.cor, method = 'number', tl.cex = 0.7, number.cex = 0.7)
corrplot.mixed(B.cor, lower.col = "black", tl.cex = 0.7, number.cex = 0.7)
corrplot.mixed(B.cor, tl.cex = 0.75, number.cex = 0.75)
Correlograma variables PCA
B %>% ggpairs(upper = list(continuous = wrap("cor", method = "spearman")),
title=NULL,
axisLabels = "none")
3. A N A L I S I S D E C O M P O N E N T E S P R I N C I P A L E S ( P C A )
PCA
names(B)## [1] "pH" "EA" "EC" "CEC" "OM" "Pav" "K" "Ca" "Mg"
## [10] "S" "Na" "Cu" "B" "BS" "Alsat" "Fe/Mn" "Fe/Zn"res.pca <- PCA(B, graph = TRUE)
Extract eigenvalues/variances
get_eig(res.pca)## eigenvalue variance.percent cumulative.variance.percent
## Dim.1 8.915122903 52.44189943 52.44190
## Dim.2 2.519799554 14.82235032 67.26425
## Dim.3 1.466089258 8.62405446 75.88830
## Dim.4 1.181438451 6.94963795 82.83794
## Dim.5 0.704201824 4.14236367 86.98031
## Dim.6 0.646763495 3.80449114 90.78480
## Dim.7 0.364586493 2.14462643 92.92942
## Dim.8 0.329931352 1.94077266 94.87020
## Dim.9 0.247262942 1.45448790 96.32468
## Dim.10 0.201001040 1.18235906 97.50704
## Dim.11 0.182388333 1.07287255 98.57992
## Dim.12 0.093995889 0.55291700 99.13283
## Dim.13 0.076491209 0.44994829 99.58278
## Dim.14 0.041969143 0.24687731 99.82966
## Dim.15 0.022175955 0.13044680 99.96010
## Dim.16 0.005908239 0.03475435 99.99486
## Dim.17 0.000873919 0.00514070 100.00000Visualize eigenvalues/variances
factoextra::fviz_screeplot(res.pca, addlabels = TRUE, ylim = c(0, 80), title=NULL)
Extract the results for variables
var <- get_pca_var(res.pca); var## Principal Component Analysis Results for variables
## ===================================================
## Name Description
## 1 "$coord" "Coordinates for the variables"
## 2 "$cor" "Correlations between variables and dimensions"
## 3 "$cos2" "Cos2 for the variables"
## 4 "$contrib" "contributions of the variables"Coordinates of variables
var$coord## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
## pH 0.7419105 -0.36736588 0.25719364 -0.08160693 -0.37213143
## EA -0.8046118 0.44023356 -0.15775541 0.27545817 -0.13255914
## EC 0.8056863 0.46739049 -0.04420906 0.11011718 -0.03403508
## CEC 0.8839086 0.26979692 0.09467178 0.07830395 -0.15503749
## OM 0.1469020 0.03634646 -0.90217765 -0.20331487 0.12910589
## Pav 0.3917497 -0.37514333 -0.19195279 0.74463753 -0.17290012
## K 0.7357174 0.25073089 -0.16314165 0.02685990 -0.01057288
## Ca 0.9265724 0.01155272 0.12819764 0.03425314 -0.15389105
## Mg 0.8151191 0.41167057 0.20194558 -0.14353063 0.11518739
## S 0.6946738 0.54683874 0.02093514 0.06426773 0.11977685
## Na 0.6985322 0.53764569 0.11362089 0.04541197 -0.14434488
## Cu 0.4175670 -0.60851571 0.10160619 0.46560213 0.23478927
## B 0.7858379 0.23715347 -0.23620779 0.25389241 0.26027618
## BS 0.8347073 -0.40006822 0.15647160 -0.21948103 0.18530463
## Alsat -0.8267000 0.41416468 -0.12414747 0.24414820 -0.20195772
## Fe/Mn -0.5350788 0.32061470 0.43698172 0.25787541 0.44505549
## Fe/Zn -0.7865636 0.27614661 0.34528191 0.02464159 -0.08987313Contribution of variables
var$contrib## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
## pH 6.174128 5.355890037 4.51190590 0.56369337 19.66507293
## EA 7.261820 7.691309706 1.69749344 6.42244257 2.49529665
## EC 7.281227 8.669494071 0.13330980 1.02635853 0.16449643
## CEC 8.763697 2.888736922 0.61133704 0.51898675 3.41331469
## OM 0.242063 0.052427383 55.51670945 3.49886485 2.36698210
## Pav 1.721432 5.585068065 2.51320815 46.93304646 4.24515419
## K 6.071481 2.494880119 1.81538722 0.06106575 0.01587412
## Ca 9.630113 0.005296661 1.12098457 0.09930922 3.36302101
## Mg 7.452721 6.725640352 2.78168718 1.74372525 1.88413802
## S 5.412956 11.867317347 0.02989449 0.34960273 2.03727007
## Na 5.473253 11.471662156 0.88055400 0.17455394 2.95873197
## Cu 1.955803 14.695270943 0.70417386 18.34927135 7.82815387
## B 6.926895 2.231993810 3.80564272 5.45617541 9.61992548
## BS 7.815218 6.351877455 1.66997754 4.07739596 4.87613134
## Alsat 7.665995 6.807381881 1.05127254 5.04540388 5.79193611
## Fe/Mn 3.211502 4.079442903 13.02465188 5.62870833 28.12750311
## Fe/Zn 6.939695 3.026310189 8.13181022 0.05139563 1.14699791Graph of variables: default plot
## Control variable colors using their contributions
fviz_pca_var(res.pca, col.var="contrib",
gradient.cols = c("darkred", "#E7B800", "darkgreen"),
repel = TRUE # Avoid text overlapping
)
Contributions of variables to PC1
fviz_contrib(res.pca, choice = "var", axes = 1, top = 20, title=F)
Contributions of variables to PC2
fviz_contrib(res.pca, choice = "var", axes = 2, top = 20, title=F)
Contributions of variables to PC3
fviz_contrib(res.pca, choice = "var", axes = 3, top = 20, title=F)
Extract the results for individuals
ind <- get_pca_ind(res.pca); ind## Principal Component Analysis Results for individuals
## ===================================================
## Name Description
## 1 "$coord" "Coordinates for the individuals"
## 2 "$cos2" "Cos2 for the individuals"
## 3 "$contrib" "contributions of the individuals"Coordinates of individuals
ind$coord## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
## BOGOTA 0.8963707 -0.02386459 -0.07983630 -0.18879900 -1.030030771
## CAJICA 0.2405183 0.59034094 0.62866413 -0.31659591 -0.003737490
## CARMEN DE CARUPA -3.0476672 0.62852122 -1.17586604 0.11477304 -0.084946725
## CHIA 4.9735875 -0.19867426 1.26827272 -0.40311595 -0.480989119
## CHIPAQUE 1.6556870 -1.65922762 1.32155623 1.18026056 -0.637140313
## CHOACHI -4.4512718 1.28180138 0.75961933 0.04648728 -0.731014841
## CHOCONTA -1.4795314 -0.07477085 -0.79147747 -0.13881540 -0.400419121
## COGUA -0.3315197 -0.60551593 0.56807578 0.06962268 0.900793620
## COTA 1.5883740 2.73037132 0.02187361 1.27901605 1.657377210
## CUCUNUBA -0.0637253 -0.79893171 1.75439371 -1.20132950 0.793593221
## EL ROSAL 3.3009806 -0.24209657 -2.23967595 -1.89701918 -0.043093714
## FACATATIVA 5.2865932 0.63991504 -0.43817044 -0.98440391 -1.111805799
## FOSCA -0.4950377 -3.21253303 0.44483005 1.77164416 0.035337218
## FUNZA 4.9094861 0.88803141 -0.04707146 0.26925492 0.662845174
## FUQUENE -0.4544724 0.47938935 1.08063272 -0.47324593 1.700056506
## GRANADA 2.6262992 -0.28570878 -1.39358859 0.18408520 1.534371766
## GUACHETA -0.3553489 0.06607123 1.45668973 -0.60254469 1.392693160
## GUASCA -0.7537495 -1.66757300 0.37109032 -1.57064792 0.340308034
## GUATAVITA -3.8364866 0.04791734 -0.59225462 -0.35359724 -0.018078681
## GUTIERREZ -1.1025289 -1.19978721 -0.44514252 0.46109907 -0.803830467
## JUNIN -5.3613499 1.40065131 0.39088343 0.18390813 -1.147889902
## LA CALERA -1.2183203 -0.14046941 -0.04632304 -0.52119443 -0.553469015
## LENGUAZAQUE -4.4597668 2.23931809 0.89220746 2.28548795 1.053683545
## MACHETA -5.3979361 2.09238180 1.36848508 0.21199181 -1.879577963
## MADRID 7.8992104 3.10358074 -0.65472634 0.42656903 -0.674204890
## MANTA 0.9108057 -1.81143774 2.05364986 -0.22417929 -0.219761790
## PACHO -3.5706459 -0.22054401 -0.52108890 -1.09108929 -0.409987834
## PASCA 2.4888318 -3.11255953 0.13157944 1.93737836 0.414218345
## SAN BERNARDO -1.3008572 -0.72053952 0.21565051 0.93795502 -0.369159807
## SAN CAYETANO -4.2101777 0.96013737 0.65994718 -0.29081509 0.106230304
## SESQUILE -2.9796176 -0.12975671 -0.86160256 -0.07846817 -0.379506844
## SIBATE 0.2196304 -0.44867550 -1.81450226 0.36144077 1.050397330
## SIMIJACA 0.8448797 3.58849804 0.99114938 0.20501015 1.474421889
## SOACHA 5.9507090 3.69041534 0.94458021 0.11850904 -1.450305417
## SOPO -0.4159880 0.21276016 2.37877094 -1.27878455 0.289914400
## SUBACHOQUE -0.2845036 -0.35651333 -2.63685317 -1.33081256 0.500143566
## SUESCA -0.9376978 -1.46607055 0.53016309 -1.32566661 0.009587139
## SUSA -2.8454173 0.55780440 -0.41691101 0.06941793 0.352194521
## SUTATAUSA 2.9918804 0.29949164 0.94428191 0.03272111 -1.217614545
## TABIO 0.0973267 -0.20850608 -0.06387439 -1.11368531 0.100692779
## TAUSA -2.5568567 1.60545935 -3.82070783 1.78613135 -0.694943684
## TENJO 3.3040450 -0.04843465 -1.25373812 0.56863390 0.020030402
## UBAQUE 0.4164768 -2.25268193 -1.32196480 0.29781642 -0.327547664
## UNE 2.5119098 -3.47096217 0.92686699 2.94771206 -0.716860545
## VILLAPINZON -1.3071825 -0.63222177 -0.80298261 0.61830409 0.569545806
## ZIPACON 0.2531064 -1.88988410 -0.52619312 -2.91868135 -0.347997643
## ZIPAQUIRA -0.1490518 -0.22491690 -0.15936228 -0.06173881 0.775478648Graph of individuals
# 1. Use repel = TRUE to avoid overplotting
# 2. Control automatically the color of individuals using the cos2
# cos2 = the quality of the individuals on the factor map
# Use points only
# 3. Use gradient color
fviz_pca_ind(res.pca, col.ind = "contrib",
gradient.cols = c("darkred", "#E7B800", "darkgreen"),
title=NULL,
repel = TRUE) # Avoid text overlapping (slow if many points)
Biplot of individuals and variables
fviz_pca_biplot(res.pca,
repel = TRUE,
labelsize=3,
col.ind="#666666",
col.var="cos2",
gradient.cols = c("darkred", "#E7B800", "darkgreen"),
title = NULL,
select.ind = list(contrib = 37) )+
theme_classic()
4. A G R U P A M I E N T O
Analisis de componentes principales y agrupamiento usando el IDPM
GR.CUN$IDPM <- as.factor(GR.CUN$IDPM)
# Visualize
# Use habillage to specify groups for coloring
fviz_pca_ind(res.pca,
label = "none", # hide individual labels
habillage = GR.CUN$IDPM, # color by groups
palette = c("blue2", "yellow3", "red3"),
addEllipses = TRUE # Concentration ellipses
)
Definimos paleta de color
color.pal <- c("#00AFBB","#e79e00","#9E1FDE", "#9D9FDE","#FC4E07")
color.pal2 <- c("#9E1FDE","#00AFBB" ,"#e79e00", "#9D9FDE","#FC4E07")Agrupamiento por por HCPC
set.seed(123)
# Fuente: Kassambara, 2017
# Compute PCA with ncp = 3
res.pca.HCPC <- PCA(B, ncp = 3, graph = FALSE)
# Compute hierarchical clustering on principal components
res.hcpc <- HCPC(res.pca.HCPC, nb.clust = -1, graph = FALSE)
#nb.clust: un número entero que especifica el número de clústeres. Los valores posibles son:
#0: el árbol se corta al nivel en el que el usuario hace clic
#-1: el árbol se corta automáticamente al nivel sugerido
#Any positive integer: el árbol se corta con racimos nb.clusters
# Dendogram
fviz_dend(res.hcpc,
cex = 0.6, # Label size
palette = color.pal2, # Color palette see ?ggpubr::ggpar
rect = TRUE,
rect_border = "#666666",
labels_track_height = 0.7 # Augment the room for labels
)
# Generacion dataframe con clusters
Gr.clus.hcpc <- data.frame(MUN = rownames(res.hcpc$data.clust),
Grupo.hcpc = res.hcpc$data.clust[, 18])Determinacion del numero de clusters, método jerarquico
es <- as.matrix(scale(B))
# "euclidean","maximum","manhattan","canberra","binary","minkowski","pearson","spearman","kendall"
dist.method <- "kendall"
# kmeans, cluster::pam, cluster::clara, cluster::fanny, hcut
# Grafico de codo
#fviz_nbclust(x=es, FUNcluster=hcut, method="wss", diss=get_dist(es, method=dist.method) )
# Gap-stat
#fviz_nbclust(x=es, FUNcluster=hcut, method="gap_stat", diss=get_dist(es, method=dist.method) )
# silhouette
#fviz_nbclust(x=es, FUNcluster=hcut, method="silhouette", diss=get_dist(es, method=dist.method))# Instala e importa los paquetes necesarios
library(factoextra)
library(cluster) # Necesario para hcut
# Escala los datos
es <- scale(B)
# Método de distancia
dist.method <- "kendall"
# Obtiene la matriz de distancia
dist_matrix <- get_dist(es, method = dist.method)
# Función de corte (hcut)
# hierarchical_clusters <- hcut(dist_matrix, k = 2:10, hc_method = "complete")
# Grafico de codo
#fviz_nbclust(x = es, FUNcluster = hcut, method = "wss", diss = dist_matrix)
# Gap-stat
#fviz_nbclust(x = es, FUNcluster = hcut, method = "gap_stat", diss = dist_matrix)
# Silhouette
#fviz_nbclust(x = es, FUNcluster = hcut, method = "silhouette", diss = dist_matrix)Determinacion del numero de clusters, método k-means
# Grafico de codo
fviz_nbclust(x=es, FUNcluster=kmeans, method="wss", diss=get_dist(es, method=dist.method) )
# Gap-stat
fviz_nbclust(x=es, FUNcluster=kmeans, method="gap_stat", diss=get_dist(es, method=dist.method) )
# silhouette
fviz_nbclust(x=es, FUNcluster=kmeans, method="silhouette", diss=get_dist(es, method=dist.method) )
Agrupamiento jerarquico
#h.res <- hcut(es, k = 4, stand = TRUE)
# Visualize
#fviz_dend(h.res, rect = TRUE, cex = 0.5, k_colors = color.pal)Agrupamiento k-means
# 2. Compute k-means
set.seed(123)
km.res <- kmeans(es, 3, nstart = 10)
# 3. Visualize
library("factoextra")
fviz_cluster(km.res, data = es, palette = c("#00AFBB","#e79e00" ,"#9E1FDE"), ggtheme = theme_classic(),
main = "Partitioning Clustering Plot", repel = TRUE)
km_clusters <- eclust(x = es, FUNcluster = "hclust", k = 3, seed = 123,
hc_metric = "euclidean", nstart = 50, graph = FALSE)
fviz_silhouette(sil.obj = km_clusters, print.summary = F, palette = c("#e79e00","#00AFBB","#9E1FDE","#2A1FDB" ),
ggtheme = theme_classic() ) 
Generacion de dataframes a partir de los cluster
# Generacion dataframe con clusters
Gr.clus <- data.frame(Gr.km = km.res$cluster,
Gr.hcpc = Gr.clus.hcpc$Grupo.hcpc) %>%
mutate(., Gr.km=as.factor(Gr.km), Gr.hcpc=as.factor(Gr.hcpc) ) %>%
tibble::rownames_to_column(., "MUN")
head(Gr.clus, n = 10)## MUN Gr.km Gr.hcpc
## 1 BOGOTA 2 2
## 2 CAJICA 2 2
## 3 CARMEN DE CARUPA 1 1
## 4 CHIA 3 3
## 5 CHIPAQUE 2 2
## 6 CHOACHI 1 1
## 7 CHOCONTA 2 2
## 8 COGUA 2 2
## 9 COTA 3 3
## 10 CUCUNUBA 2 2# Dataframe sin agrupar por municipios con grupos
DF.cluster.tot <- RAS.CUN %>% left_join(., y=Gr.clus) ## Joining with `by = join_by(MUN)`head(DF.cluster.tot, n = 10)## # A tibble: 10 × 59
## MUN CUL Acos1 Acos2 Acos3 ALT Asem1 Asem2 Asem3 IDPM Prod1 Prod2
## <chr> <fct> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 FUNZA Uchu… 408. 586. NA 2548 492. 651. NA 1 8161. 11779.
## 2 BOGOTA Papa… NA NA NA 2625 NA NA NA 3 NA NA
## 3 CHIA Frut… 118. 195. NA 2600 125. 189. NA 2 2035 3465.
## 4 CARMEN D… Papa… 855. 2083. 69.4 2600 945. 2276. 77.6 1 15280. 39427
## 5 CARMEN D… Papa… 855. 2083. 69.4 2600 945. 2276. 77.6 1 15280. 39427
## 6 CARMEN D… Papa… 855. 2083. 69.4 2600 945. 2276. 77.6 1 15280. 39427
## 7 CARMEN D… Papa… 855. 2083. 69.4 2600 945. 2276. 77.6 1 15280. 39427
## 8 CARMEN D… Arve… 855. 2083. 69.4 2600 945. 2276. 77.6 1 15280. 39427
## 9 CARMEN D… Past… 855. 2083. 69.4 2600 945. 2276. 77.6 1 15280. 39427
## 10 CARMEN D… Past… 855. 2083. 69.4 2600 945. 2276. 77.6 1 15280. 39427
## # ℹ 47 more variables: Prod3 <dbl>, Rmax1 <dbl>, Rmax2 <dbl>, Rmax3 <dbl>,
## # Rmea1 <dbl>, Rmea2 <dbl>, Rmea3 <dbl>, Rmed1 <dbl>, Rmed2 <dbl>,
## # Rmed3 <dbl>, TEM <dbl>, pH <dbl>, CEC <dbl>, EC <dbl>, OM <dbl>, P <dbl>,
## # K <dbl>, Ca <dbl>, S <dbl>, Mg <dbl>, Na <dbl>, Fe <dbl>, Cu <dbl>,
## # Mn <dbl>, Zn <dbl>, B <dbl>, EA <dbl>, Al <dbl>, Pav <dbl>, N.tot <dbl>,
## # N.dis <dbl>, BS <dbl>, Alsat <dbl>, S.Ca <dbl>, S.Mg <dbl>, S.K <dbl>,
## # S.Na <dbl>, Ca_Mg <dbl>, Mg_K <dbl>, Ca_K <dbl>, Ca.Mg_K <dbl>, …# Dataframe agrupado por municipios con grupos
A1 <- tibble::rownames_to_column(A, "MUN")
A2 <- tibble::rownames_to_column(GR.CUN[2:3], "MUN")
DF.cluster.tot <- A1 %>% full_join(., A2, by = "MUN") %>% left_join(., y=Gr.clus)## Joining with `by = join_by(MUN)`head(DF.cluster.tot, n=10)## MUN Rmed1 ALT TEM pH EA Al EC
## 1 BOGOTA NA 2625 13.1 5.559103 1.48842672 1.14813493 0.4744674
## 2 CAJICA 20.00 2558 14.0 5.316000 1.10836635 0.87816355 0.4607279
## 3 CARMEN DE CARUPA 15.00 2600 12.0 5.062155 2.06117455 1.63019473 0.3139291
## 4 CHIA 17.00 2600 14.0 5.918421 0.08288328 0.05627741 0.8380108
## 5 CHIPAQUE 14.50 2400 13.0 5.771273 0.76391403 0.61289557 0.5017144
## 6 CHOACHI 15.00 1923 18.0 5.151642 2.73579100 2.20641876 0.2094300
## 7 CHOCONTA 20.00 2689 10.0 5.114296 1.57266228 1.31436399 0.3881219
## 8 COGUA 19.67 2600 14.0 5.452424 1.13583515 0.94955086 0.3825780
## 9 COTA 23.00 2566 14.0 5.014000 1.68571680 1.16897724 1.2191779
## 10 CUCUNUBA 15.24 2590 14.0 5.555217 0.50128910 0.36834854 0.3394613
## CEC OM N.tot N.dis P Pav K
## 1 12.685202 9.778200 0.4889100 73.33650 55.97627 128.18565 0.8372324
## 2 10.384645 8.094967 0.4047484 60.71225 29.65147 67.90187 0.6229017
## 3 6.417927 12.586212 0.6293106 94.39659 50.34010 115.27883 0.6813839
## 4 14.984885 7.087336 0.3543668 53.15502 59.40440 136.03607 1.3264158
## 5 12.069595 6.234734 0.3117367 46.76050 106.19090 243.17715 0.6996663
## 6 7.449557 7.019542 0.3509771 52.64657 20.21929 46.30217 0.4339747
## 7 7.185019 9.893829 0.4946914 74.20372 51.97780 119.02915 0.7592882
## 8 8.582600 10.649070 0.5324535 79.86803 48.10996 110.17181 0.5279060
## 9 11.399998 8.561286 0.4280643 64.20965 37.84977 86.67598 1.3810860
## 10 7.801094 6.126479 0.3063240 45.94859 21.26958 48.70733 0.8593057
## Ca Mg S Na Fe Cu Mn
## 1 7.609837 1.7951393 10.507277 0.12319214 569.1982 2.330431 7.651463
## 2 6.828402 1.5094093 13.102876 0.30956664 668.5700 2.600000 7.291400
## 3 2.703054 0.7181425 8.190727 0.08674512 532.3260 1.506336 5.883974
## 4 10.384599 2.5320848 31.323270 0.27498420 488.5985 3.551789 6.985933
## 5 8.785649 1.4130209 12.127141 0.10290401 669.0530 4.207126 7.183408
## 6 2.978155 0.8782449 6.669945 0.07743290 793.1288 1.318992 6.973152
## 7 3.426372 0.9684186 8.993108 0.17511828 520.5033 2.564419 9.027619
## 8 5.137754 1.4431673 15.964870 0.13715575 763.5126 3.872457 6.770941
## 9 5.196928 2.7790936 20.586415 0.35717320 795.6744 4.328800 10.535400
## 10 4.451287 1.5760503 11.428284 0.12308086 624.6912 2.737838 5.479770
## Zn B BS Alsat S.Ca S.Mg S.K S.Na
## 1 5.463535 0.3275570 80.05188 14.7452755 55.16836 13.56820 7.243108 1.320546
## 2 6.398600 0.3558244 80.67935 15.4288655 57.30971 13.82068 6.534641 3.014313
## 3 1.920793 0.2525678 66.48783 26.3535645 41.91207 11.13792 10.871190 1.392300
## 4 11.472134 0.4177221 99.31008 0.4654782 68.81959 16.91423 8.629590 1.929716
## 5 5.334073 0.3241120 89.93421 7.9391507 69.18755 12.30076 6.428829 0.959038
## 6 3.039539 0.1272001 59.04389 31.9970562 39.42074 11.50512 5.908447 1.127415
## 7 2.880188 0.2669020 76.23524 19.5441738 47.32400 13.26507 10.551657 2.568885
## 8 5.599856 0.2983707 85.37701 12.1882056 59.12271 17.29694 6.084893 1.773415
## 9 28.979400 0.4398446 72.98395 19.4199499 39.26746 19.91903 11.360041 2.437411
## 10 2.904986 0.2615827 92.33921 5.3468453 56.74571 20.73359 11.023580 1.583335
## Ca.Mg Mg.K Ca.K Ca.Mg.K Ca.B Fe.Mn Fe.Zn Rmed2
## 1 4.207091 2.613154 10.679175 13.431593 22.91193 76.11250 159.48352 NA
## 2 4.650541 4.041579 23.226207 27.267786 22.10321 99.79597 224.82995 19.510
## 3 3.952252 1.244362 4.686444 5.946376 13.35808 101.88795 348.51124 19.040
## 4 4.150889 2.791340 13.325539 16.116879 28.28309 82.06782 44.59346 18.000
## 5 6.157040 2.107519 13.378390 15.714203 31.35667 106.15040 192.95755 15.670
## 6 3.411689 2.389101 8.001158 10.422192 44.79601 127.84451 469.15464 15.460
## 7 4.155828 1.600995 5.741427 7.352359 15.30352 64.50734 279.89002 28.945
## 8 4.019314 3.544209 11.898310 15.463167 21.59959 127.99761 231.20179 20.020
## 9 2.313750 1.995038 4.270404 6.265442 10.59602 137.53683 58.32030 26.000
## 10 3.124588 2.223052 5.858726 8.184432 23.36961 124.07392 283.95300 15.685
## Rmed3 Gr.km Gr.hcpc
## 1 NA 2 2
## 2 11.670 2 2
## 3 13.950 1 1
## 4 NA 3 3
## 5 11.690 2 2
## 6 14.470 1 1
## 7 17.885 2 2
## 8 17.890 2 2
## 9 15.725 3 3
## 10 16.660 2 25. E X P L O R A T O R I O E N T R E G R U P O S
5.1. Graficos con datos agrupados por municipio
Annual potato yield
#
DF.cluster.tot.pa <- DF.cluster.tot %>% filter(., !is.na(Rmed2) )
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot.pa$Rmed2, DF.cluster.tot.pa$Gr.hcpc)
kw.l <- kruskal(y = DF.cluster.tot.pa$Rmed2, trt = DF.cluster.tot.pa$Gr.hcpc,
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot.pa$Rmed2 ~ DF.cluster.tot.pa$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 2.013408
## Degrees of freedom: 2
## Pvalue Chisq : 0.3654214
##
## DF.cluster.tot.pa$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.pa.Rmed2 r
## 1 19.09091 11
## 2 25.91667 24
## 3 22.63636 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.491261
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot.pa$Rmed2 groups
## 2 25.91667 a
## 3 22.63636 a
## 1 19.09091 akw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.pa.Rmed2, Kruskall=groups) %>%
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot.pa$Rmed2) * 1.3) {
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(which(!is.na(y))), "\n",
"Mean =", round(mean(y, na.rm=TRUE), 2), "\n",
"Median =", round(median(y, na.rm=TRUE), 2), "\n"
)
))
}
#
ggplot(DF.cluster.tot.pa, aes(x=Gr.hcpc, y=Rmed2, colour=Gr.hcpc) ) +
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot.pa$Rmed2)*0.4,
label = Kruskall), size=5)+
ylim(0, 40) + # Limites eje y para visualizar mejor
geom_hline(yintercept = 22, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo
ylab("Annual potato yield (t/ha)") # Etiquetas de eje, con unidades
Rendimiento papa criolla
#
DF.cluster.tot.pc <- DF.cluster.tot %>% filter(., !is.na(Rmed3) )
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot.pc$Rmed3, DF.cluster.tot.pc$Gr.hcpc)
kw.l <- kruskal(y = DF.cluster.tot.pc$Rmed3, trt = DF.cluster.tot.pc$Gr.hcpc,
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot.pc$Rmed3 ~ DF.cluster.tot.pc$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 2.67051
## Degrees of freedom: 2
## Pvalue Chisq : 0.2630911
##
## DF.cluster.tot.pc$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.pc.Rmed3 r
## 1 14.25000 10
## 2 20.57143 21
## 3 21.41667 6
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.518259
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot.pc$Rmed3 groups
## 3 21.41667 a
## 2 20.57143 a
## 1 14.25000 akw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.pc.Rmed3, Kruskall=groups) %>%
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot.pc$Rmed3) * 1.7) {
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(which(!is.na(y))), "\n",
"Mean =", round(mean(y, na.rm=TRUE), 2), "\n",
"Median =", round(median(y, na.rm=TRUE), 2), "\n"
)
))
}
#
ggplot(DF.cluster.tot.pc, aes(x=Gr.hcpc, y=Rmed3, colour=Gr.hcpc) ) +
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot.pc$Rmed3)*0.4,
label = Kruskall), size=5)+
ylim(0, 40) + # Limites eje y para visualizar mejor
geom_hline(yintercept = 14, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo
ylab("Criolla potato yield (t/ha)") # Etiquetas de eje, con unidades
Altitud
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$ALT, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$ALT, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$ALT ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 0.9909717
## Degrees of freedom: 2
## Pvalue Chisq : 0.6092748
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.ALT r
## 1 23.68182 11
## 2 25.60000 25
## 3 20.68182 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$ALT groups
## 2 25.60000 a
## 1 23.68182 a
## 3 20.68182 akw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.ALT, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$ALT) * 1.6) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 0), "\n",
"Median =", round(median(y), 0), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=ALT, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.2) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$ALT)*0.8,###
label = Kruskall), size=5)+
ylim(0, 4500) + # Limites eje y para visualizar mejor###
geom_hline(yintercept =2300, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept =3200, linetype="dashed", color = "red3", size=0.6) + # Valor optimo máximo###
ylab("Altitude (masl)")# Etiquetas de eje, con unidades###
pH
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$pH, DF.cluster.tot$Gr.hcpc)
kw.l <- kruskal(y = DF.cluster.tot$pH, trt = DF.cluster.tot$Gr.hcpc,
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$pH ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 16.91482
## Degrees of freedom: 2
## Pvalue Chisq : 0.0002123217
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.pH r
## 1 9.636364 11
## 2 26.760000 25
## 3 32.090909 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$pH groups
## 3 32.090909 a
## 2 26.760000 a
## 1 9.636364 bkw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.pH, Kruskall=groups) %>%
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$pH) * 1.2) {
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 1), "\n",
"Median =", round(median(y), 1), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=pH, colour=Gr.hcpc) ) +
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$pH)*0.9,
label = Kruskall), size=5)+
ylim(4, 7) + # Limites eje y para visualizar mejor
geom_hline(yintercept = 5.1, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo
geom_hline(yintercept = 6.0, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo
ylab("pH") # Etiquetas de eje, con unidades
Calcio Ca
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$Ca, DF.cluster.tot$Gr.hcpc)
kw.l <- kruskal(y = DF.cluster.tot$Ca, trt = DF.cluster.tot$Gr.hcpc,
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$Ca ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 32.79691
## Degrees of freedom: 2
## Pvalue Chisq : 7.55514e-08
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.Ca r
## 1 6.00000 11
## 2 25.20000 25
## 3 39.27273 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$Ca groups
## 3 39.27273 a
## 2 25.20000 b
## 1 6.00000 ckw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.Ca, Kruskall=groups) %>%
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$Ca) * 1.3) {
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 2), "\n",
"Median =", round(median(y), 2), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=Ca, colour=Gr.hcpc) ) +
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$Ca)*0.5,
label = Kruskall), size=5)+
ylim(0, 15) + # Limites eje y para visualizar mejor
geom_hline(yintercept = 3, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo
geom_hline(yintercept = 10, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo
ylab("Ca (cmol(+)/kg)") # Etiquetas de eje, con unidades
Materia organica OM
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$OM, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$OM, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$OM ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 2.097795
## Degrees of freedom: 2
## Pvalue Chisq : 0.3503238
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.OM r
## 1 24.90909 11
## 2 21.56000 25
## 3 28.63636 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$OM groups
## 3 28.63636 a
## 1 24.90909 a
## 2 21.56000 akw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.OM, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$OM) * 1.25) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 1), "\n",
"Median =", round(median(y), 1), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=OM, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$OM)*0.4,###
label = Kruskall), size=5)+
ylim(0, 25) + # Limites eje y para visualizar mejor###
geom_hline(yintercept = 5.0, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept = 10, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo###
ylab("Organic matter (%)") # Etiquetas de eje, con unidades###
Capacidad de intercambio cationico CEC
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$CEC, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$CEC, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$CEC ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 28.06963
## Degrees of freedom: 2
## Pvalue Chisq : 8.030762e-07
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.CEC r
## 1 9.727273 11
## 2 23.000000 25
## 3 40.545455 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$CEC groups
## 3 40.545455 a
## 2 23.000000 b
## 1 9.727273 ckw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.CEC, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$CEC) * 1.6) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 1), "\n",
"Median =", round(median(y), 1), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=CEC, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$CEC)*0.4,###
label = Kruskall), size=5)+
ylim(0, 28) + # Limites eje y para visualizar mejor###
geom_hline(yintercept = 20, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept = 10, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo###
ylab("Cation exchange capacity (meq/100g)") # Etiquetas de eje, con unidades###
Conductividad electrica EC
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$EC, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$EC, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$EC ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 31.36944
## Degrees of freedom: 2
## Pvalue Chisq : 1.542456e-07
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.EC r
## 1 9.636364 11
## 2 22.400000 25
## 3 42.000000 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$EC groups
## 3 42.000000 a
## 2 22.400000 b
## 1 9.636364 ckw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.EC, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$EC) * 1.6) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 1), "\n",
"Median =", round(median(y), 1), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=EC, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$EC)*0.4,###
label = Kruskall), size=5)+
ylim(0, 3.4) + # Limites eje y para visualizar mejor###
geom_hline(yintercept = 2, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept = 0.5, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo###
ylab("Electrical conductivity (dS/m)") # Etiquetas de eje, con unidades###
Sodio Na
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$Na, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$Na, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$Na ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 24.11501
## Degrees of freedom: 2
## Pvalue Chisq : 5.800857e-06
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.Na r
## 1 11.90909 11
## 2 22.24000 25
## 3 40.09091 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$Na groups
## 3 40.09091 a
## 2 22.24000 b
## 1 11.90909 ckw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.Na, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$Na) * 1.3) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 2), "\n",
"Median =", round(median(y), 2), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=Na, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$Na)*0.4,###
label = Kruskall), size=5)+
ylim(0, 1) + # Limites eje y para visualizar mejor###
geom_hline(yintercept = 0.5, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept = 0.1, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo###
ylab("Na (meq/100g)") # Etiquetas de eje, con unidades###
Potasio K
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$K, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$K, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$K ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 17.88789
## Degrees of freedom: 2
## Pvalue Chisq : 0.000130525
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.K r
## 1 12.45455 11
## 2 23.32000 25
## 3 37.09091 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$K groups
## 3 37.09091 a
## 2 23.32000 b
## 1 12.45455 ckw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.K, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$K) * 1.35) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 2), "\n",
"Median =", round(median(y), 2), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=K, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$K)*0.4,###
label = Kruskall), size=5)+
ylim(0, 2.6) + # Limites eje y para visualizar mejor###
geom_hline(yintercept = 0.4, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept = 0.2, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo###
ylab("K (cmol(+)/kg)") # Etiquetas de eje, con unidades###
Magnesio Mg
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$Mg, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$Mg, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$Mg ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 30.47072
## Degrees of freedom: 2
## Pvalue Chisq : 2.417509e-07
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.Mg r
## 1 7.909091 11
## 2 23.960000 25
## 3 40.181818 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$Mg groups
## 3 40.181818 a
## 2 23.960000 b
## 1 7.909091 ckw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.Mg, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$Mg) * 1.5) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 2), "\n",
"Median =", round(median(y), 2), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=Mg, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$Mg)*0.4,###
label = Kruskall), size=5)+
ylim(0, 5.5) + # Limites eje y para visualizar mejor###
geom_hline(yintercept = 3.0, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept = 1.5, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo###
ylab("Mg (cmol(+)/kg)") # Etiquetas de eje, con unidades###
Saturacion de aluminio Alsat
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$Alsat, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$Alsat, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$Alsat ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 28.01733
## Degrees of freedom: 2
## Pvalue Chisq : 8.243543e-07
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.Alsat r
## 1 42.00000 11
## 2 21.24000 25
## 3 12.27273 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$Alsat groups
## 1 42.00000 a
## 2 21.24000 b
## 3 12.27273 ckw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.Alsat, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$Alsat) * 1.5) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 1), "\n",
"Median =", round(median(y), 1), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=Alsat, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$Alsat)*0.4,###
label = Kruskall), size=5)+
ylim(0, 60) + # Limites eje y para visualizar mejor###
geom_hline(yintercept = 20, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo###
ylab("Aluminum saturation (%)") # Etiquetas de eje, con unidades###
Azufre S
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$S, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$S, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$S ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 25.94197
## Degrees of freedom: 2
## Pvalue Chisq : 2.32687e-06
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.S r
## 1 13.27273 11
## 2 21.00000 25
## 3 41.54545 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$S groups
## 3 41.54545 a
## 2 21.00000 b
## 1 13.27273 bkw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.S, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$S) * 1.3) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 2), "\n",
"Median =", round(median(y), 2), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=S, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$S)*0.4,###
label = Kruskall), size=5)+
ylim(0, 100) + # Limites eje y para visualizar mejor###
geom_hline(yintercept =20, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept =10, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo###
ylab("S (ppm)")# Etiquetas de eje, con unidades###
Aluminio Al
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$Al, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$Al, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$Al ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 25.89215
## Degrees of freedom: 2
## Pvalue Chisq : 2.385568e-06
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.Al r
## 1 41.45455 11
## 2 21.12000 25
## 3 13.09091 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$Al groups
## 1 41.45455 a
## 2 21.12000 b
## 3 13.09091 bkw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.Al, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$Al) * 1.4) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 2), "\n",
"Median =", round(median(y), 2), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=Al, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$Al)*0.4,###
label = Kruskall), size=5)+
ylim(0, 4) + # Limites eje y para visualizar mejor###
geom_hline(yintercept =1, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo###
ylab("Al (meq/100g)")# Etiquetas de eje, con unidades###
Hierro Fe
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$Fe, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$Fe, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$Fe ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 6.126422
## Degrees of freedom: 2
## Pvalue Chisq : 0.04673739
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.Fe r
## 1 28.18182 11
## 2 26.04000 25
## 3 15.18182 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$Fe groups
## 1 28.18182 a
## 2 26.04000 a
## 3 15.18182 akw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.Fe, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$Fe) * 1.1) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 1), "\n",
"Median =", round(median(y), 1), "\n"
) ))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=Fe, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.1) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$Fe)*0.1,###
label = Kruskall), size=5)+
ylim(0, 1400) + # Limites eje y para visualizar mejor###
geom_hline(yintercept =100, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept =30, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo###
ylab("Fe (mg/kg)")# Etiquetas de eje, con unidades###
Manganeso Mn
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$Mn, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$Mn, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$Mn ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 3.96325
## Degrees of freedom: 2
## Pvalue Chisq : 0.1378451
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.Mn r
## 1 17.72727 11
## 2 27.40000 25
## 3 22.54545 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$Mn groups
## 2 27.40000 a
## 3 22.54545 a
## 1 17.72727 akw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.Mn, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$Mn) * 1.9) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 2), "\n",
"Median =", round(median(y), 2), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=Mn, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$Al)*0.4,###
label = Kruskall), size=5)+
ylim(0, 27) + # Limites eje y para visualizar mejor###
geom_hline(yintercept =20, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept =10, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo###
ylab("Mn (mg/kg)")# Etiquetas de eje, con unidades###
Boro B
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$B, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$B, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$B ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 27.67188
## Degrees of freedom: 2
## Pvalue Chisq : 9.797797e-07
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.B r
## 1 11.27273 11
## 2 21.96000 25
## 3 41.36364 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$B groups
## 3 41.36364 a
## 2 21.96000 b
## 1 11.27273 ckw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.B, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$B) * 1.7) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 2), "\n",
"Median =", round(median(y), 2), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=B, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$Cu)*0.04,###
label = Kruskall), size=5)+
ylim(0, 1.0) + # Limites eje y para visualizar mejor###
geom_hline(yintercept =0.6, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept =0.3, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo###
ylab("B (mg/kg)")# Etiquetas de eje, con unidades###
Zinc Zn
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$Zn, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$Zn, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$Zn ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 29.11652
## Degrees of freedom: 2
## Pvalue Chisq : 4.758043e-07
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.Zn r
## 1 8.363636 11
## 2 23.880000 25
## 3 39.909091 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$Zn groups
## 3 39.909091 a
## 2 23.880000 b
## 1 8.363636 ckw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.Zn, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$Zn) * 1.25) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 2), "\n",
"Median =", round(median(y), 2), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=Zn, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$Zn)*0.4,###
label = Kruskall), size=5)+
ylim(0, 36) + # Limites eje y para visualizar mejor###
geom_hline(yintercept =4, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept =3, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo###
ylab("Zn (mg/kg)")# Etiquetas de eje, con unidades###
### Cobre Cu
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$Cu, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$Cu, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$Cu ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 16.48549
## Degrees of freedom: 2
## Pvalue Chisq : 0.0002631605
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.Cu r
## 1 9.363636 11
## 2 29.000000 25
## 3 27.272727 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$Cu groups
## 2 29.000000 a
## 3 27.272727 a
## 1 9.363636 bkw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.Cu, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$Cu) * 1.25) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 2), "\n",
"Median =", round(median(y), 2), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=Cu, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$Zn)*0.1,###
label = Kruskall), size=5)+
ylim(0, 11) + # Limites eje y para visualizar mejor###
geom_hline(yintercept =3, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept =1, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo###
ylab("Cu (mg/kg)")# Etiquetas de eje, con unidades###
### Saturación de potasio S.K
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$S.K, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$S.K, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$S.K ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 1.41501
## Degrees of freedom: 2
## Pvalue Chisq : 0.4928725
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.S.K r
## 1 21.09091 11
## 2 23.56000 25
## 3 27.90909 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$S.K groups
## 3 27.90909 a
## 2 23.56000 a
## 1 21.09091 akw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.S.K, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$S.K) * 1.5) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 1), "\n",
"Median =", round(median(y), 1), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=S.K, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$S.K)*0.4,###
label = Kruskall), size=5)+
ylim(0, 20) + # Limites eje y para visualizar mejor###
geom_hline(yintercept =6, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept =3, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo###
ylab("Potassium saturation (%)")# Etiquetas de eje, con unidades###
Saturacion de calcio S.Ca
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$S.Ca, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$S.Ca, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$S.Ca ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 23.68503
## Degrees of freedom: 2
## Pvalue Chisq : 7.192192e-06
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.S.Ca r
## 1 6.545455 11
## 2 28.360000 25
## 3 31.545455 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$S.Ca groups
## 3 31.545455 a
## 2 28.360000 a
## 1 6.545455 bkw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.S.Ca, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$S.Ca) * 1.4) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 1), "\n",
"Median =", round(median(y), 1), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=S.Ca, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$S.Ca)*0.4,###
label = Kruskall), size=5)+
ylim(0, 100) + # Limites eje y para visualizar mejor###
geom_hline(yintercept =70, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept =50, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo###
ylab("Calcium saturation (%)")# Etiquetas de eje, con unidades###
Saturacion de sodio S.Na
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$S.Na, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$S.Na, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$S.Na ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 12.19041
## Degrees of freedom: 2
## Pvalue Chisq : 0.002253652
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.S.Na r
## 1 18.18182 11
## 2 21.08000 25
## 3 36.45455 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$S.Na groups
## 3 36.45455 a
## 2 21.08000 b
## 1 18.18182 bkw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.S.Na, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$S.Na) * 2) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 1), "\n",
"Median =", round(median(y), 1), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=S.Na, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$S.Ca)*0.02,###
label = Kruskall), size=5)+
ylim(0, 11) + # Limites eje y para visualizar mejor###
geom_hline(yintercept = 7, linetype="dashed", color = "red3", size=0.6) +
geom_hline(yintercept = 5, linetype="dashed", color = "blue3", size=0.6) +# Valor optimo mínimo###
ylab("Sodium saturation (%)")# Etiquetas de eje, con unidades###
Saturacion de magnesio S.Mg
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$S.Mg, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$S.Mg, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$S.Mg ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 13.51636
## Degrees of freedom: 2
## Pvalue Chisq : 0.001161339
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.S.Mg r
## 1 12.45455 11
## 2 24.76000 25
## 3 33.81818 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$S.Mg groups
## 3 33.81818 a
## 2 24.76000 a
## 1 12.45455 bkw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.S.Mg, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$S.Mg) *1.6) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 1), "\n",
"Median =", round(median(y), 1), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=S.Mg, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$S.Ca)*0.04,###
label = Kruskall), size=5)+
ylim(0, 38) + # Limites eje y para visualizar mejor###
geom_hline(yintercept =25, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept =10, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo###
ylab("Magnesium saturation (%)")# Etiquetas de eje, con unidades###
Saturacion de bases BS
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$BS, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$BS, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$BS ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 28.11342
## Degrees of freedom: 2
## Pvalue Chisq : 7.856835e-07
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.BS r
## 1 6.00000 11
## 2 26.72000 25
## 3 35.81818 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$BS groups
## 3 35.81818 a
## 2 26.72000 b
## 1 6.00000 ckw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.BS, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$BS) * 1.3) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 1), "\n",
"Median =", round(median(y), 1), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=BS, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$BS)*0.4,###
label = Kruskall), size=5)+
ylim(0, 130) + # Limites eje y para visualizar mejor###
geom_hline(yintercept =50, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept =35, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo###
ylab("Base saturation (%)")# Etiquetas de eje, con unidades###
Relacion Calcio Magnesio Ca.Mg
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$Ca.Mg, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$Ca.Mg, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$Ca.Mg ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 6.062785
## Degrees of freedom: 2
## Pvalue Chisq : 0.0482484
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.Ca.Mg r
## 1 15.54545 11
## 2 27.76000 25
## 3 23.90909 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$Ca.Mg groups
## 2 27.76000 a
## 3 23.90909 ab
## 1 15.54545 bkw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.Ca.Mg, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$Ca.Mg) * 1.4) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 1), "\n",
"Median =", round(median(y), 1), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=Ca.Mg, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$Ca.Mg)*0.4,###
label = Kruskall), size=5)+
ylim(0, 12) + # Limites eje y para visualizar mejor###
geom_hline(yintercept =5, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept =3, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo###
ylab("Ca/Mg ratio")# Etiquetas de eje, con unidades###
Relacion magnesio potasio Mg.K
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$Mg.K, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$Mg.K, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$Mg.K ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 5.120232
## Degrees of freedom: 2
## Pvalue Chisq : 0.07729577
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.Mg.K r
## 1 17.81818 11
## 2 23.64000 25
## 3 31.00000 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$Mg.K groups
## 3 31.00000 a
## 2 23.64000 a
## 1 17.81818 akw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.Mg.K, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$Mg.K) * 1.8) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 1), "\n",
"Median =", round(median(y), 1), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=Mg.K, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$Mg.K)*0.04,###
label = Kruskall), size=5)+
ylim(0, 11) + # Limites eje y para visualizar mejor###
geom_hline(yintercept =8, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept =6, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo###
ylab("Mg/K Ratio")# Etiquetas de eje, con unidades###
Relacion calcio potasio Ca.K
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$Ca.K, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$Ca.K, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$Ca.K ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 8.059265
## Degrees of freedom: 2
## Pvalue Chisq : 0.01778086
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.Ca.K r
## 1 13.90909 11
## 2 26.28000 25
## 3 28.90909 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$Ca.K groups
## 3 28.90909 a
## 2 26.28000 a
## 1 13.90909 bkw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.Ca.K, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$Ca.K) * 0.8) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 1), "\n",
"Median =", round(median(y), 1), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=Mg.K, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$Ca.K)*0.04,###
label = Kruskall), size=5)+
ylim(0, 24) + # Limites eje y para visualizar mejor###
geom_hline(yintercept =18, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept =12, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo###
ylab("Ca/K Ratio")# Etiquetas de eje, con unidades###
Relacion calcio magnesio sobre potasio Ca.Mg.K
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$Ca.Mg.K, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$Ca.Mg.K, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$Ca.Mg.K ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 8.346692
## Degrees of freedom: 2
## Pvalue Chisq : 0.01540064
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.Ca.Mg.K r
## 1 13.81818 11
## 2 26.12000 25
## 3 29.36364 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$Ca.Mg.K groups
## 3 29.36364 a
## 2 26.12000 a
## 1 13.81818 bkw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.Ca.Mg.K, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$Ca.Mg.K) * 1.3) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 1), "\n",
"Median =", round(median(y), 1), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=Ca.Mg.K, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$Ca.Mg.K)*0.04,###
label = Kruskall), size=5)+
ylim(0, 43) + # Limites eje y para visualizar mejor###
geom_hline(yintercept =20, linetype="dashed", color = "red3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept =12, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo máximo###
ylab("(Ca+Mg)/K Ratio")# Etiquetas de eje, con unidades###
Acidez EA
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$EA, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$EA, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$EA ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 26.853
## Degrees of freedom: 2
## Pvalue Chisq : 1.475522e-06
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.EA r
## 1 41.81818 11
## 2 21.00000 25
## 3 13.00000 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$EA groups
## 1 41.81818 a
## 2 21.00000 b
## 3 13.00000 bkw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.EA, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$EA) * 1.35) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 2), "\n",
"Median =", round(median(y), 2), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=EA, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.5) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$EA)*0.04,###
label = Kruskall), size=5)+
ylim(0, 4.5) + # Limites eje y para visualizar mejor###
geom_hline(yintercept =2, linetype="dashed", color = "red", size=0.6) + # Valor optimo máximo###
ylab("Exchangeable acidity (meq/100g)")# Etiquetas de eje, con unidades###
Relacion Hierro Manganeso Fe.Mn
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$Fe.Mn, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$Fe.Mn, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$Fe.Mn ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 9.575164
## Degrees of freedom: 2
## Pvalue Chisq : 0.00833258
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.Fe.Mn r
## 1 33.09091 11
## 2 23.96000 25
## 3 15.00000 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$Fe.Mn groups
## 1 33.09091 a
## 2 23.96000 ab
## 3 15.00000 bkw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.Fe.Mn, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$Fe.Mn) * 1.25) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 1), "\n",
"Median =", round(median(y), 1), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=Fe.Mn, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.2) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$Fe.Mn)*0.04,###
label = Kruskall), size=5)+
ylim(0, 240) + # Limites eje y para visualizar mejor###
geom_hline(yintercept =8, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept =10, linetype="dashed", color = "red3", size=0.6) + # Valor optimo máximo###
ylab("Fe/Mn Ratio")# Etiquetas de eje, con unidades###
Relacion hierro zinc Fe.Zn
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$Fe.Zn, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$Fe.Zn, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$Fe.Zn ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 23.12217
## Degrees of freedom: 2
## Pvalue Chisq : 9.529835e-06
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.Fe.Zn r
## 1 37.636364 11
## 2 24.360000 25
## 3 9.545455 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$Fe.Zn groups
## 1 37.636364 a
## 2 24.360000 b
## 3 9.545455 ckw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.Fe.Zn, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$Fe.Zn) * 1.25) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 2), "\n",
"Median =", round(median(y), 2), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=Fe.Zn, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.2) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$Fe.Zn)*0.04,###
label = Kruskall), size=5)+
ylim(0, 950) + # Limites eje y para visualizar mejor###
geom_hline(yintercept =15, linetype="dashed", color = "red3", size=0.6) + # Valor optimo máximo###
ylab("Fe/Zn Ratio")# Etiquetas de eje, con unidades###
<<<### Cobre Cu
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$Cu, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$Cu, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$Cu ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 16.48549
## Degrees of freedom: 2
## Pvalue Chisq : 0.0002631605
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.Cu r
## 1 9.363636 11
## 2 29.000000 25
## 3 27.272727 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$Cu groups
## 2 29.000000 a
## 3 27.272727 a
## 1 9.363636 bkw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.Cu, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$Cu) * 1.3) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 2), "\n",
"Median =", round(median(y), 2), "\n"
)
))
}<
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=Cu, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.2) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$Cu)*0.04,###
label = Kruskall), size=5)+
ylim(0, 11.5) + # Limites eje y para visualizar mejor###
geom_hline(yintercept =2, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept =3, linetype="dashed", color = "red3", size=0.6) + # Valor optimo máximo###
ylab("Cu (ppm)")# Etiquetas de eje, con unidades###Fosforo
#
# Kruskall wallis test
kw.t <- kruskal.test(DF.cluster.tot$P, DF.cluster.tot$Gr.hcpc)###
kw.l <- kruskal(y = DF.cluster.tot$P, trt = DF.cluster.tot$Gr.hcpc,###
console=T, p.adj="bonferroni", group=T, alpha = 0.05)##
## Study: DF.cluster.tot$P ~ DF.cluster.tot$Gr.hcpc
## Kruskal-Wallis test's
## Ties or no Ties
##
## Critical Value: 4.644178
## Degrees of freedom: 2
## Pvalue Chisq : 0.09806851
##
## DF.cluster.tot$Gr.hcpc, means of the ranks
##
## DF.cluster.tot.P r
## 1 16.54545 11
## 2 25.32000 25
## 3 28.45455 11
##
## Post Hoc Analysis
##
## P value adjustment method: bonferroni
## t-Student: 2.488968
## Alpha : 0.05
## Groups according to probability of treatment differences and alpha level.
##
## Treatments with the same letter are not significantly different.
##
## DF.cluster.tot$P groups
## 3 28.45455 a
## 2 25.32000 a
## 1 16.54545 akw.etiqu <- data.frame(rownames_to_column(kw.l$groups)) %>%
mutate(., rowname=as.numeric(rowname) ) %>%
arrange(., rowname) %>%
rename(., Gr.hcpc=rowname, Promedio=DF.cluster.tot.P, Kruskall=groups) %>% ###
mutate(., Gr.hcpc=as.factor(Gr.hcpc) )
# Gráfico
# Referencia: https://appsilon.com/ggplot2-boxplots/
get_box_stats <- function(y, upper_limit = max(DF.cluster.tot$P) * 0.9) {###
return(data.frame(
y = 0.95 * upper_limit,
label = paste(
"N =", length(y), "\n",
"Mean =", round(mean(y), 2), "\n",
"Median =", round(median(y), 2), "\n"
)
))
}
#
ggplot(DF.cluster.tot, aes(x=Gr.hcpc, y=P, colour=Gr.hcpc) ) +###
geom_jitter(alpha=0.3, width = 0.0)+
geom_boxplot(outlier.shape = T, alpha=0.8)+
scale_fill_manual(values = color.pal) +
scale_color_manual(values = color.pal)+
theme_classic()+
theme(legend.position = "none")+
xlab("Groups") +
stat_summary(fun.data = get_box_stats, geom = "text",
hjust = 0.5, vjust = 0.9, size=3.2) +
geom_text(data = kw.etiqu,
aes(x = Gr.hcpc, y = min(DF.cluster.tot$P)*0.04,###
label = Kruskall), size=5)+
ylim(0, 160) + # Limites eje y para visualizar mejor###
geom_hline(yintercept =15, linetype="dashed", color = "blue3", size=0.6) + # Valor optimo mínimo###
geom_hline(yintercept =40, linetype="dashed", color = "red3", size=0.6) + # Valor optimo máximo###
ylab("P (ppm)")# Etiquetas de eje, con unidades###