Simular el teorema del limite central con un conjunto de datos.
Con una base de datos y librerías adecuadas, simular el valor de la media muestral comparado con el valor de la media poblacional asociando con ello con el teorema de limite central.
library(dplyr)
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
library(mosaic)
## Registered S3 method overwritten by 'mosaic':
## method from
## fortify.SpatialPolygonsDataFrame ggplot2
##
## The 'mosaic' package masks several functions from core packages in order to add
## additional features. The original behavior of these functions should not be affected by this.
##
## Attaching package: 'mosaic'
## The following object is masked from 'package:Matrix':
##
## mean
## The following object is masked from 'package:ggplot2':
##
## stat
## The following objects are masked from 'package:dplyr':
##
## count, do, tally
## The following objects are masked from 'package:stats':
##
## binom.test, cor, cor.test, cov, fivenum, IQR, median, prop.test,
## quantile, sd, t.test, var
## The following objects are masked from 'package:base':
##
## max, mean, min, prod, range, sample, sum
library(readr)
library(ggplot2)
library(knitr)
library(fdth)
##
## Attaching package: 'fdth'
## The following objects are masked from 'package:mosaic':
##
## sd, var
## The following objects are masked from 'package:stats':
##
## sd, var
library(gtools)
##
## Attaching package: 'gtools'
## The following object is masked from 'package:mosaic':
##
## logit
library(Rmpfr)
## Loading required package: gmp
##
## Attaching package: 'gmp'
## The following object is masked from 'package:mosaic':
##
## factorize
## The following objects are masked from 'package:Matrix':
##
## crossprod, tcrossprod
## The following objects are masked from 'package:base':
##
## %*%, apply, crossprod, matrix, tcrossprod
## C code of R package 'Rmpfr': GMP using 64 bits per limb
##
## Attaching package: 'Rmpfr'
## The following object is masked from 'package:gmp':
##
## outer
## The following objects are masked from 'package:stats':
##
## dbinom, dgamma, dnbinom, dnorm, dpois, pnorm
## The following objects are masked from 'package:base':
##
## cbind, pmax, pmin, rbind
N <- 1000000;
edad.poblacion <- round(rnorm(N, mean = 35, sd = 5), 0)
summary(edad.poblacion)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 11.00 32.00 35.00 35.01 38.00 60.00
paste("El valor de edad de una pobacion. Los primeros cincuenta valores ...")
## [1] "El valor de edad de una pobacion. Los primeros cincuenta valores ..."
head(edad.poblacion, 50)
## [1] 33 27 33 46 37 40 30 44 29 39 27 39 37 38 33 38 36 35 43 29 30 27 38 30 32
## [26] 34 37 31 33 34 38 39 37 39 39 31 33 37 33 47 26 24 43 31 29 39 37 40 33 34
paste("El valor de edad de una pobacion. Los últimos cincuenta valores ...")
## [1] "El valor de edad de una pobacion. Los últimos cincuenta valores ..."
tail(edad.poblacion, 50)
## [1] 32 35 28 32 30 40 42 39 25 33 30 30 36 36 46 35 38 37 30 35 38 40 39 38 37
## [26] 30 38 31 32 34 37 34 35 39 34 39 42 35 35 35 33 43 36 44 39 32 33 35 33 40
media.pob <- mean(edad.poblacion)
desv.std <- sd(edad.poblacion)
paste("Los parámetros de la media y desviación estándard de la población")
## [1] "Los parámetros de la media y desviación estándard de la población"
media.pob; desv.std
## [1] 35.00537
## [1] 5.010417
n <- 500
muestras <- data.frame(m1=sample(edad.poblacion, n),
m2=sample(edad.poblacion, n),
m3=sample(edad.poblacion, n),
m4=sample(edad.poblacion, n),
m5=sample(edad.poblacion, n))
summary(muestras)
## m1 m2 m3 m4 m5
## Min. :21.00 Min. :19.00 Min. :19.00 Min. :19.00 Min. :18.0
## 1st Qu.:32.00 1st Qu.:32.00 1st Qu.:31.00 1st Qu.:31.00 1st Qu.:32.0
## Median :35.00 Median :35.50 Median :35.00 Median :35.00 Median :35.0
## Mean :34.97 Mean :35.35 Mean :34.93 Mean :34.68 Mean :35.2
## 3rd Qu.:38.00 3rd Qu.:39.00 3rd Qu.:38.00 3rd Qu.:38.00 3rd Qu.:39.0
## Max. :48.00 Max. :50.00 Max. :51.00 Max. :52.00 Max. :49.0
kable(head(muestras, 10), caption = "Muestras de la población. Los primeros diez de 500 registros")
| m1 | m2 | m3 | m4 | m5 |
|---|---|---|---|---|
| 33 | 42 | 40 | 32 | 36 |
| 40 | 33 | 30 | 29 | 36 |
| 37 | 34 | 35 | 33 | 41 |
| 31 | 40 | 36 | 33 | 34 |
| 25 | 38 | 37 | 36 | 33 |
| 41 | 40 | 37 | 37 | 29 |
| 33 | 35 | 29 | 33 | 43 |
| 33 | 37 | 37 | 29 | 38 |
| 28 | 34 | 21 | 19 | 29 |
| 32 | 40 | 32 | 32 | 35 |
kable(head(muestras, 10), caption = "Muestras de la población. Los últimos diez de 500 registros")
| m1 | m2 | m3 | m4 | m5 |
|---|---|---|---|---|
| 33 | 42 | 40 | 32 | 36 |
| 40 | 33 | 30 | 29 | 36 |
| 37 | 34 | 35 | 33 | 41 |
| 31 | 40 | 36 | 33 | 34 |
| 25 | 38 | 37 | 36 | 33 |
| 41 | 40 | 37 | 37 | 29 |
| 33 | 35 | 29 | 33 | 43 |
| 33 | 37 | 37 | 29 | 38 |
| 28 | 34 | 21 | 19 | 29 |
| 32 | 40 | 32 | 32 | 35 |
medias <- 0
error <- 0
for(i in 1:5) {
medias[i] <- mean(muestras[,i])
error[i] <- medias[i] - media.pob
}
error.muestreo <- data.frame(Media.Poblacion = media.pob, Media.Muestras = medias, Errores = error)
kable(error.muestreo, caption = "Error de media de edad de cada muestra con respecto a la media de la población")
| Media.Poblacion | Media.Muestras | Errores |
|---|---|---|
| 35.00537 | 34.968 | -0.03737 |
| 35.00537 | 35.350 | 0.34463 |
| 35.00537 | 34.928 | -0.07737 |
| 35.00537 | 34.680 | -0.32537 |
| 35.00537 | 35.198 | 0.19263 |
hist(edad.poblacion, main = "Histrograma de la edad de la población")
hist(muestras$m1, main = "Histrograma de la edad de la muestra 1", ylab = "Edades", xlab="Observacaiones")
hist(muestras$m2, main = "Histrograma de la edad de la muestra 2", ylab = "Edades", xlab="Observacaiones")
hist(muestras$m3, main = "Histrograma de la edad de la muestra 3", ylab = "Edades", xlab="Observacaiones")
hist(muestras$m4, main = "Histrograma de la edad de la muestra 4", ylab = "Edades", xlab="Observacaiones")
hist(muestras$m5, main = "Histrograma de la edad de la muestra 5", ylab = "Edades", xlab="Observacaiones")
Con el ejercicio anterior, anterior se encontró el error de muestreo y se presentaron los resultados de comparar un estadístico para una muestra (como la media de la muestra) con la media de la población; bajo este contexto, cuando se usa la media muestral para estudiar la media de la población, ¿cómo se determina la exactitud de la estimación?, es decir, como saber si la media de la muestra es un estimador real con respecto a la población.
Se recapitulan los datos iniciales:
options(scipen = 999) # Para mostrar notación normal y no científica en el valor de N: 1e+06
N; n;
## [1] 1000000
## [1] 500
options(scipen = 0) # Regresa a notación numérica normal
summary(edad.poblacion)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 11.00 32.00 35.00 35.01 38.00 60.00
media.pob; desv.std
## [1] 35.00537
## [1] 5.010417
\[C(N/(n))=N!/n!⋅(N−n)!\] Por lo anterior el número de muestra con una población de un millón (1000000) en grupo de 500 es demasiado grande para tratarlo.
Para ejemplificar el teorema de límite central se reduce la población a 10 con muestras de 2 personas.
Simulando una nueva población con los mismos valores de media de edad igual a 35 y desviación de 5.
Se determinan los nuevos parámetros de medias y desviaciones estándar de la población. * N = 10, tamaño de la población * n = 2, tamaño de la muestra
N <- 10; n <- 2
edad.poblacion <- round(rnorm(N, mean = 35, sd = 5), 0)
edad.poblacion
## [1] 39 35 34 49 30 34 37 36 34 23
media.pob <- mean(edad.poblacion)
desv.std <- sd(edad.poblacion)
media.pob; desv.std
## [1] 35.1
## [1] 6.573516
De acuerdo al concepto distribución muestral de la media ¿cuál es el número de muestras que hay que determinar en grupos de 2 para una población de 10?
n.combinaciones <- factorialMpfr(N) / (factorialMpfr(n) * (factorialMpfr(N-n)))
as.integer(n.combinaciones)
## [1] 45
Ahora bien ¿cuál es el valor estadístico de la media de la edad de la primera muestra, de la segunda, de la tercera y de la 45 ava muestra.
muestras <- cbind(1:as.integer(n.combinaciones))
muestras <- cbind(muestras, combinations(N, n, 1:N))
muestras <- cbind(muestras, edad.poblacion[muestras[,2]], edad.poblacion[muestras[,3]])
medias <- 0
error <- 0
for(i in 1:as.integer(n.combinaciones)) {
medias[i] <- mean(muestras[i,c(4,5)])
error[i] <- medias[i] - media.pob
}
muestras <- cbind(muestras, medias)
muestras <- cbind(muestras, media.pob)
muestras <- cbind(muestras, error)
muestras <- data.frame(muestras)
colnames(muestras) <- c("Muestra", "Pos.1", "Pos.2", "Valor.1", "Valor.2", "Media muestra", "Media pob.", "Error")
kable(muestras, caption = "Las muestras")
| Muestra | Pos.1 | Pos.2 | Valor.1 | Valor.2 | Media muestra | Media pob. | Error |
|---|---|---|---|---|---|---|---|
| 1 | 1 | 2 | 39 | 35 | 37.0 | 35.1 | 1.9 |
| 2 | 1 | 3 | 39 | 34 | 36.5 | 35.1 | 1.4 |
| 3 | 1 | 4 | 39 | 49 | 44.0 | 35.1 | 8.9 |
| 4 | 1 | 5 | 39 | 30 | 34.5 | 35.1 | -0.6 |
| 5 | 1 | 6 | 39 | 34 | 36.5 | 35.1 | 1.4 |
| 6 | 1 | 7 | 39 | 37 | 38.0 | 35.1 | 2.9 |
| 7 | 1 | 8 | 39 | 36 | 37.5 | 35.1 | 2.4 |
| 8 | 1 | 9 | 39 | 34 | 36.5 | 35.1 | 1.4 |
| 9 | 1 | 10 | 39 | 23 | 31.0 | 35.1 | -4.1 |
| 10 | 2 | 3 | 35 | 34 | 34.5 | 35.1 | -0.6 |
| 11 | 2 | 4 | 35 | 49 | 42.0 | 35.1 | 6.9 |
| 12 | 2 | 5 | 35 | 30 | 32.5 | 35.1 | -2.6 |
| 13 | 2 | 6 | 35 | 34 | 34.5 | 35.1 | -0.6 |
| 14 | 2 | 7 | 35 | 37 | 36.0 | 35.1 | 0.9 |
| 15 | 2 | 8 | 35 | 36 | 35.5 | 35.1 | 0.4 |
| 16 | 2 | 9 | 35 | 34 | 34.5 | 35.1 | -0.6 |
| 17 | 2 | 10 | 35 | 23 | 29.0 | 35.1 | -6.1 |
| 18 | 3 | 4 | 34 | 49 | 41.5 | 35.1 | 6.4 |
| 19 | 3 | 5 | 34 | 30 | 32.0 | 35.1 | -3.1 |
| 20 | 3 | 6 | 34 | 34 | 34.0 | 35.1 | -1.1 |
| 21 | 3 | 7 | 34 | 37 | 35.5 | 35.1 | 0.4 |
| 22 | 3 | 8 | 34 | 36 | 35.0 | 35.1 | -0.1 |
| 23 | 3 | 9 | 34 | 34 | 34.0 | 35.1 | -1.1 |
| 24 | 3 | 10 | 34 | 23 | 28.5 | 35.1 | -6.6 |
| 25 | 4 | 5 | 49 | 30 | 39.5 | 35.1 | 4.4 |
| 26 | 4 | 6 | 49 | 34 | 41.5 | 35.1 | 6.4 |
| 27 | 4 | 7 | 49 | 37 | 43.0 | 35.1 | 7.9 |
| 28 | 4 | 8 | 49 | 36 | 42.5 | 35.1 | 7.4 |
| 29 | 4 | 9 | 49 | 34 | 41.5 | 35.1 | 6.4 |
| 30 | 4 | 10 | 49 | 23 | 36.0 | 35.1 | 0.9 |
| 31 | 5 | 6 | 30 | 34 | 32.0 | 35.1 | -3.1 |
| 32 | 5 | 7 | 30 | 37 | 33.5 | 35.1 | -1.6 |
| 33 | 5 | 8 | 30 | 36 | 33.0 | 35.1 | -2.1 |
| 34 | 5 | 9 | 30 | 34 | 32.0 | 35.1 | -3.1 |
| 35 | 5 | 10 | 30 | 23 | 26.5 | 35.1 | -8.6 |
| 36 | 6 | 7 | 34 | 37 | 35.5 | 35.1 | 0.4 |
| 37 | 6 | 8 | 34 | 36 | 35.0 | 35.1 | -0.1 |
| 38 | 6 | 9 | 34 | 34 | 34.0 | 35.1 | -1.1 |
| 39 | 6 | 10 | 34 | 23 | 28.5 | 35.1 | -6.6 |
| 40 | 7 | 8 | 37 | 36 | 36.5 | 35.1 | 1.4 |
| 41 | 7 | 9 | 37 | 34 | 35.5 | 35.1 | 0.4 |
| 42 | 7 | 10 | 37 | 23 | 30.0 | 35.1 | -5.1 |
| 43 | 8 | 9 | 36 | 34 | 35.0 | 35.1 | -0.1 |
| 44 | 8 | 10 | 36 | 23 | 29.5 | 35.1 | -5.6 |
| 45 | 9 | 10 | 34 | 23 | 28.5 | 35.1 | -6.6 |
La media de la distribución muestral de la media se obtiene al sumar las medias muestrales y dividir el resultado entre el número de muestras. La media de todas las medias muestrales se representa mediante μx¯
Entonces … la media de la distribución muestral comparado con la media poblaciónal
paste("La media poblacional es: ", media.pob, " y la media de la edad de la distribución muestral es: ", mean(muestras$`Media muestra`))
## [1] "La media poblacional es: 35.1 y la media de la edad de la distribución muestral es: 35.1"