paquetes <- c(“ggplot2”, “dplyr”, “tidyr”, “car”, “scales”)
for (p in paquetes) { if (!requireNamespace(p, quietly = TRUE)) { install.packages(p, repos = “https://cloud.r-project.org”) } }
library(ggplot2) library(dplyr) library(tidyr) library(car) library(scales)
tema_agricola <- theme_minimal(base_size = 13) + theme( plot.title = element_text( face = “bold”, size = 17, hjust = 0.5, color = “#24543D” ), plot.subtitle = element_text( size = 11, hjust = 0.5 ), axis.title = element_text( face = “bold” ), legend.position = “bottom”, panel.grid.minor = element_blank() )
cat(“”) cat(“============================================================”) cat(” ANALISIS AGROMETEOROLOGICO DE PANAMA“) cat(”============================================================“) cat(”Cultivos: Arroz y Maiz“) cat(”Provincias: Herrera y Los Santos“) cat(”Periodo: 2014-2015 a 2016-2017“) cat(”Analisis: ANOVA, pruebas estadisticas y correlaciones“) cat(”============================================================“)
datos <- data.frame( Anio = c( “2014-2015”, “2014-2015”, “2015-2016”, “2015-2016”, “2016-2017”, “2016-2017” ),
Provincia = c( “Herrera”, “Los Santos”, “Herrera”, “Los Santos”, “Herrera”, “Los Santos” ),
Precipitacion = c( 940.5, 704.0, 860.6, 644.3, 1371.3, 1232.1 ),
Maiz = c( 182924, 1378443, 440500, 1302800, 470000, 1701100 ),
Arroz = c( 152922, 702391, 66412, 648100, 116953, 936960 ) )
datos\(Anio <- factor( datos\)Anio, levels = c(“2014-2015”, “2015-2016”, “2016-2017”) )
datos\(Provincia <- factor( datos\)Provincia, levels = c(“Herrera”, “Los Santos”) )
cat(“ORIGINAL DE DATOS”) print(datos)
datos_largos <- datos %>% tidyr::pivot_longer( cols = c(Maiz, Arroz), names_to = “Cultivo”, values_to = “Produccion” ) %>% dplyr::mutate( Cultivo = base::ifelse( Cultivo == “Maiz”, “Maiz”, “Arroz” ) )
cat(“EN FORMATO LARGO”) print(datos_largos)
cat(“============================================================”) cat(“ESTADISTICA DESCRIPTIVA DE LA PRODUCCION”) cat(“============================================================”)
descriptivos <- datos_largos %>% dplyr::group_by(Provincia, Cultivo) %>% dplyr::summarise( N = dplyr::n(), Media = mean(Produccion), Desviacion = sd(Produccion), Minimo = min(Produccion), Maximo = max(Produccion), .groups = “drop” )
print(descriptivos)
cat(“DESCRIPTIVA DE PRECIPITACION”)
descriptivos_pp <- datos %>% dplyr::group_by(Provincia) %>% dplyr::summarise( N = dplyr::n(), Media = mean(Precipitacion), Desviacion = sd(Precipitacion), Minimo = min(Precipitacion), Maximo = max(Precipitacion), .groups = “drop” )
print(descriptivos_pp)
cat(“============================================================”) cat(“ANOVA PARA LA PRODUCCION DE MAIZ”) cat(“============================================================”)
modelo_maiz <- aov(Maiz ~ Provincia, data = datos)
resultado_maiz <- summary(modelo_maiz)
print(resultado_maiz)
p_maiz <- resultado_maiz[[1]][[“Pr(>F)”]][1]
cat(“p del ANOVA de maiz:”, p_maiz, “”)
if (p_maiz < 0.05) { cat(“Se rechaza H0 al nivel de significacion del 5%.”) } else { cat(“No se rechaza H0 al nivel de significacion del 5%.”) }
cat(“============================================================”) cat(“ANOVA PARA LA PRODUCCION DE ARROZ”) cat(“============================================================”)
modelo_arroz <- aov(Arroz ~ Provincia, data = datos)
resultado_arroz <- summary(modelo_arroz)
print(resultado_arroz)
p_arroz <- resultado_arroz[[1]][[“Pr(>F)”]][1]
cat(“p del ANOVA de arroz:”, p_arroz, “”)
if (p_arroz < 0.05) { cat(“Se rechaza H0 al nivel de significacion del 5%.”) } else { cat(“No se rechaza H0 al nivel de significacion del 5%.”) }
cat(“============================================================”) cat(“PRUEBAS DE NORMALIDAD SHAPIRO-WILK”) cat(“============================================================”)
cat(“:”) normalidad_maiz <- shapiro.test(residuals(modelo_maiz)) print(normalidad_maiz)
cat(“:”) normalidad_arroz <- shapiro.test(residuals(modelo_arroz)) print(normalidad_arroz)
cat(“============================================================”) cat(“PRUEBAS DE HOMOGENEIDAD DE VARIANZAS”) cat(“============================================================”)
cat(“:”) levene_maiz <- car::leveneTest( Maiz ~ Provincia, data = datos ) print(levene_maiz)
cat(“:”) levene_arroz <- car::leveneTest( Arroz ~ Provincia, data = datos ) print(levene_arroz)
cat(“============================================================”) cat(“PRUEBAS NO PARAMETRICAS DE KRUSKAL-WALLIS”) cat(“============================================================”)
cat(“:”) kruskal_maiz <- kruskal.test( Maiz ~ Provincia, data = datos ) print(kruskal_maiz)
cat(“:”) kruskal_arroz <- kruskal.test( Arroz ~ Provincia, data = datos ) print(kruskal_arroz)
cat(“============================================================”) cat(“CORRELACION DE SPEARMAN: PRECIPITACION Y MAIZ”) cat(“============================================================”)
cor_maiz <- cor.test( datos\(Precipitacion, datos\)Maiz, method = “spearman”, exact = FALSE )
print(cor_maiz)
cat(“============================================================”) cat(“CORRELACION DE SPEARMAN: PRECIPITACION Y ARROZ”) cat(“============================================================”)
cor_arroz <- cor.test( datos\(Precipitacion, datos\)Arroz, method = “spearman”, exact = FALSE )
print(cor_arroz)
grafico_pp <- ggplot( datos, aes( x = Anio, y = Precipitacion, color = Provincia, group = Provincia ) ) + geom_line(linewidth = 1.2) + geom_point(size = 4) + geom_text( aes(label = round(Precipitacion, 1)), vjust = -1, size = 3.5, show.legend = FALSE ) + scale_color_manual( values = c( “Herrera” = “#D99A18”, “Los Santos” = “#218C55” ) ) + labs( title = “Precipitacion anual por provincia”, subtitle = “Anios agricolas 2014-2015 a 2016-2017”, x = “Anio agricola”, y = “Precipitacion (mm)”, color = “Provincia” ) + tema_agricola
print(grafico_pp)
grafico_maiz <- ggplot( datos, aes( x = Anio, y = Maiz, fill = Provincia ) ) + geom_col( position = position_dodge(width = 0.8), width = 0.7 ) + scale_fill_manual( values = c( “Herrera” = “#E4AD31”, “Los Santos” = “#39865B” ) ) + scale_y_continuous( labels = scales::label_number(big.mark = “,”, decimal.mark = “.”) ) + labs( title = “Produccion de maiz por provincia”, x = “Anio agricola”, y = “Produccion (quintales)”, fill = “Provincia” ) + tema_agricola
print(grafico_maiz)
grafico_arroz <- ggplot( datos, aes( x = Anio, y = Arroz, fill = Provincia ) ) + geom_col( position = position_dodge(width = 0.8), width = 0.7 ) + scale_fill_manual( values = c( “Herrera” = “#E4AD31”, “Los Santos” = “#39865B” ) ) + scale_y_continuous( labels = scales::label_number(big.mark = “,”, decimal.mark = “.”) ) + labs( title = “Produccion de arroz por provincia”, x = “Anio agricola”, y = “Produccion (quintales)”, fill = “Provincia” ) + tema_agricola
print(grafico_arroz)
grafico_cor_maiz <- ggplot( datos, aes( x = Precipitacion, y = Maiz, color = Provincia, label = Anio ) ) + geom_point(size = 4) + geom_text( vjust = -1, size = 3, show.legend = FALSE ) + geom_smooth( method = “lm”, formula = y ~ x, se = FALSE, color = “#555555”, linetype = “dashed” ) + scale_color_manual( values = c( “Herrera” = “#D99A18”, “Los Santos” = “#218C55” ) ) + scale_y_continuous( labels = scales::label_number(big.mark = “,”, decimal.mark = “.”) ) + labs( title = “Relacion entre precipitacion y produccion de maiz”, x = “Precipitacion anual (mm)”, y = “Produccion de maiz (quintales)”, color = “Provincia” ) + tema_agricola
print(grafico_cor_maiz)
grafico_cor_arroz <- ggplot( datos, aes( x = Precipitacion, y = Arroz, color = Provincia, label = Anio ) ) + geom_point(size = 4) + geom_text( vjust = -1, size = 3, show.legend = FALSE ) + geom_smooth( method = “lm”, formula = y ~ x, se = FALSE, color = “#555555”, linetype = “dashed” ) + scale_color_manual( values = c( “Herrera” = “#D99A18”, “Los Santos” = “#218C55” ) ) + scale_y_continuous( labels = scales::label_number(big.mark = “,”, decimal.mark = “.”) ) + labs( title = “Relacion entre precipitacion y produccion de arroz”, x = “Precipitacion anual (mm)”, y = “Produccion de arroz (quintales)”, color = “Provincia” ) + tema_agricola
print(grafico_cor_arroz)
cat(“============================================================”) cat(“RESUMEN DE RESULTADOS”) cat(“============================================================”)
cat(“DE MAIZ: p =”, p_maiz) cat(“DE ARROZ: p =”, p_arroz)
cat(“PARA MAIZ:”) cat(“Rho =”, unname(cor_maiz\(estimate)) cat("\np =", cor_maiz\)p.value)
cat(“PARA ARROZ:”) cat(“Rho =”, unname(cor_arroz\(estimate)) cat("\np =", cor_arroz\)p.value)
cat(“del analisis.”)
write.csv( datos, “datos_agrometeorologicos.csv”, row.names = FALSE )
write.csv( descriptivos, “estadistica_descriptiva.csv”, row.names = FALSE )
ggsave( “01_precipitacion.png”, grafico_pp, width = 10, height = 6, dpi = 300 )
ggsave( “02_produccion_maiz.png”, grafico_maiz, width = 10, height = 6, dpi = 300 )
ggsave( “03_produccion_arroz.png”, grafico_arroz, width = 10, height = 6, dpi = 300 )
ggsave( “04_correlacion_maiz.png”, grafico_cor_maiz, width = 10, height = 6, dpi = 300 )
ggsave( “05_correlacion_arroz.png”, grafico_cor_arroz, width = 10, height = 6, dpi = 300 )
cat(“exportaron las tablas y los cinco graficos.”)