#UNIVERSIDAD CENTRAL DEL ECUADOR
#Facultad de Ingeniería en Geología, Minas, Petroleos y Ambiental
#Ingeniería Ambiental
#Autor: Grupo 4
#fecha:
#Cargar datos
setwd("C:/Users/LENOVO/OneDrive/Escritorio/ESTADISTICA")
datos <- read.csv("water_pollution_disease.csv")
#extraccion variable cuantitativa continua
Nivel_nitrato <- datos$Nitrate.Level..mg.L.
#Tabla de distribución de frecuencia
#Manualmente
min <- min(Nivel_nitrato)
max <- max(Nivel_nitrato)
R <- max-min
K <- floor(1+3.33*log10(length(Nivel_nitrato)))
A <- R/K
lim_inf <- round(seq(from=min,to=max-A,by=A),2)
lim_sup <- round(seq(from=min+A,to=max,by=A),2)
MC <- (lim_inf+lim_sup)/2
ni <- c()
for (i in 1:K) {
if (i < K) {
ni[i] <- length(subset(Nivel_nitrato, Nivel_nitrato >= lim_inf[i] & Nivel_nitrato < lim_sup[i]))
} else {
ni[i] <- length(subset(Nivel_nitrato, Nivel_nitrato >= lim_inf[i] & Nivel_nitrato <= lim_sup[i]))
}
}
sum(ni)
## [1] 3000
hi <- ni/sum(ni)*100
sum(hi)
## [1] 100
Ni_asc <- cumsum(ni)
Hi_asc <- cumsum(hi)
Ni_desc <- rev(cumsum(rev(ni)))
Hi_desc <- rev(cumsum(rev(hi)))
TDF_nitrato <- data.frame(lim_inf,
lim_sup,
MC,ni,
round(hi,2),
Ni_asc,
Ni_desc,
round(Hi_asc,2),
round(Hi_desc,2))
colnames(TDF_nitrato) <- c("Lim inf","Lim sup","MC","ni","hi(%)",
"Ni asc","Ni desc","Hi asc(%)","Hi desc(%)")
# crear de fila de totales
totales <- c( lim_inf= "TOTAL",
lim_sup= "-",
MC= "-",
ni= sum(ni),
hi= sum(hi),
Ni_asc= "-",
Ni_des= "-",
Hi_asc= "-",
Hi_des= "-")
TDF_nitrato_total <- rbind(TDF_nitrato, totales)
# Tabla mas estetica
library(knitr)
## Warning: package 'knitr' was built under R version 4.5.1
library(kableExtra)
## Warning: package 'kableExtra' was built under R version 4.5.1
kable(TDF_nitrato_total, align = 'c',
caption = "Tabla de Distribucion de Frecuencia del Nivel de Concentración
de Nitrato (mg/L) en los países analizados") %>%
kable_styling(full_width = FALSE, position = "center",
bootstrap_options = c("striped", "hover", "condensed"))
Tabla de Distribucion de Frecuencia del Nivel de Concentración de
Nitrato (mg/L) en los países analizados
|
Lim inf
|
Lim sup
|
MC
|
ni
|
hi(%)
|
Ni asc
|
Ni desc
|
Hi asc(%)
|
Hi desc(%)
|
|
0.05
|
4.21
|
2.13
|
235
|
7.83
|
235
|
3000
|
7.83
|
100
|
|
4.21
|
8.37
|
6.29
|
248
|
8.27
|
483
|
2765
|
16.1
|
92.17
|
|
8.37
|
12.54
|
10.455
|
268
|
8.93
|
751
|
2517
|
25.03
|
83.9
|
|
12.54
|
16.7
|
14.62
|
271
|
9.03
|
1022
|
2249
|
34.07
|
74.97
|
|
16.7
|
20.86
|
18.78
|
239
|
7.97
|
1261
|
1978
|
42.03
|
65.93
|
|
20.86
|
25.02
|
22.94
|
251
|
8.37
|
1512
|
1739
|
50.4
|
57.97
|
|
25.02
|
29.18
|
27.1
|
237
|
7.9
|
1749
|
1488
|
58.3
|
49.6
|
|
29.18
|
33.34
|
31.26
|
248
|
8.27
|
1997
|
1251
|
66.57
|
41.7
|
|
33.34
|
37.51
|
35.425
|
234
|
7.8
|
2231
|
1003
|
74.37
|
33.43
|
|
37.51
|
41.67
|
39.59
|
250
|
8.33
|
2481
|
769
|
82.7
|
25.63
|
|
41.67
|
45.83
|
43.75
|
252
|
8.4
|
2733
|
519
|
91.1
|
17.3
|
|
45.83
|
49.99
|
47.91
|
267
|
8.9
|
3000
|
267
|
100
|
8.9
|
|
TOTAL
|
|
|
3000
|
100
|
|
|
|
|
#Simplificación con el histograma
Hist_Nitrato <- hist(Nivel_nitrato,breaks = 10,plot=F)
k <- length(Hist_Nitrato$breaks)
Li <- Hist_Nitrato$breaks[1:(length(Hist_Nitrato$breaks)-1)]
Ls <- Hist_Nitrato$breaks[2:length(Hist_Nitrato$breaks)]
ni <- Hist_Nitrato$counts
sum(ni)
## [1] 3000
MC <- Hist_Nitrato$mids
hi <- (ni/sum(ni))
sum(hi)
## [1] 1
Ni_asc <- cumsum(ni)
Hi_asc <- cumsum(hi)
Ni_desc <- rev(cumsum(rev(ni)))
Hi_desc <- rev(cumsum(rev(hi)))
TDF_nitrato <- data.frame(Li = round(Li, 2),
Ls = round(Ls, 2),
MC = round(MC, 2),
ni = ni,
hi = round(hi * 100, 2),
Ni_asc = Ni_asc,
Ni_desc = Ni_desc,
Hi_asc = round(Hi_asc * 100, 2),
Hi_desc = round(Hi_desc * 100, 2))
colnames(TDF_nitrato) <- c("Lim inf","Lim sup","MC","ni","hi(%)",
"Ni asc","Ni desc","Hi asc(%)","Hi desc(%)")
# crear de fila de totales
totales <- c( lim_inf= "TOTAL",
lim_sup= "-",
MC= "-",
ni= sum(ni),
hi= sum(hi*100),
Ni_asc= "-",
Ni_des= "-",
Hi_asc= "-",
Hi_des= "-")
library(knitr)
library(kableExtra)
kable(TDF_nitrato_total, align = 'c',
caption = "Tabla de Distribución de frecuencias del nivel de concentración
de nitrato (mg/L) en los países analizados") %>%
kable_styling(full_width = FALSE, position = "center",
bootstrap_options = c("striped", "hover", "condensed"))
Tabla de Distribución de frecuencias del nivel de concentración de
nitrato (mg/L) en los países analizados
|
Lim inf
|
Lim sup
|
MC
|
ni
|
hi(%)
|
Ni asc
|
Ni desc
|
Hi asc(%)
|
Hi desc(%)
|
|
0.05
|
4.21
|
2.13
|
235
|
7.83
|
235
|
3000
|
7.83
|
100
|
|
4.21
|
8.37
|
6.29
|
248
|
8.27
|
483
|
2765
|
16.1
|
92.17
|
|
8.37
|
12.54
|
10.455
|
268
|
8.93
|
751
|
2517
|
25.03
|
83.9
|
|
12.54
|
16.7
|
14.62
|
271
|
9.03
|
1022
|
2249
|
34.07
|
74.97
|
|
16.7
|
20.86
|
18.78
|
239
|
7.97
|
1261
|
1978
|
42.03
|
65.93
|
|
20.86
|
25.02
|
22.94
|
251
|
8.37
|
1512
|
1739
|
50.4
|
57.97
|
|
25.02
|
29.18
|
27.1
|
237
|
7.9
|
1749
|
1488
|
58.3
|
49.6
|
|
29.18
|
33.34
|
31.26
|
248
|
8.27
|
1997
|
1251
|
66.57
|
41.7
|
|
33.34
|
37.51
|
35.425
|
234
|
7.8
|
2231
|
1003
|
74.37
|
33.43
|
|
37.51
|
41.67
|
39.59
|
250
|
8.33
|
2481
|
769
|
82.7
|
25.63
|
|
41.67
|
45.83
|
43.75
|
252
|
8.4
|
2733
|
519
|
91.1
|
17.3
|
|
45.83
|
49.99
|
47.91
|
267
|
8.9
|
3000
|
267
|
100
|
8.9
|
|
TOTAL
|
|
|
3000
|
100
|
|
|
|
|
# GRAFICAS
# Histograma
hist(Nivel_nitrato, breaks = 11,
main = "Gráfica N°1: Distribución de la concentración de nitratos (mg/L)
presente en el estudio sobre contaminación del agua y enfermedades asociadas",
xlab = "Niveles de Nitrato (mg/L)",
ylab = "Cantidad",
ylim = c(0, max(ni)),
col = "purple",
cex.main = 0.9,
cex.lab = 1,
cex.axis = 0.9,
xaxt = "n")
axis(1, at = Hist_Nitrato$breaks,
labels = Hist_Nitrato$breaks, las = 1,
cex.axis = 0.9)
# Cargar librerías necesarias
library(fitdistrplus)
## Warning: package 'fitdistrplus' was built under R version 4.5.1
## Cargando paquete requerido: MASS
## Cargando paquete requerido: survival
# Supón que ya tienes tus datos de nitrato:
# Nivel_nitrato <- tus datos aquí
# Filtrar valores positivos
Nivel_nitrato_pos <- Nivel_nitrato[Nivel_nitrato > 0]
# Puntos del centro de las barras que quieres usar (de la 3 a la 7 aprox)
x_puntos <- TDF_nitrato$MC[3:7]
y_puntos <- TDF_nitrato$ni[3:7] + 10 # un poco más arriba de la parte superior de las barras
# Generamos una secuencia suave para el eje x
x_vals <- seq(min(x_puntos), max(x_puntos), length.out = 200)
# Ajustamos curva suave (loess) a los puntos elegidos
modelo_loess <- loess(y_puntos ~ x_puntos, span = 1) # puedes ajustar span para más o menos suavidad
y_suave <- predict(modelo_loess, newdata = data.frame(x_puntos = x_vals))
# Obtener posiciones reales del gráfico
bp <- barplot(
TDF_nitrato$ni,
space = 0,
col = "purple",
border = "black",
names.arg = TDF_nitrato$MC,
main = "Gráfica N°1: Distribución de la concentración de nitratos (mg/L)\npresente en el estudio sobre contaminación del agua y enfermedades asociadas",
xlab = "Niveles de Nitrato (mg/L)",
ylab = "Cantidad",
ylim = c(0, max(TDF_nitrato$ni) + 50)
)
# Interpolamos las posiciones para el eje x
interp <- approxfun(TDF_nitrato$MC, bp)
x_grafico <- interp(x_vals)
# Dibujamos la curva suavizada
lines(x_grafico, y_suave, col = "gold", lwd = 4)
# Subconjunto de la curva (de la barra 3 a la 7)
Nivel_nitrato_pos <- Nivel_nitrato[Nivel_nitrato > 0]
Histograma <- hist(Nivel_nitrato_pos, breaks = 10, plot = FALSE)
lambda <- 1 / mean(Nivel_nitrato_pos) # parámetro del modelo exponencial
# Frecuencias observadas y esperadas SOLO de la barra 3 a la 7
indices <- 3:7
FO <- Histograma$counts[indices]
FE <- c()
for (i in indices) {
P <- pexp(Histograma$breaks[i + 1], rate = lambda) - pexp(Histograma$breaks[i], rate = lambda)
FE <- c(FE, P * length(Nivel_nitrato_pos))
}
# Test de Pearson
pearson <- cor(FO, FE)
pearson
## [1] 0.9870695
# Test de Chi-cuadrado (aprueba si X2 < chi)
X2 <- sum((FO - FE)^2 / FE)
X2
## [1] 133.0127
gl <- length(FO) - 1 # grados de libertad
chi <- qchisq(0.99, df = gl)
X2 > chi # TRUE = no pasa
## [1] TRUE
# Probabilidad entre 5 y 8 mg/L
P_5_8 <- (pexp(8, rate = lambda) - pexp(5, rate = lambda)) * 100
P_5_8
## [1] 9.236193
# Gráfico con curva y área sombreada entre 5 y 8
x <- seq(0, max(Nivel_nitrato_pos), length.out = 1000)
plot(x, dexp(x, rate = lambda), type = "l", col = "red", lwd = 2,
main = "Distribución Exponencial Ajustada: Nivel de Nitrato (mg/L)",
xlab = "Nivel de Nitrato (mg/L)", ylab = "Densidad de Probabilidad",
ylim = c(0, max(dexp(x, rate = lambda)) * 1.1))
x_area <- seq(5, 8, length.out = 100)
polygon(c(5, x_area, 8), c(0, dexp(x_area, rate = lambda), 0),
col = rgb(0, 0, 1, 0.3), border = NA)
legend("topright",
legend = c("Modelo Exponencial", "Área entre 5 y 8 mg/L"),
col = c("red", rgb(0, 0, 1, 0.3)),
lwd = 2,
pch = c(NA, 15),
pt.cex = 2,
bty = "n")
# Intervalo de confianza (95%) de la media
media <- mean(Nivel_nitrato_pos)
error <- sd(Nivel_nitrato_pos) / sqrt(length(Nivel_nitrato_pos))
lim_inf <- media - 2 * error
lim_sup <- media + 2 * error
# Mostrar resultados
cat("Test de Pearson:", round(pearson, 3), "\n")
## Test de Pearson: 0.987
cat("Test Chi-cuadrado:", round(X2, 2), "vs", round(chi, 2), "→ Aprobado?", X2 < chi, "\n")
## Test Chi-cuadrado: 133.01 vs 13.28 → Aprobado? FALSE
cat("Probabilidad entre 5 y 8 mg/L:", round(P_5_8, 2), "%\n")
## Probabilidad entre 5 y 8 mg/L: 9.24 %
cat("Intervalo de Confianza (95%) de la media:", round(lim_inf, 2), "-", round(lim_sup, 2), "mg/L\n")
## Intervalo de Confianza (95%) de la media: 24.55 - 25.61 mg/L
# Parámetro lambda
lambda <- 1 / mean(Nivel_nitrato)
# Rango de interés
a <- 3
b <- 5
# Probabilidad de que Nivel_nitrato esté entre a y b
P1 <- (pexp(b, rate = lambda) - pexp(a, rate = lambda)) * 100
P1 # Esto da la probabilidad en porcentaje
## [1] 6.800619
# Estadísticos
media <- mean(Nivel_nitrato)
desviacion <- sd(Nivel_nitrato)
n <- length(Nivel_nitrato)
# Intervalo 95%
error <- 2 * (desviacion / sqrt(n))
lim_inf <- round(media - error, 3)
lim_sup <- round(media + error, 3)
plot(1, type = "n", axes = FALSE, xlab = "", ylab = "") # Gráfico vacío para texto
text(x = 1, y = 1,
labels = paste(
"Cálculos de Probabilidad y Estimación de Parámetros\n",
"\n• Probabilidad de que el nivel de nitrato esté entre 3 y 5 mg/L:",
sprintf("%.2f %%", P1),
"\n\n• Intervalo de confianza del 95%% para la media poblacional del
Nivel de Nitrato:",
sprintf("\nEntre %.2f mg/L y %.2f mg/L", lim_inf, lim_sup),
"\n\nInterpretación:\nCon un 95%% de confianza, se espera que la media
real del nivel de nitrato en los países analizados esté dentro de este
rango.\nTambién, hay aproximadamente un",
sprintf(" %.2f %%", P1),
"de probabilidad de que un valor observado esté entre 3 y 5 mg/L, según
el modelo exponencial ajustado."),
cex = 1.3,
col = "blue",
font = 6)
# Crear un resumen del test de bondad de ajuste
resumen_test <- data.frame(
Test = c("Coeficiente de Pearson", "Chi-cuadrado"),
Estadístico = c(round(pearson, 3), round(X2, 2)),
Valor_critico = c(NA, round(chi, 2)),
Aprobado = c("Sí", ifelse(X2 < chi, "Sí", "No"))
)
# Mostrar resumen con kable
library(knitr)
library(kableExtra)
kable(resumen_test,
caption = "Resumen del Test de Bondad de Ajuste para el Modelo Exponencial") %>%
kable_styling(bootstrap_options = c("striped", "hover"),
full_width = FALSE, position = "center")
Resumen del Test de Bondad de Ajuste para el Modelo Exponencial
|
Test
|
Estadístico
|
Valor_critico
|
Aprobado
|
|
Coeficiente de Pearson
|
0.987
|
NA
|
Sí
|
|
Chi-cuadrado
|
133.010
|
13.28
|
No
|