SAA Práctica 1: Variable Aleatoria, Probabilidad y teorema central del límite

Author

Jesús Turpín

Published

October 24, 2023

Requisitos previos:

Teoría Temas 1.1 y 1.2
Haber realizado la práctica 0.

library(tidyverse)

1. Experimento aleatorio. Juegos con 2 dados de 6 caras

1.1 Define el espacio muestral cuando realizamos tiradas con dos dados donde obtenemos como resultado de cada experimento aleatorio la suma de los números de cada dado.

S = {2,3,4,5,6,7,8,9,10,11,12}

caras_dado <- 6
combinaciones <- expand.grid(
  dado_1 = 1:caras_dado,
  dado_2 = 1:caras_dado)

combinaciones <- combinaciones %>%
  mutate(suma = dado_1 + dado_2)
combinaciones

   dado_1 dado_2 suma
1       1      1    2
2       2      1    3
3       3      1    4
4       4      1    5
5       5      1    6
6       6      1    7
7       1      2    3
8       2      2    4
9       3      2    5
10      4      2    6
11      5      2    7
12      6      2    8
13      1      3    4
14      2      3    5
15      3      3    6
16      4      3    7
17      5      3    8
18      6      3    9
19      1      4    5
20      2      4    6
21      3      4    7
22      4      4    8
23      5      4    9
24      6      4   10
25      1      5    6
26      2      5    7
27      3      5    8
28      4      5    9
29      5      5   10
30      6      5   11
31      1      6    7
32      2      6    8
33      3      6    9
34      4      6   10
35      5      6   11
36      6      6   12

s <- combinaciones %>%
  select(suma) %>%
  distinct()
s

s <- s$suma
s

 [1]  2  3  4  5  6  7  8  9 10 11 12

1.2 Calcula las probabilidades de cada uno de los elementos del conjunto espacio muestral

probs <- combinaciones %>%
  group_by(suma) %>%
  summarise(n = n()) %>%
  mutate(p = n/nrow(combinaciones))

1.3 Representa gráficamente, mediante un diagrama de columnas las probabilidades anteriores (eje y), donde aparezcan los resultados del espacio muestral en el eje x

probs %>%
  ggplot(aes(x = suma, y = p)) +
  geom_col() +
  scale_x_continuous(breaks = 2:12)

1.4 ¿Los ensayos son dependientes o independientes?

Los ensayos son independientes. Ahora si lanzamos primero un dado y calculamos las probabilidades de obtener suma antes de tirar el otro, las probabilidades de obtener una suma dependen de lo que haya salido en el primer dado.

1.5 Calcula la probabilidad de obtener tres veces consecutivas un 7 en las tres primeras tiradas.

p_7 <- probs %>%
  filter(suma == 7) %>%
  select(p) %>% pull()

p_7**3

[1] 0.00462963

2. Simulación experimento aleatorio con recompensa: dado y ruleta

2.1 Realiza el siguiente experimento utilizando código R y las librerías dplyr (sample_n y operador %>%) y ggplot2: Simula el lanzamiento de un dado. Cada vez que salga un 2, recibes 5 € (el tuyo y 4 de ganancia neta), cada vez que no aciertes, pierdes 1. Calcula la esperanza matemática teórica y mediante un experimento con 10000 simulaciones y una semilla = “777”.

Representa gráficamente en una línea, tu ganancia acumulada (eje y) sobre el número de simulaciones (eje x) y agrega una línea a la gráfica con la esperanza matemática, para medir el grado de suerte. Aplica ahora la semilla 999. ¿Con qué semilla ha habido más suerte?

set.seed(777)
sims <- 10000
dado <- data.frame(n=1:6, ev=rep(-1/6,times=6))
  dado %>%
  sample_n(sims, replace=TRUE) %>%
  mutate(jugador = case_when((n == 2) ~ 4, TRUE ~ -1)) %>%
  mutate(intento = row_number()) %>%
  mutate(ganancia = cumsum(jugador)) %>%
  mutate(ev = -1/6) %>%
  mutate(ev = cumsum(ev)) %>%
  ggplot(aes(intento, ganancia)) +
  geom_line()+
  geom_line(aes(x=intento, y=ev, color="red")) +
    labs(title = "Tirada dado. 10000 simulaciones. Semilla 777") +
    theme_bw()

set.seed(999)
sims <- 10000
dado <- data.frame(n=1:6, ev=rep(-1/6,times=6))
  dado %>%
  sample_n(sims, replace=TRUE) %>%
  mutate(jugador = case_when((n == 2) ~ 4, TRUE ~ -1)) %>%
  mutate(intento = row_number()) %>%
  mutate(ganancia = cumsum(jugador)) %>%
  mutate(ev = -1/6) %>%
  mutate(ev = cumsum(ev)) %>%
  ggplot(aes(intento, ganancia)) +
  geom_line()+
  geom_line(aes(x=intento, y=ev, color="red")) +
    labs(title = "Tirada dado. 10000 simulaciones. Semilla 999") +
    theme_bw()

2.2 Realiza el experimento, apostando al 7 en la ruleta, 1.000.000 veces, con 35 € de ganancia neta cada vez que aciertes y un euro de pérdidas cada vez que falles. Calcula la esperanza matemática teórica y la esperanza matemática según el experimento. Representa visualmente la línea de la esperanza matemática (si en cada tirada ganas Ex). Utiliza como semillas 1 y 123. En base a los resultados obtenidos ¿Qué opinas sobre la varianza y la esperanza matemática de este juego?

set.seed(1)
sims <- 1000000
ruleta <- data.frame(n=0:36, ev=rep(-1/37, times=37))
  ruleta %>%
  sample_n(sims, replace=TRUE) %>%
  mutate(jugador = case_when(n == 7 ~ 35, TRUE ~ -1)) %>%
  mutate(intento = row_number()) %>%
  mutate(ganancia = cumsum(jugador)) %>%
  mutate(ev = -1/37) %>%
  mutate(ev = cumsum(ev)) %>%
  ggplot(aes(intento, ganancia)) +
  geom_line()+
  geom_line(aes(x=intento, y=ev, color="red"))

set.seed(123)
sims <- 1000000
ruleta <- data.frame(n=0:36, ev=rep(-1/37, times=37))
  ruleta %>%
  sample_n(sims, replace=TRUE) %>%
  mutate(jugador = case_when(n == 7 ~ 35, TRUE ~ -1)) %>%
  mutate(intento = row_number()) %>%
  mutate(ganancia = cumsum(jugador)) %>%
  mutate(ev = -1/37) %>%
  mutate(ev = cumsum(ev)) %>%
  ggplot(aes(intento, ganancia)) +
  geom_line()+
  geom_line(aes(x=intento, y=ev, color="red"))

3. El teorema central del límite: lqsa

Utilizando el dataset “lqsa”, simula 1000 veces la extracción de muestras aleatorias de tamaño n = 15 de las 30 cartas. Calcula la media muestral para cada extracción y construye una matriz o data frame con las 5 variables numéricas de la baraja como columnas y en cada fila la media muestral del valor de cada extracción.

Calcula la media y la desviación típica de las mil filas para cada variable.

Representa gráficamente la distribución de los datos para cada variable. Utiliza geom_density() en lugar de barras. Ejemplo del experimento con la variable “locura”.

lqsa <- read.csv(
  "https://raw.githubusercontent.com/jesusturpin/curintel2324/main/data/lqsa.csv"
  )
lqsa_n <- lqsa %>%
  select_if(is.numeric)

simular_experimento <- function(df, n) {
  muestra <- sample_n(df, n)
  medias <- colMeans(muestra)
  return(medias)
}

sims <- 1000
n <- 15
resultados <- replicate(sims, simular_experimento(lqsa_n,n))

df_resultados <- as.data.frame(t(resultados))
df_resultados %>% head()

     Poder Convivencia   Liante Atractivo   Locura
1 4.133333    9.933333 81.20000  12.40000 70.86667
2 6.133333    9.800000 75.73333  11.33333 66.40000
3 4.933333    9.200000 69.06667  12.40000 69.40000
4 3.000000   10.666667 68.26667  13.93333 69.13333
5 3.800000    8.933333 74.93333  12.06667 78.20000
6 4.400000   10.866667 65.86667  16.33333 62.53333

df_resultados %>%
  ggplot(aes(x = Locura)) +
  geom_histogram(bins = 39)

df_resultados %>%
  ggplot(aes(x = Locura)) +
  geom_density()

df_resultados %>%
  ggplot(aes(x = Locura)) +
  geom_histogram(aes(y = after_stat(density)), bins = 39, fill = "blue", color = "black") +
  geom_density(color = "red") +
  scale_y_continuous(sec.axis = sec_axis(~ . * length(df_resultados), name = "Probabilidad"))+
  theme_bw()

df_resultados %>%
  ggplot(aes(x = Locura)) +
  geom_density()

Investiga por tu cuenta el tamaño de la muestra necesario para que se cumpla el teorema.

4. El problema del cumpleaños

Prueba el código del problema del cumpleaños variando los parámetros. Investiga, busca ejemplos y practica el funcionamiento de las funciones sapply, replicate y any

mismo_dia_cumple <- function(n) {
  dias <- sample(1:365, n, replace = TRUE)
  any(duplicated(dias))
}

calcular_prob_sim <- function(n, B=10000) {
  resultados <- replicate(B, mismo_dia_cumple(n))
  mean(resultados)
}

calcular_prob_exac <- function(n){
  prob_unique <- seq(365,365-n+1)/365 
  1 - prod(prob_unique)
}

n <- seq(1,60)
prob_sim <- sapply(n, calcular_prob_sim, B=500)
prob_exac <- sapply(n, calcular_prob_exac)
ggplot(data.frame(n, prob_sim, prob_exac), aes(n, prob_sim))+
  geom_point()+
  geom_line(aes(n, prob_exac), color="red")