library(ggplot2)
library(dplyr)
my_tamamu = 1e5 #tamaño muestral

Motivación: las 3 principales distribuciones

La distribución binomial: lanzar un dado 90 veces

Todo experimento consiste en realizar una actividad y, en el curso de ella, medir “algo”. Si repetimos el experimento muchas veces (o si lo realizamos sobre muchos individuos) vamos obteniendo una muestra de resultados. A este número de veces que se repite el experimento se le llama “tamaño muestral”, que es nuestro caso es 10^{5}.

En este caso:

Actividad: lanzar un dado 90 veces.
Medida: contar cuantas veces sale el “6”. Cómo mínimo, serán 0 veces y, como máximo, 90 veces.
Repeticiones del experimento: 10^{5} veces.

La función “rbinom” realiza todo esto por nosotros. Sus parámetros son:

size: número de veces que lanzamos el dado, que es 90.
prob: probabilidad teórica de que salga el “6” al lanzar el dado una sóla vez, que es 1/6 porque el dado tiene 6 caras.
n: número de veces que repetimos el experimento, que es 10^{5}.

my_n = my_tamamu
my_size = 90
my_prob = 1/6
muestra <- rbinom(n=my_n, size=my_size, prob=my_prob)

La muestra de resultados es:

head(muestra, 30)

##  [1] 20 17 14 12 17 14 20 15 15 14 19 16 12 19 16 15 11 15 12 22 18 16 12
## [24] 19 21 13 11 16 12 13

ggplot() + 
    geom_histogram(aes(x = muestra, y = ..density..),
                   color = "white", 
                   fill = "orange", 
                   binwidth = 1) +
    labs(title = "Histograma")

Podemos recoger en una tabla cuantas veces ha salido cada resultado:

df <- data.frame(nro_seises = muestra) %>% 
    group_by(nro_seises) %>% 
    summarise(veces = n()) %>%
    arrange(-veces)
df

## # A tibble: 31 × 2
##    nro_seises veces
##         <int> <int>
## 1          15 11217
## 2          14 10982
## 3          16 10457
## 4          13 10241
## 5          17  8941
## 6          12  8416
## 7          18  7656
## 8          11  6350
## 9          19  5516
## 10         10  4399
## # ... with 21 more rows

Vemos que el resultado más repetido es 15, que ha ocurrido 11217 veces. Como el máximo de veces que podía repetirse era 10^{5}, esto supone una proporción de 0.1122 o un porcentaje del 11.22%. A esta proporción se le llama frecuencia del resultado 15.

Podemos recoger en una tabla las frecuencias de todos los resultados:

distri_bin <- mutate(df, frecuencia = round(veces / my_n,4))
distri_bin <- select(distri_bin, nro_seises, frecuencia)
distri_bin

## # A tibble: 31 × 2
##    nro_seises frecuencia
##         <int>      <dbl>
## 1          15     0.1122
## 2          14     0.1098
## 3          16     0.1046
## 4          13     0.1024
## 5          17     0.0894
## 6          12     0.0842
## 7          18     0.0766
## 8          11     0.0635
## 9          19     0.0552
## 10         10     0.0440
## # ... with 21 more rows

Abstracción

Cualquier experimento que se haga para medir “algo” acabará siempre como el anterior, es decir, con una muestra de resultados obtenidos y una tabla de resultados y frecuencias.

Si la tabla de resultados y frecuencias coincide con la que acabamos de obtener (mismos resultados y misma frecuencia para cada resultado) diremos que el experimento indica que lo que se está midiendo sigue una distribución binomial de tamaño 90 y probabilidad 1/6.

Por ejemplo, si medimos la edad de las personas ingresadas en un hospital podría ser que saliera la tabla anterior (los niños de 15 años serían, en particular, los más habituales). Diríamos que la edad de las personas internadas en ese hospital sigue una distribución binomial de tamaño 90 y probabilidad 1/6.

Obviamente, podemos decir que el número de “6s”" que salen al tirar 90 veces un dado sigue una distribución binomial de tamaño 90 y probabilidad 1/6.

Otros ejemplos serían:

El número de “caras” que salen al lanzar una moneda 100 veces sigue una distribución binomial de tamaño 100 y probabilidad 1/2.
El número de “ases de corazones” que salen al sacar una carta 100 veces de una baraja de poker (la carta se devuelve tras cada extracción y el mazo se mezcla) sigue una distribución binomial de tamaño 100 y probabilidad 1/52.

La distribución normal: medir una longitud de 5m con un aparato cuyo margen de error es 25cm

En este caso:

Actividad: Medir una barra cuya longitud sabemos que es 5 metros con un aparato cuya precisión es de 25 cm, es decir, 0.25 metros.
Medida: la implícita en la actividad.
Repeticiones del experimento:10^{5} veces.

La función “rnorm” realiza el experimento por nosotros. Sus parámetros son:

mean: la longitud verdadera de la barra, que es 5
sd: el margen de error del aparato de medida, que es 0.25
n: número de veces que repetimos el proceso de medición, que es 10^{5}.

my_n = my_tamamu
my_mean = 5
my_sd = 0.25
muestra <- round( rnorm(n=my_n, mean=my_mean, sd=my_sd) ,2)

La muestra de resultados:

head(muestra, 30)

##  [1] 4.83 5.27 5.04 4.54 5.03 4.70 5.13 5.08 4.74 4.81 5.12 5.07 4.87 4.95
## [15] 4.96 4.48 5.46 5.05 4.50 4.87 4.79 4.73 4.83 5.08 5.25 5.36 5.52 5.10
## [29] 5.14 4.64

ggplot() + 
    geom_histogram(aes(x = muestra, y = ..density..),
                   color = "white", 
                   fill = "orange", 
                   binwidth = 0.01) +
    labs(title = "Histograma")

Podemos recoger en una tabla cuantas veces ha salido cada resultado:

df <- data.frame(medidas = muestra) %>% 
    group_by(medidas) %>% 
    summarise(veces = n()) %>%
    arrange(-veces)
df

## # A tibble: 197 × 2
##    medidas veces
##      <dbl> <int>
## 1     5.00  1678
## 2     5.01  1647
## 3     5.02  1637
## 4     4.96  1603
## 5     4.95  1599
## 6     4.99  1596
## 7     5.05  1584
## 8     4.97  1574
## 9     5.04  1560
## 10    4.92  1556
## # ... with 187 more rows

Vemos que el resultado más repetido es 5, que ha ocurrido 1678 veces. Como el máximo de veces que podía repetirse era 10^{5}, esto supone una proporción de 0.0168 o un porcentaje del 1.68%. A esta proporción se le llama frecuencia del resultado 5.

Podemos recoger en una tabla las frecuencias de todas los resultados:

distri_nor <- mutate(df, frecuencia = round(veces / my_n,4))
distri_nor <- select(distri_nor, medidas, frecuencia)
distri_nor

## # A tibble: 197 × 2
##    medidas frecuencia
##      <dbl>      <dbl>
## 1     5.00     0.0168
## 2     5.01     0.0165
## 3     5.02     0.0164
## 4     4.96     0.0160
## 5     4.95     0.0160
## 6     4.99     0.0160
## 7     5.05     0.0158
## 8     4.97     0.0157
## 9     5.04     0.0156
## 10    4.92     0.0156
## # ... with 187 more rows

Abstracción

Cualquier experimento que se haga para medir “algo” acabará siempre como el anterior, es decir, con una muestra de resultados obtenidos y una tabla de resultados y frecuencias.

Por ejemplo, si medimos los pesos de niños de 3 meses de edad que han pasado por una clínica pediátrica, podría ser que saliera la tabla anterior (en particular, el peso más repetido sería 5 kg). En ese caso diríamos que los pesos de los niños de 3 meses de edad que han pasado por esa clínica pediátrica siguen una distribución normal de media 5 kg y desviación típica 0.25 kg.

Obviamente, podemos decir que las medidas de un barra de 5 metros con una aparato cuya precisión es de 0.25 metros siguen una distribución normal de media 5 metros y desviación típica 0.25 metros.

Otros ejemplos serían:

Cualquier medida de características naturales (altura, coeficiente de inteligencia, tensión arterial…) suele seguir una distribución normal.

La distribución de Poisson: esperar 2 horas un autobús que suele pasar 4 veces cada hora

En este caso:

Actividad: Esperar 2 horas en una parada de un autobús que tiene una frecuencia habitual de 4 autobuses por hora.
Medida: contar cuántos autobuses realmente pasan.
Repeticiones del experimento: 10^{5} veces.

La función “rpois” realiza el experimento por nosotros. Sus parámetros son:

\(\lambda\): resultado de multiplicar la frecuencia horaria habitual por el número de horas que vamos a esperar, es decir: 4 * 2 = 8.
n: número de veces que vamos a esperar 2 horas contando autobuses, que es 10^{5}.

my_ = my_tamamu
my_lambda = 4 * 2
muestra <- rpois(n=my_n, lambda=my_lambda)

La muestra de resultados es:

head(muestra, 30)

##  [1]  4  5  5 10 10  5  5 10  8  9  7  7  7  2 11  5  5  3 10  9  4  8  6
## [24] 10 13 10  6 10  8  5

ggplot() + 
    geom_histogram(aes(x = muestra, y = ..density..),
                   color = "white", 
                   fill = "orange", 
                   binwidth = 1) +
    labs(title = "Histograma")

Podemos recoger en una tabla cuantas veces ha salido cada resultado:

df <- data.frame(nro_autobuses = muestra) %>% 
    group_by(nro_autobuses) %>% 
    summarise(veces = n()) %>%
    arrange(-veces)
df

## # A tibble: 24 × 2
##    nro_autobuses veces
##            <int> <int>
## 1              8 13861
## 2              7 13812
## 3              9 12392
## 4              6 12325
## 5             10  9964
## 6              5  9221
## 7             11  7243
## 8              4  5824
## 9             12  4813
## 10            13  2946
## # ... with 14 more rows

Vemos que el resultado más repetido es 8, que ha ocurrido 13861 veces. Como el máximo de veces que podía repetirse era 10^{5}, esto supone una proporción de 0.1386 o un porcentaje del 13.86%. A esta proporción se le llama frecuencia del resultado 8.

Podemos recoger en una tabla las frecuencias de todos los resultados:

distri_poi <- mutate(df, frecuencia = round(veces / my_n,4))
distri_poi <- select(distri_poi, nro_autobuses, frecuencia)
distri_poi

## # A tibble: 24 × 2
##    nro_autobuses frecuencia
##            <int>      <dbl>
## 1              8     0.1386
## 2              7     0.1381
## 3              9     0.1239
## 4              6     0.1232
## 5             10     0.0996
## 6              5     0.0922
## 7             11     0.0724
## 8              4     0.0582
## 9             12     0.0481
## 10            13     0.0295
## # ... with 14 more rows

Abstracción

Cualquier experimento que se haga para medir “algo” acabará siempre como el anterior, es decir, con una muestra de resultados obtenidos y una tabla de resultados y frecuencias.

Por ejemplo, si medimos la edad de las personas ingresadas en un hospital podría ser que saliera la tabla anterior (los niños de 8 años serían, en particular, los más habituales). Diríamos que la edad de las personas internadas en ese hospital sigue una distribución de Poisson de parámetro \(\lambda\) = 4*2 = 8.

Obviamente, podemos decir que el recuento durante 2 horas de los autobuses que llegan a una parada cuya frecuencia habitual es de 4 autobuses por hora sigue una distribución de Poisson de parámetro \(\lambda\) = 4*2 = 8.

Otros ejemplos serían:

El número de usuarios que accederán a una web de 9:00 a 17:00, si lo habitual es 15 usuarios por hora, sigue una distribución de Poisson de parámetro \(\lambda\) = 15 * 8 = 120.
El número de errores que cometerá una mecanógrafa en las próximas 15 páginas que escriba, si habitualmente comete 2 errores cada 5 páginas, sigue una distribución de Poisson de parámetro \(\lambda\) = 2/5 * 15 = 6.
El número de camisetas que vende una tienda, si lo habitual es 12 camisetas cada día, sigue una distribución de Poisson de parámetro \(\lambda\) = 15*1 = 15.

La ley de los grandes números

Lo que se puede ver en los experimentos anteriores es lo que dice la Ley de los grandes números: a partir de un tamaño muestral, posiblemente muy grande, es decir, a partir de cierto número de repeticiones del experimento, las tablas de resultados se estabilizan; mismos resultados y mismas frecuencias.

En efecto, rehaciendo este documento una y otra vez se ve que los resultados son prácticamente los mismos y las frecuencias coinciden hasta el tercer o cuarto decimal, apróximadamente.

distri_bin

## # A tibble: 31 × 2
##    nro_seises frecuencia
##         <int>      <dbl>
## 1          15     0.1122
## 2          14     0.1098
## 3          16     0.1046
## 4          13     0.1024
## 5          17     0.0894
## 6          12     0.0842
## 7          18     0.0766
## 8          11     0.0635
## 9          19     0.0552
## 10         10     0.0440
## # ... with 21 more rows

distri_nor

## # A tibble: 197 × 2
##    medidas frecuencia
##      <dbl>      <dbl>
## 1     5.00     0.0168
## 2     5.01     0.0165
## 3     5.02     0.0164
## 4     4.96     0.0160
## 5     4.95     0.0160
## 6     4.99     0.0160
## 7     5.05     0.0158
## 8     4.97     0.0157
## 9     5.04     0.0156
## 10    4.92     0.0156
## # ... with 187 more rows

distri_poi

## # A tibble: 24 × 2
##    nro_autobuses frecuencia
##            <int>      <dbl>
## 1              8     0.1386
## 2              7     0.1381
## 3              9     0.1239
## 4              6     0.1232
## 5             10     0.0996
## 6              5     0.0922
## 7             11     0.0724
## 8              4     0.0582
## 9             12     0.0481
## 10            13     0.0295
## # ... with 14 more rows

Esto se interpreta como que existen unos resultados teóricos del experimento y una frecuencia teórica para cada resultado. R los tiene almacenados en funciones. Por desgracia, estas funciones lo que nos dan es la frecuencia “acumulada” y es necesario algo de trabajo para extraer la frecuencia “a secas”.

Frecuencias acumuladas y percentiles de la distribución binomial (funciones pbinom y qbinom)

La función pbinom nos da la frecuencia acumulada, es decir, la suma de las frecuencias de todos los valores hasta llegar al que le indiquemos (incluido). Por ejemplo, para nuestra binomial de tamaño 90 y probabilidad 1/6, podemos preguntarnos cuánto suman las frecuencias hasta llegar al valor más repetido (el 15):

pbinom(15, size=90, prob=1/6)

## [1] 0.5683735

La frecuencia “a secas” del valor 15 se obtendría restando a su frecuencia acumulada la del valor anterior:

pbinom(15, size=90, prob=1/6) - pbinom(15 - 1, size=90, prob=1/6)

## [1] 0.1121921

Es fácil construir una función fbinom que automatice esto:

fbinom <- function(valor, size, prob) {
    pbinom(valor, size=size, prob=prob) - 
        pbinom(valor - 1, size=size, prob=prob)
}

fbinom(15, size=90, prob=1/6)

## [1] 0.1121921

Veamos como suben las frecuencias acumuladas hasta llegar al valor más repetido (15):

valores <- seq(12,15)
frecuencias_acumuladas <- pbinom(valores, size=90, prob=1/6)
binomial <- data.frame(valores, frecuencias_acumuladas)
binomial

##   valores frecuencias_acumuladas
## 1      12              0.2448148
## 2      13              0.3454656
## 3      14              0.4561814
## 4      15              0.5683735

En el valor 15 superamos, entonces, por primera vez la frecuencia acumulada del 50% (frecuencia acumulada, 0.5). Esta misma información (en qué valor igualamos o superamos la frecuencia acumulada 0.5) nos la da directamente la función qbinom. A este valor se le llama percentil 50:

percentil_50 <- qbinom(0.5, size=90, prob=1/6)
percentil_50

## [1] 15

Por construcción, la suma de todas las frecuencias tiene que ser 1. Entonces, una forma de saber en qué valor empieza a acumularse algo de frecuencia y en qué valor ya está casi toda acumulada, es preguntar a qbinom en qué valor igualamos o superamos las frecuencias acumuladas 10^{-5} y 0.99999. Estos valores se llaman, respectivamente, percentil mínimo y percentil máximo:

percentil_min <- qbinom(frecuencia_min, size=90, prob=1/6)
percentil_min

## [1] 2

percentil_max <- qbinom(frecuencia_max, size=90, prob=1/6)
percentil_max

## [1] 32

Sabiendo el rango de valores que manejamos, podemos construir una tabla de valores y frecuencias práctica. A esta tabla es a lo que se llama, en la práctica, la distribución binomial de tamaño 90 y probabilidad 1/6.

valores <- round(seq(percentil_min, percentil_max, by=1) ,0)
frecuencias <- fbinom(valores, size=90, prob=1/6)
frecuencias_acumuladas <- pbinom(valores, size=90, prob=1/6)
frecuencias <- round(frecuencias, cifras)
frecuencias_acumuladas <- round(frecuencias_acumuladas, cifras)
distri_bin <- data.frame(valores, frecuencias, frecuencias_acumuladas)

curva_bin <- ggplot(data = distri_bin, 
                    aes(x=valores, y=frecuencias)) + 
    geom_line(color = "blue", size=1) + 
    labs(title = "Distribución binomial de tamaño 90 y probabilidad 1/6")
curva_bin

Finalmente, una muestra “perfecta” de resultados que daría justo las frecuencias teóricas, se puede construir:

muestra_bin <- rep(distri_bin$valores, 
                     round(distri_bin$frecuencias * my_tamamu,0))

head(muestra_bin, 200)

##   [1] 2 3 3 3 3 3 3 3 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4 4
##  [36] 4 4 4 4 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5
##  [71] 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5
## [106] 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5
## [141] 5 5 5 5 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6
## [176] 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6

ggplot() + 
    geom_histogram(aes(x = muestra_bin, y = ..density..),
                   color = "white", 
                   fill = "orange", 
                   binwidth = 01) +
    labs(title = "Histograma")

Frecuencias acumuladas y percentiles de la distribución normal (funciones pnorm y qnorm)

La función pnorm nos da la frecuencia acumulada, es decir, la suma de las frecuencias de todos los valores hasta llegar al que le indiquemos (incluido). Por ejemplo, para nuestra normal de media 5 y desviación típica 0.25, podemos preguntarnos cuánto suman las frecuencias hasta llegar al valor más repetido (el 5):

pnorm(5, mean=5, sd=0.25)

## [1] 0.5

La frecuencia “a secas” del valor 5 se obtendría restando a su frecuencia acumulada la del valor anterior. En la distribución normal, el concepto de “valor anterior” depende de la precisión con que medimos los valores. En nuestro caso, ha sido hasta el centímetro:

pnorm(5, mean=5, sd=0.25)- pnorm(5 - 0.01, mean=5, sd=0.25)

## [1] 0.01595344

Es fácil construir una función fnorm que automatice esto:

fnorm <- function(valor, mean, sd) {
    pnorm(valor, mean=mean, sd=sd) - 
        pnorm(valor - 0.01, mean=mean, sd=sd)
}

fnorm(5, mean=5, sd=0.25)

## [1] 0.01595344

Veamos como suben las frecuencias acumuladas hasta llegar al valor más repetido (5):

valores <- seq(4,5, by=0.25)
frecuencias_acumuladas <- round( pnorm(valores, mean=5, sd=0.25) ,4)
normal <- data.frame(valores, frecuencias_acumuladas)
normal

##   valores frecuencias_acumuladas
## 1    4.00                 0.0000
## 2    4.25                 0.0013
## 3    4.50                 0.0228
## 4    4.75                 0.1587
## 5    5.00                 0.5000

En el valor 5 llegamos, entonces, por primera vez a la frecuencia acumulada del 50% (frecuencia acumulada, 0.5). Esta misma información (en qué valor igualamos o superamos la frecuencia acumulada 0.5) nos la da directamente la función qnorm. A este valor se le llama percentil 50:

percentil_50 <- qnorm(0.5, mean=5, sd=0.25)
percentil_50

## [1] 5

Por construcción, la suma de todas las frecuencias tiene que ser 1. Entonces, una forma de saber en qué valor empieza a acumularse algo de frecuencia y en qué valor ya está casi toda acumulada, es preguntar a qnorm en qué valor igualamos o superamos las frecuencias acumuladas 10^{-5} y 0.99999. Estos valores se llaman, respectivamente, percentil mínimo y percentil máximo:

percentil_min <- qnorm(frecuencia_min, mean=5, sd=0.25)
percentil_min

## [1] 3.933777

percentil_max <- qnorm(frecuencia_max, mean=5, sd=0.25)
percentil_max

## [1] 6.066223

Sabiendo el rango de valores que manejamos, podemos construir una tabla de valores y frecuencias práctica. A esta tabla es a lo que se llama, en la práctica, la distribución normal de media 5 y desviación típica 0.25.

valores <- round(seq(percentil_min,percentil_max, by=0.01), 2)
frecuencias <- fnorm(valores, mean=5, sd=0.25)
frecuencias_acumuladas <- pnorm(valores, mean=5, sd=0.25)
frecuencias <- round(frecuencias, cifras)
frecuencias_acumuladas <- round(frecuencias_acumuladas, cifras)
distri_nor <- data.frame(valores, frecuencias, frecuencias_acumuladas)

curva_nor <- ggplot(data = distri_nor, 
                    aes(x=valores, y=frecuencias)) + 
    geom_line(color = "blue", size=1) + 
    labs(title = "Distribución normal de media 5 y desviación típica 0.25")
curva_nor

Finalmente, una muestra “perfecta” de resultados que daría justo las frecuencias teóricas, se puede construir:

muestra_nor <- rep(distri_nor$valores, 
                     round(distri_nor$frecuencias * my_tamamu, 0))

head(muestra_nor, 200)

##   [1] 4.01 4.02 4.03 4.04 4.05 4.06 4.07 4.08 4.08 4.09 4.09 4.10 4.10 4.11
##  [15] 4.11 4.11 4.12 4.12 4.12 4.13 4.13 4.13 4.14 4.14 4.14 4.14 4.15 4.15
##  [29] 4.15 4.15 4.15 4.16 4.16 4.16 4.16 4.16 4.17 4.17 4.17 4.17 4.17 4.17
##  [43] 4.18 4.18 4.18 4.18 4.18 4.18 4.18 4.19 4.19 4.19 4.19 4.19 4.19 4.19
##  [57] 4.19 4.20 4.20 4.20 4.20 4.20 4.20 4.20 4.20 4.20 4.21 4.21 4.21 4.21
##  [71] 4.21 4.21 4.21 4.21 4.21 4.21 4.22 4.22 4.22 4.22 4.22 4.22 4.22 4.22
##  [85] 4.22 4.22 4.22 4.22 4.23 4.23 4.23 4.23 4.23 4.23 4.23 4.23 4.23 4.23
##  [99] 4.23 4.23 4.23 4.24 4.24 4.24 4.24 4.24 4.24 4.24 4.24 4.24 4.24 4.24
## [113] 4.24 4.24 4.24 4.24 4.25 4.25 4.25 4.25 4.25 4.25 4.25 4.25 4.25 4.25
## [127] 4.25 4.25 4.25 4.25 4.25 4.25 4.25 4.26 4.26 4.26 4.26 4.26 4.26 4.26
## [141] 4.26 4.26 4.26 4.26 4.26 4.26 4.26 4.26 4.26 4.26 4.26 4.26 4.27 4.27
## [155] 4.27 4.27 4.27 4.27 4.27 4.27 4.27 4.27 4.27 4.27 4.27 4.27 4.27 4.27
## [169] 4.27 4.27 4.27 4.27 4.27 4.28 4.28 4.28 4.28 4.28 4.28 4.28 4.28 4.28
## [183] 4.28 4.28 4.28 4.28 4.28 4.28 4.28 4.28 4.28 4.28 4.28 4.28 4.28 4.28
## [197] 4.28 4.29 4.29 4.29

ggplot() + 
    geom_histogram(aes(x = muestra_nor, y = ..density..),
                   color = "white", 
                   fill = "orange", 
                   binwidth = 0.01) +
    labs(title = "Histograma")

Frecuencias acumuladas y percentiles de la distribución de Poisson (funciones ppois y qpois)

La función ppois nos da la frecuencia acumulada, es decir, la suma de las frecuencias de todos los valores hasta llegar al que le indiquemos (incluido). Por ejemplo, para nuestra poisson de \(\lambda = 4*2 = 8\), podemos preguntarnos cuánto suman las frecuencias hasta llegar al valor más repetido (el 8):

ppois(8, lambda=4*2)

## [1] 0.5925473

La frecuencia “a secas” del valor 8 se obtendría restando a su frecuencia acumulada la del valor anterior:

ppois(8, lambda=4*2) - ppois(8 - 1, lambda=4*2)

## [1] 0.1395865

Es fácil construir una función fpois que automatice esto:

fpois <- function(valor, lambda) {
    ppois(valor, lambda = lambda) - 
        ppois(valor - 1, lambda = lambda)
}

fpois(8, lambda=4*2)

## [1] 0.1395865

Veamos como suben las frecuencias acumuladas hasta llegar al valor más repetido (8):

valores <- seq(4,8)
frecuencias_acumuladas <- round( ppois(valores, lambda=4*2) ,4)
poisson <- data.frame(valores, frecuencias_acumuladas)
poisson

##   valores frecuencias_acumuladas
## 1       4                 0.0996
## 2       5                 0.1912
## 3       6                 0.3134
## 4       7                 0.4530
## 5       8                 0.5925

En el valor 8 superamos, entonces, por primera vez la frecuencia acumulada del 50% (frecuencia acumulada, 0.5). Esta misma información (en qué valor igualamos o superamos la frecuencia acumulada 0.5) nos la da directamente la función qpois. A este valor se le llama percentil 50:

percentil_50 <- qpois(0.5, lambda=4*2)
percentil_50

## [1] 8

Por construcción, la suma de todas las frecuencias tiene que ser 1. Entonces, una forma de saber en qué valor empieza a acumularse algo de frecuencia y en qué valor ya está casi toda acumulada, es preguntar a qpois en qué valor igualamos o superamos las frecuencias acumuladas 10^{-5} y 0.99999. Estos valores se llaman, respectivamente, percentil mínimo y percentil máximo:

percentil_min <- qpois(frecuencia_min, lambda=4*2)
percentil_min

## [1] 0

percentil_max <- qpois(frecuencia_max, lambda=4*2)
percentil_max

## [1] 23

Sabiendo el rango de valores que manejamos, podemos construir una tabla de valores y frecuencias práctica. A esta tabla es a lo que se llama, en la práctica, la distribución de Poisson de \(\lambda=4*2\).

valores <- round(seq(percentil_min, percentil_max, by=1) ,0)
frecuencias <- fpois(valores, lambda=4*2)
frecuencias_acumuladas <- ppois(valores, lambda=4*2)
frecuencias <- round(frecuencias, cifras)
frecuencias_acumuladas <- round(frecuencias_acumuladas, cifras)
distri_poi <- data.frame(valores, frecuencias, frecuencias_acumuladas)

curva_poi <- ggplot(data = distri_poi, 
                    aes(x=valores, y=frecuencias)) + 
    geom_line(color = "blue", size=1) + 
    labs(title = "Distribución de Poisson de lambda = 4*2")
curva_poi

Finalmente, una muestra “perfecta” de resultados que daría justo las frecuencias teóricas, se puede construir:

muestra_poi <- rep(distri_poi$valores, 
                     round(distri_poi$frecuencias * my_tamamu, 0))

head(muestra_poi, 200)

##   [1] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1
##  [36] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
##  [71] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [106] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [141] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [176] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1

ggplot() + 
    geom_histogram(aes(x = muestra_poi, y = ..density..),
                   color = "white", 
                   fill = "orange", 
                   binwidth = 1) +
    labs(title = "Histograma")

Estadística descriptiva de la distribución binomial

Utilizamos nuestra binomial de tamaño 90 y probabilidad 1/6 en los ejemplos:

pbinom(valor, size=90, prob=1/6)
fbinom(valor, size=90, prob=1/6)
qbinom(valor, size=90, prob=1/6)
distri_bin
curva_bin
muestra_bin

La media (función mean)

La media es un indicador de la tendencia central, es decir, de cuál es el valor más céntrico de la distribución, del cual emanan los demás valores a modo, podríamos decir, de desviaciones:

mean(muestra_bin)

## [1] 14.99978

Se puede demostrar que coincide siempre con: \(size \times prob\).

La mediana (función median)

La mediana es el valor en el que se iguala o supera por primera vez la frecuencia acumulado de 0.5:

mean(muestra_bin)

## [1] 14.99978

Se puede demostrar que coincide siempre la media.

Percentiles (función quantile)

El percentil dd es el valor en el que se iguala o supera por primera vez la frecuencia 0.dd:

quantile(x=muestra_bin, probs=0.05) #percentil 5

## 5% 
##  9

quantile(x=muestra_bin, probs=0.95) #percentil 95

## 95% 
##  21

quantile(x=muestra_bin, probs=0.50) #percentil 50

## 50% 
##  15

El percentil 50 es, por definición, la mediana.

Sumario (función summary y geometría boxplot)

El sumario es un compendio de 6 indicadores:

Mínimo (valor más pequeño con algo de frecuencia)
percentil 25, también llamado “primer cuartil”
percentil 50, también llamado “segundo cuartil” o “mediana”
media
percentil 75, también llamado “tercer cuartil”
Máximo (valor más grande con algo de frecuencia)

summary(muestra_bin)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##       2      13      15      15      17      32

El sumario queda muy bien acompañado de un diagrama de cajas. El diagrama de cajas muestra una caja cuya base es el primer cuartil, la tapa es el tercer cuartil, la mediana (que cae dentro de la caja) se marca claramente y se añaden dos líneas perpendiculares a la caja (llamadas “bigotes de la caja”) con longitud igual a 1,5 veces la altura de la caja.

También se dibujan aquelos valores que quedan fuera de los bigotes de la caja. Estos valores se consideran atípicos (“outliers”).

ggplot(data = NULL, aes(x="", y = muestra_bin)) +
    geom_boxplot()

Histograma (geometría geom_histogram)

El histograma agrupa los valores en intervalos y para cada intervalo levanta una barra proporcional a la frecuencia de los valores que contiene. La anchura del intervalo es calculada automáticamente por R pero se puede sustituir con el parámetro “binwidth”. Es habitual acompañarlo con una línea vertical que indica la media:

ggplot() + 
    geom_histogram(aes(x = muestra_bin, y = ..density..),
                   color = "white", 
                   fill = "orange", binwidth = 1) +
    geom_vline(aes(xintercept = mean(muestra_bin)),
               size=1,
               color="black")

Varianza y desviación típica (funciones var y sd)

La varianza es la media de los cuadrados de las desviaciones con respecto a la media. Como la idea es que la media es “el valor más céntrico de la distribución, del cual emanan los demás valores a modo de desviaciones”, la varianza valora la magnitud de estas desviaciones o, mejor dicho, su cuadrado.

var(muestra_bin)

## [1] 12.49739

Se puede demostrar que coincide siempre con \(size \times prob \times (1 - prob)\).

Como la varianza valora el cuadrado de las desviaciones, es más significativa (porque tiene las unidades de medida de la media) su raiz cuadrada, llamada desviación típica:

sqrt(var(muestra_bin))

## [1] 3.535165

sd(muestra_bin)

## [1] 3.535165

Coeficiente de asimetría (función skewness del paquete moments)

El coeficiente de asimetría indica el grado de simetría que tiene el histograma con respecto al eje que supone la media. “Simetría” quiere decir que los valores a izquierda y derecha del eje son prácticamente los mismos y sus frecuencias también. Otro forma de verlo: que las barras del histograma se “reflejan” en el eje de la media.

library(moments)
skewness(muestra_bin)

## [1] 0.1879125

El coeficiente de asimetría es positivo así que hay asimetría “hacia la derecha”. Esto se ve en el histograma: a la derecha del eje de la media hay más barras. Esto no quiere decir que haya más frecuencia acumulada porque, como la media coincide con la mediana, a ambos lados del eje de la media hay justo un 50% de la frecuencia. Lo que pasa es que hay más valores, más barras. Por eso, son más bajas; de lo contrario, habría más frecuencia acumulada en el lado derecho.

También se suele decir que “la cola derecha es más larga que la izquierda”.

Otra “pista” que teníamos en este sentido es que hasta la media hay 14 valores (2, 3, …15) y desde la media 18 valores (15, 16, …32).

Coeficiente de curtosis (funcion kurtosis del paquete moments)

El coeficiente de curtosis indica si en los extremos del histograma hay una cantidad despreciable de frecuencia acumulada (“colas ligeras”), una cantidad digna de tener en cuenta (“colas pesadas”) o una cantidad “normal”.

En R, el coeficiente de curtosis gira alrededor del valor “3”:

curtosis menor que 3: colas ligeras
curtosis mayor que 3: colas pesadas

library(moments)
kurtosis(muestra_bin)

## [1] 3.009978

En este caso, la curtosis es prácticamente 3 luego las colas no son ni ligeras ni pesadas.

Estadística descriptiva de la distribución normal

Utilizamos nuestra normal de media 5 y desviación típica 0.25

pnorm(valor, mean=5, sd=0.25)
fnorm(valor, mean=5, sd=0.25)
qnorm(valor, mean=5, sd=0.25)
distri_nor
curva_nor
muestra_nor

La media (función mean)

La media es un indicador de la tendencia central, es decir, de cuál es el valor más céntrico de la distribución, del cual emanan los demás valores a modo, podríamos decir, de desviaciones:

mean(muestra_nor)

## [1] 5.005

Se puede demostrar que coincide siempre con el parámetro mean.

La mediana (función median)

La mediana es el valor en el que se iguala o supera por primera vez la frecuencia acumulado de 0.5:

median(muestra_nor)

## [1] 5.005

Se puede demostrar que coincide siempre con la media.

Percentiles (función quantile)

El percentil dd es el valor en el que se iguala o supera por primera vez la frecuencia 0.dd:

quantile(x=muestra_nor, probs=0.05) #percentil 5

##   5% 
## 4.59

quantile(x=muestra_nor, probs=0.95) #percentil 95

##  95% 
## 5.42

quantile(x=muestra_nor, probs=0.50) #percentil 50

##   50% 
## 5.005

El percentil 50 es, por definición, la mediana.

Sumario (función summary y geometría boxplot)

El sumario es un compendio de 6 indicadores:

Mínimo (valor más pequeño con algo de frecuencia)
percentil 25, también llamado “primer cuartil”
percentil 50, también llamado “segundo cuartil” o “mediana”
media
percentil 75, también llamado “tercer cuartil”
Máximo (valor más grande con algo de frecuencia)

summary(muestra_nor)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   4.010   4.840   5.005   5.005   5.170   6.000

También se dibujan aquelos valores que quedan fuera de los bigotes de la caja. Estos valores se consideran atípicos (“outliers”).

ggplot(data = NULL, aes(x="", y = muestra_nor)) +
    geom_boxplot()

Histograma (geometría geom_histogram)

ggplot() + 
    geom_histogram(aes(x = muestra_nor, y = ..density..),
                   color = "white", 
                   fill = "orange", binwidth = 0.01) +
    geom_vline(aes(xintercept = mean(muestra_nor)),
               size=1,
               color="black")

Varianza y desviación típica (funciones var y sd)

var(muestra_nor)

## [1] 0.06244198

Se puede demostrar que coincide siempre con el cuadrado del parámetro sd.

Como la varianza valora el cuadrado de las desviaciones, es más significativa (porque tiene las unidades de medida de la media) su raiz cuadrada, llamada desviación típica:

sqrt(var(muestra_nor))

## [1] 0.2498839

sd(muestra_nor)

## [1] 0.2498839

Se puede demostrar que coincide siempre con el parámetro sd.

Coeficiente de asimetría (función skewness del paquete moments)

library(moments)
skewness(muestra_nor)

## [1] 1.317306e-15

Se puede demostrar que es cero así que la distribución normal es simétrica.

Coeficiente de curtosis (funcion kurtosis del paquete moments)

En R, el coeficiente de curtosis gira alrededor del valor “3”:

curtosis menor que 3: colas ligeras
curtosis mayor que 3: colas pesadas

library(moments)
kurtosis(muestra_nor)

## [1] 2.988044

Se puede demostrar que es exactamente 3 luego las colas no son ni ligeras ni pesadas.

Estadística descriptiva de la distribución de Poisson

Utilizamos nuestra poisson de \(\lambda = 4*2\):

ppois(valor, lambda = 4*2)
fpois(valor, lambda = 4*2)
qpois(valor, lambda = 4*2)
distri_poi
curva_poi
muestra_poi

La media (función mean)

La media es un indicador de la tendencia central, es decir, de cuál es el valor más céntrico de la distribución, del cual emanan los demás valores a modo, podríamos decir, de desviaciones:

mean(muestra_poi)

## [1] 7.99993

Se puede demostrar que coincide siempre con el parámetro \(\lambda\).

La mediana (función median)

La mediana es el valor en el que se iguala o supera por primera vez la frecuencia acumulado de 0.5:

mean(muestra_poi)

## [1] 7.99993

Percentiles (función quantile)

El percentil dd es el valor en el que se iguala o supera por primera vez la frecuencia 0.dd:

quantile(x=muestra_poi, probs=0.05) #percentil 5

## 5% 
##  4

quantile(x=muestra_poi, probs=0.95) #percentil 95

## 95% 
##  13

quantile(x=muestra_poi, probs=0.50) #percentil 50

## 50% 
##   8

El percentil 50 es, por definición, la mediana.

Sumario (función summary y geometría boxplot)

El sumario es un compendio de 6 indicadores:

Mínimo (valor más pequeño con algo de frecuencia)
percentil 25, también llamado “primer cuartil”
percentil 50, también llamado “segundo cuartil” o “mediana”
media
percentil 75, también llamado “tercer cuartil”
Máximo (valor más grande con algo de frecuencia)

summary(muestra_poi)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##       0       6       8       8      10      23

También se dibujan aquelos valores que quedan fuera de los bigotes de la caja. Estos valores se consideran atípicos (“outliers”).

ggplot(data = NULL, aes(x="", y = muestra_poi)) +
    geom_boxplot()

Histograma (geometría geom_histogram)

ggplot() + 
    geom_histogram(aes(x = muestra_poi, y = ..density..),
                   color = "white", 
                   fill = "orange", binwidth = 1) +
    geom_vline(aes(xintercept = mean(muestra_poi)),
               size=1,
               color="black")

Varianza y desviación típica (funciones var y sd)

var(muestra_poi)

## [1] 7.99885

Se puede demostrar que coincide siempre con el parámetro \(\lambda\).

Como la varianza valora el cuadrado de las desviaciones, es más significativa (porque tiene las unidades de medida de la media) su raiz cuadrada, llamada desviación típica:

sqrt(var(muestra_poi))

## [1] 2.828224

sd(muestra_poi)

## [1] 2.828224

Se puede demostrar que coincide siempre con \(\sqrt{\lambda}\).

Coeficiente de asimetría (función skewness del paquete moments)

library(moments)
skewness(muestra_poi)

## [1] 0.3529308

El coeficiente de asimetría es positivo y bastante grande así que hay mucha asimetría “hacia la derecha”. Se puede demostrar que coincide siempre con \(1/\sqrt{\lambda}\).

Esta asimetría se ve en el histograma: a la derecha del eje de la media hay más barras. También se suele decir que “la cola derecha es más larga que la izquierda”.

Otra “pista” que teníamos en este sentido es que hasta la media hay 9 valores (0, 1, …8) y desde la media 16 valores (8, 9, …23).

Coeficiente de curtosis (funcion kurtosis del paquete moments)

En R, el coeficiente de curtosis gira alrededor del valor “3”:

curtosis menor que 3: colas ligeras
curtosis mayor que 3: colas pesadas

library(moments)
kurtosis(muestra_poi)

## [1] 3.121995

En este caso, la curtosis es un poco mayor que 3, luego las colas son algo pesadas.

Se puede demostrar que coincide siempre con \(3 + 1/\lambda\).