Modelos Probabilisticos

Paquetes en R

Código

library(kableExtra)
library(tidyverse)
library(ggplot2)
library(MASS)
library(flextable)
library(vembedr)

[@kableExtra-2; @tidyverse-2; @ggplot2-2; @MASS-2; @flextable-2; @vembedr]

1 Modelos Probabilisticos

Los modelos probabilísticos son herramientas matemáticas que se utilizan para representar y analizar situaciones que involucran incertidumbre. Estos modelos se basan en la teoría de la probabilidad y permiten estimar la probabilidad de ocurrencia de eventos o resultados futuros.

En un modelo probabilístico, se definen variables aleatorias que representan los eventos o resultados que estamos interesados en estudiar. Estas variables aleatorias pueden tomar diferentes valores con ciertas probabilidades asociadas.

Existen dos tipos principales de modelos probabilísticos: los modelos de probabilidad discreta y los modelos de probabilidad continua.

2 Modelos probabilísticos de variable discreta

Se utilizan cuando las variables aleatorias pueden tomar un conjunto finito o numerable de valores posibles. Estos modelos se basan en distribuciones de probabilidad discreta, donde la probabilidad de cada valor se puede asignar de manera precisa. Un ejemplo de esto es el lanzamiento de una moneda, lanzamiento de un dado.

Función masa:

\(f(x)=P(X=x)=\frac{1}{n}\)

donde:

\(n=x_{1},x_{2},x_{3},…,x_{n}\)

Media: \(\mu=E(X)=\frac{n+1}{2}\)

Varianza: \(\sigma ^{2}=Var(X)=\frac{n^{2}-1}{12}\)

Modelos discretos:

2.1 Distribución bernoulli:

Es una distribución de probabilidad discreta que describe el resultado de un experimento o ensayo aleatorio con dos posibles resultados: éxito (generalmente etiquetado como 1) o fracaso (generalmente etiquetado como 0).

Los parámetros clave de la distribución de Bernoulli son:

• Probabilidad de éxito (p): Es la probabilidad de que ocurra un éxito en un solo ensayo. Toma valores entre 0 y 1.

La función de probabilidad de la distribución de Bernoulli está dada por:

\(P(X = k) = p^k * (1 - p)^{1 - k}\)

donde X es la variable aleatoria que representa el resultado del ensayo, k es el valor del resultado (0 o 1) y p es la probabilidad de éxito.

La distribución de Bernoulli es útil para modelar situaciones en las que solo se espera un resultado binario. Algunos ejemplos comunes incluyen el lanzamiento de una moneda (éxito = cara, fracaso = cruz), el éxito o fracaso en una prueba, o la ocurrencia o no de un evento específico.

Parámetros:

Media: \(\mu=E(X)=p\)

Varianza: \(\sigma ^{2}=Var(X)=p(1-p)\)

Ejemplo:

En este problema, la variable \(X\) indica si el dado cae seis \((X = 1)\) o no \((X = 0)\). Sabemos que la probabilidad de que salga un seis es de \(1/6\), por lo tanto \(p = 1/6\). Sustituyendo en la fórmula de la distribución de Bernoulli:

\(P(X=0) = (1-1/6)^{1-0} = (5/6)^1 = 5/6\)

\(P(X=1) = (1/6)^1(1-1/6)^{1-1} = 1/6\)

Por lo tanto, la distribución de la variable \(X\) es:

\(P(X=0) = \frac{5}{6}\)

\(P(X=1) = \frac{1}{6}\)

\(f(x)=P(X=x)=\left\{\begin{matrix}\frac{5}{6}\: ;\: si\: x=0\\ \\\frac{1}{6}\:;\: si\: x=1 \end{matrix}\right.\)

\(f(x)=P(X=x)=\left ( \frac{1}{6} \right )^{x}\left ( \frac{5}{6} \right )^{1-x}\; ;x=0;1\: \rightarrow X\sim Be\left ( \frac{1}{6} \right )\)

Solución en R

Solución en R

Código

# Probabilidad de obtener un seis en un lanzamiento de dado
p <- 1/6

# Generar valores aleatorios de la distribución de Bernoulli
n <- 1000  # número de lanzamientos
x <- rbinom(n, size = 1, prob = p)

# Contar la frecuencia de ocurrencia de cada valor
frecuencia <- table(x)

# Calcular las probabilidades teóricas
probabilidad_0 <- (1 - p)
probabilidad_1 <- p

# Mostrar los resultados
cat("Distribución de Bernoulli para el lanzamiento de un dado:\n")

Distribución de Bernoulli para el lanzamiento de un dado:

Código

cat("X = 0 con probabilidad:", probabilidad_0, "\n")

X = 0 con probabilidad: 0.8333333 

Código

cat("X = 1 con probabilidad:", probabilidad_1, "\n\n")

X = 1 con probabilidad: 0.1666667 

Grafica

Código

library(ggplot2)

# Probabilidad de obtener un seis en un lanzamiento de dado
p <- 1/6

# Generar valores aleatorios de la distribución de Bernoulli
n <- 1000  # número de lanzamientos
x <- rbinom(n, size = 1, prob = p)

# Crear un data frame con los valores observados
data <- data.frame(x = factor(x), count = 1)

# Calcular las probabilidades teóricas
probabilidad_0 <- round((1 - p), 4)
probabilidad_1 <- round(p, 4)

# Graficar la distribución de Bernoulli
ggplot(data, aes(x = x, y = count)) +
  geom_bar(stat = "identity", fill = "steelblue", width = 0.1) +
  geom_text(aes(label = count), vjust = -0.2, color = "black") +
  xlab("Valor de X") +
  ylab("Frecuencia") +
  ggtitle("Distribución de Bernoulli para el lanzamiento de un dado") +
  theme_minimal() +
  annotate("text", x = c("0", "1"), y = c(0, 0), label = c(paste0("Prob. 0: ", probabilidad_0), paste0("Prob. 1: ", probabilidad_1)), vjust = -1.5) +
  scale_x_discrete(labels = c("0", "1"), limits = c("0", "1"))

2.2 Distribución binomial:

Es una generalización del modelo de Bernoulli. Bernoulli es un experimento estadístico que solo puede asumir de dos maneras diferentes, mutuamente excluyentes e independientes conocidos como éxito y fracaso\(X\sim Bin(n,p)\). en este modelo en la mayoría de casos el tamaño de la población (N) es infinito (desconocido), cuando este tamaño de poblacional es finito(conocido) la selección de la muestra se debe hacer con reposición garantizando que la probabilidad de éxito \(p\) permanezca constante en todo el experimento.

Función masa:

\(f(x)=P(X=x)=\binom{n}{x}p^{x}(1-p)^{n-x}\)

donde:

• \(n=1,2,3,…\)

• \(0\leq p\leq 1\)

• \(x=0,1,…,n\)

  Parámetros:

Media: \(\mu=E(X)=np\)

Varianza: \(\sigma ^{2}=Var(X)=np(1-p)\)

Ejemplo:

En una fábrica de galletas, se sabe que el 80% de las galletas producidas son de chocolate. Se selecciona al azar una muestra de 20 galletas de un lote reciente.

a)¿Cuál es la probabilidad de que exactamente 15 galletas sean de chocolate?

\(f(x)=P(X=15)=\binom{20}{15}0.8^{15}(1-0.8)^{20-15}\)

\(P(X=15)=\binom{20}{15} \cdot 0.8^{15}\cdot(0.2)^{5}\)

\(P(X=15)=\frac{20!}{15!-5!}\cdot (0.8^{15})\cdot (0.2^{5})\)

\(P(X=15)\approx 0.1745\)

La probabilidad de que exactamente 15 galletas sean de chocolate es aproximadamente 0.1745.

Resolución en R

Código

n <- 20 # número de ensayos (galletas)
x <- 15 # número de éxitos (galletas de chocolate) deseados
p <- 0.8 # probabilidad de éxito (galleta de chocolate) en un solo ensayo

prob_15_choco <- dbinom(x, size = n, prob = p)

mensaje <- paste0("La probabilidad de que exactamente 15 galletas sean de chocolate es aproximadamente ", round(prob_15_choco, 4), ".")

print(mensaje)

[1] "La probabilidad de que exactamente 15 galletas sean de chocolate es aproximadamente 0.1746."

Grafica

Código

library(ggplot2)

n <- 20  # número de ensayos (galletas)
x <- 0:n  # valores de x para graficar
p <- 0.8  # probabilidad de éxito (galleta de chocolate) en un solo ensayo

# Calcular las probabilidades para cada valor de x
probs <- dbinom(x, size = n, prob = p)

# Crear un data frame con los valores de x y las probabilidades correspondientes
data <- data.frame(x = x, prob = probs)

# Graficar la distribución binomial
ggplot(data, aes(x = as.factor(x), y = prob)) +
  geom_bar(stat = "identity", fill = ifelse(data$x == 15, "red", "steelblue")) +
  geom_text(aes(label = ifelse(data$x == 15, round(prob, 4), "")), vjust = -0.5, color = "black") +
  xlab("Número de galletas de chocolate") +
  ylab("Probabilidad") +
  ggtitle("Distribución Binomial") +
  theme_minimal()

b)¿Cuál es la probabilidad de que al menos 18 galletas sean de chocolate?

La fórmula para calcular esta probabilidad es:

\(P(al menos 18)=\sum_{x=18}^{20}\binom{20}{x}\cdot (0.8^{x})\cdot (0.2^{20-x})\)

Donde:

\(P(X=18)=\binom{20}{18}\cdot (0.8^{18})\cdot (0.2^{2})\)

\(P(X=19)=\binom{20}{19}\cdot (0.8^{19})\cdot (0.2^{1})\)

\(P(X=20)=\binom{20}{20}\cdot (0.8^{20})\cdot (0.2^{0})\)

Aplicando las formulas obtendremos:

\(P(X=18)\approx 0.1369\)

\(P(X=19)\approx 0.0576\)

\(P(X=20)\approx 0.0115\)

Por lo tanto, la probabilidad de que al menos 18 galletas sean de chocolate es la suma de estas tres probabilidades:

\(P(al menos 18)=P(18)+P(19)+P(20)\approx 0.206\)

El resultado de este cálculo es que la probabilidad de que al menos 18 galletas sean de chocolate es aproximadamente 20.6%

Resolución por R

Código

n <- 20  # número de ensayos (galletas)
p <- 0.8  # probabilidad de éxito (galleta de chocolate) en un solo ensayo

# Calcular las probabilidades individuales
prob_18 <- dbinom(18, size = n, prob = p)
prob_19 <- dbinom(19, size = n, prob = p)
prob_20 <- dbinom(20, size = n, prob = p)

# Calcular la probabilidad de al menos 18 galletas de chocolate
prob_at_least_18 <- prob_18 + prob_19 + prob_20

mensaje <- paste0("La probabilidad de que al menos 18 galletas sean de chocolate es aproximadamente ", round(prob_at_least_18, 4), ".")

print(mensaje)

[1] "La probabilidad de que al menos 18 galletas sean de chocolate es aproximadamente 0.2061."

Grafica

Grafica opción 1

Código

library(ggplot2)

n <- 20  # número de ensayos (galletas)
p <- 0.8  # probabilidad de éxito (galleta de chocolate) en un solo ensayo

# Calcular las probabilidades individuales
prob_18 <- dbinom(18, size = n, prob = p)
prob_19 <- dbinom(19, size = n, prob = p)
prob_20 <- dbinom(20, size = n, prob = p)

# Calcular la probabilidad de al menos 18 galletas de chocolate
prob_at_least_18 <- prob_18 + prob_19 + prob_20

# Crear un data frame con los valores y probabilidades correspondientes
data <- data.frame(x = c("18", "19", "20", "Mínimo 18"),
                   prob = c(prob_18, prob_19, prob_20, prob_at_least_18))

# Graficar la distribución binomial
ggplot(data, aes(x = x, y = prob)) +
  geom_bar(stat = "identity", fill = "steelblue") +
  geom_text(aes(label = round(prob, 4)), vjust = -0.5, color = "black") +
  xlab("Número de galletas de chocolate") +
  ylab("Probabilidad") +
  ggtitle("Distribución Binomial") +
  ylim(0, 0.25) +
  theme_minimal()

Grafica opción 2

Código

library(ggplot2)

n <- 20  # número de ensayos (galletas)
x <- 0:n  # valores de x para graficar
p <- 0.8  # probabilidad de éxito (galleta de chocolate) en un solo ensayo

# Calcular las probabilidades para cada valor de x
probs <- dbinom(x, size = n, prob = p)

# Crear un data frame con los valores de x y las probabilidades correspondientes
data <- data.frame(x = x, prob = probs)

# Definir el color de relleno
data$fill <- ifelse(data$x %in% c(18, 19, 20), "red", "steelblue")

# Graficar la distribución binomial
ggplot(data, aes(x = as.factor(x), y = prob, fill = fill)) +
  geom_bar(stat = "identity") +
  geom_text(data = subset(data, x %in% c(18, 19, 20)), aes(label = round(prob, 4)), vjust = -0.5, color = "black") +
  xlab("Número de galletas de chocolate") +
  ylab("Probabilidad") +
  ggtitle("Distribución Binomial") +
  scale_fill_manual(values = c("red", "steelblue")) +
  theme_minimal() +
  guides(fill = FALSE)

2.3 Distribución Poisson:

Diseñada por el matemático francés Possion en 1837. Para denotar que una variable aleatoria X sigue una distribución de Poisson, se usa \(𝑋 ~ 𝑃𝑜𝑖(𝜆𝑖 ), 𝑖 = 1, ⋯ , 𝑁\).La distribución de Poisson se usa para describir el número de ocurrencias de eventos en un espacio de observación limitado. Por ejemplo, puede describir la cantidad de defectos en el sistema mecánico de una aeronave o la cantidad de llamadas a un centro de llamadas en una hora. La distribución de Poisson se usa comúnmente para control de calidad, investigación de confiabilidad / vida útil y seguros. La ocurrencia de eventos de un cierto tipo en el tiempo o un espacio sigue un Proceso de Poisson con tasa de ocurrencia \(λ>0\) por unidad de tiempo o espacio,

La distribución Poisson se utiliza para caracterizar probabilísticamente el número de veces que ocurre un evento en relación con una unidad de medida bien definida (como una unidad de tiempo o espacio, por ejemplo), de forma que:

1. La probabilidad de que el evento ocurra en una unidad de medida dada es igual para todas las unidades.

2. El número de eventos que ocurren en una unidad de medida es independiente del número de eventos que ocurren en otras unidades.

   \(P(X=k) = \frac{{e^{-\lambda} \lambda^k}}{{k!}}\)

\(P(X=k)\) es la probabilidad de que la variable aleatoria X tome el valor k, \(( e)\) es la base del logaritmo natural (aproximadamente 2.71828), \(( \lambda )\)es el parámetro de la distribución, que representa el número esperado de ocurrencias del evento en un intervalo dado, y \(( k )\) es el número de ocurrencias del evento.

Parámetros:

• Esperanza: \(\textsf{E}[k] =\lambda\)

• Varianza: \(\textsf{V}[k] =\lambda\)

Ejemplo:

Cada año ocurre un promedio de 24 accidentes aéreos. Calcule la probabilidad de que ocurra un accidente en un mes.

Solución:

X=n° de accidentes aéreos que ocurren en 1 mes

\(X\sim Pois (2 )\)

x=1

\(\mu =2\)

\[ f(x)=P(X=x)=\frac{e^{-\mu }\cdot \mu ^{x}}{x!} \]

\[ f(x)=P(X=1)=\frac{e^{-2 }\cdot 2 ^{1}}{1!}=e^{-2}\cdot 2\]

\[ f(x)=P(X=1)=0.2707=27.07% \]

La probabilidad de que ocurra un accidente en un mes es de 27.07%

Para realizarlo en R studio:

Solución utilizando R studio

Código

dpois(1,2)

[1] 0.2706706

Código

# Promedio anual de accidentes
promedio_anual <- 24

# Cálculo del promedio mensual
promedio_mensual <- promedio_anual / 12

# Número de accidentes en un mes (x)
x <- 1

# Cálculo de la probabilidad utilizando la distribución de Poisson
probabilidad <- dpois(x, lambda = promedio_mensual)

# Imprimir la probabilidad
probabilidad

[1] 0.2706706

Grafica en ggplot2

Código

library(ggplot2)
# Promedio anual de accidentes
promedio_anual <- 24

# Cálculo del promedio mensual
promedio_mensual <- promedio_anual / 12

# Generación de valores para x
x <- 0:6

# Cálculo de las probabilidades utilizando la distribución de Poisson
probabilidades <- dpois(x, lambda = promedio_mensual)

# Creación del data frame
data <- data.frame(x = x, probabilidad = probabilidades)

# Gráfico de la distribución de Poisson
ggplot(data, aes(x = as.factor(x), y = probabilidad, fill = (x == 1))) +
  geom_col(color = "black", width = 0.1) +
  scale_fill_manual(values = c("FALSE" = "blue", "TRUE" = "red")) +
  geom_text(aes(label = round(probabilidad, 4)), vjust = -0.5, color = "black", size = 4) +
  scale_y_continuous(limits = c(0, 0.4), expand = c(0, 0)) +
  labs(title = "Diagrama de Bastones", x = "Número de accidentes en un mes", y = "Probabilidad") +
  theme(legend.position = "none")

b . Calcule la probabilidad de que ocurran máximo 5 accidentes en un año

\[ X\sim Pois (4) \]

\[ x\leq 5 \]

\[ \mu =4 \]

\[ f(x)=P(X=x)=\frac{e^{-\mu }\cdot \mu ^{x}}{x!} \]

\[ f(x)P(X\leq 5) =P(X=0)+P(X=1)+P(X=2)+P(X=3)+P(X=4)+P(X=5) \]

\(P(X\leq 5) =\frac{e^{-4 }\cdot 4 ^{0}}{0!}+\frac{e^{-4 }\cdot 4 ^{1}}{1!}+\frac{e^{-4 }\cdot 4 ^{2}}{2!}+\frac{e^{-4 }\cdot 4 ^{3}}{3!}+\frac{e^{-4 }\cdot 4 ^{4}}{4!}+\frac{e^{-4 }\cdot 4 ^{5}}{5!}\)

\(P(X\leq5)=0.0183+0.0733+0.1465+0.1954+0.1954+0.1563=0.7852=78.52%\)

La probabilidad de que ocurra un accidente en un mes es de 78.52%

Para resolver en R base

Resolver por R

Código

ppois(5,4)

[1] 0.7851304

Grafica por ggplot2

Código

library(ggplot2)

pois_sum <- function(lambda, lb, ub, col = 4, lwd = 1, ...) {
  x <- 0:(lambda + lambda * 2)
  
  if (missing(lb)) {
    lb <- min(x)
  }
  if (missing(ub)) {
    ub <- max(x)
  }
  
  data <- data.frame(x = x, y = dpois(x, lambda = lambda))
  data$y_rounded <- round(data$y, 4)
  
  data$color <- ifelse(x >= 0 & x <= 5, "red", "blue")
  
  ggplot(data, aes(x = x, y = y)) +
    geom_segment(aes(xend = x, yend = 0, color = color), size = lwd) +
    geom_point(aes(y = y, color = color), size = 3) +
    geom_text(aes(label = y_rounded), vjust = -0.5, size = 3) +
    theme_minimal() +
    labs(x = "accidentes aereos", y = "P(X = x)") +
    ylim(0, max(data$y) * 1.1) +
    scale_color_manual(values = c("blue", "red")) +
    scale_x_continuous(breaks = seq(0, 12), limits = c(0, 12)) +
    theme(legend.position = "none")
}

pois_sum(lambda = 4, ub = 5, lwd = 2,
         ylab = "P(X = x)", xlab = "accidentes aereos")

2.4 Distribución Hipergeométrica:

La distribución hipergeométrica surge a partir del “número de éxitos en \(n\) ensayos dependientes de un experimento de Bernoulli”. Un experimento hipergeométrico con parámetros n, M, y N está basado en las siguientes condiciones:

1. Se elige una muestra sin reemplazo de n elementos de un conjunto compuesto por N elementos, de los cuales M tienen una característica de interés..

2. Cada elemento se puede caracterizar como un “éxito (el elemento tiene la característica de interés) o como un”fracaso” (el elemento no tiene la característica de interés).

3. La probabilidad de exito puede variar según cada nueva observación

Se dice que una variable aleatoria X tiene una distribución hipergeométrica con parámetros n, M, y N, lo que se escribe: \(X\sim Hg(n,M,N)\). La función masa o cuantia es:

\(f(x)=P(X=x)=\frac{\binom{M}{x}\binom{N-M}{n-x}}{\binom{N}{n}}\)

\(n\) :muestra

\(N\):Población

\(M\):Numero de sucesos que cumple con la caracteristica objeto de estudio

\(X\):Número de exitos obtenidos en la muestra

\(x\) \([0,1,2,3…n]\)

Parámetros:

Media: \(\mu=E(X)=n\frac{M}{N}\)

Varianza: \(\sigma ^{2}=Var(X)=n\frac{M}{N}\left ( 1-\frac{M}{N} \right )\left ( \frac{N-n}{N-1} \right )\)

ejemplo:

De un lote de 10 misiles, se seleccionan 4 al azar y se disparan. Si el lote contiene 3 misiles defectuosos que no pueden dispararse, ¿cuál es la probabilidad de que a) los 4 puedan dispararse? b) a lo sumo fallen 2?

Solución:

a) N=10, M=7, n=4, x=4

\(P(X=4)=\frac{\binom{7}{4}\binom{10-3}{4-4}}{\binom{10}{4}}\)

\(P(X=4)=\frac{\binom{7}{4}\binom{7}{0}}{\binom{10}{4}}\)

\(P(X=4)=\frac{35}{210}\)

\(P(X=4)=\frac{1}{6}\)

La probabilidad que los 4 miseles puedan disparase es de 16.6%

b) N=10, M=3, n=4 X=0,1,2

\(P(X\leq 2)= P(X=0)+P(X=1)+P(X=2)\)

\(P(X\leq 2)=\frac{1}{6}+\frac{1}{2}+\frac{3}{10}\)

\(P(X\leq 2)=\frac{24}{25}\)

La probabilidad de que a lo sumo fallen 2 misiles es de 96.6%

Para realizarlo en R studio:

a)N=10, M=7, n=4, x=4

Para calcular la probabilidad de que los 4 misiles puedan dispararse, necesitamos calcular la probabilidad de seleccionar los 4 misiles funcionales de los 7 misiles restantes (descontando los 3 misiles defectuosos).La probabilidad se calcula como:Calcular la probabilidad de que los 4 misiles puedan dispararse

Código

# parametros
n <- 4
M <- 7
N <- 10
# P(X = 4)
# la parametrizacion de esta rutina es diferente a la presentada en la formula
sum(dhyper(x = 4, m = M, n = N-M, k = n))

[1] 0.1666667

Código

# Parámetros
n <- 4
M <- 7
N <- 10
x <- 4

# Calcular las probabilidades de la distribución
prob <- dhyper(x = 0:n, m = M, n = N - M, k = n)

# Redondear las probabilidades a cuatro decimales
prob_rounded <- round(prob, 4)

# Función de masa de probabilidad
plot(x = 0:n, y = prob_rounded, xlab = "Número de éxitos", ylab = "P(x)", pch = 15, col = "blue", ylim = c(0, 0.6), main = "Distribución Hipergeométrica")
segments(x0 = 0:n, y0 = 0, x1 = 0:n, y1 = prob_rounded, lwd = 2, col = "blue")

# Línea en x=4 de color rojo
segments(x0 = x, y0 = 0, x1 = x, y1 = prob_rounded[x+1], lwd = 2, col = "red")

# Punto en x=4 de color rojo
points(x, prob_rounded[x+1], pch = 15, col = "red")

# Etiquetas con las probabilidades redondeadas
text(x = 0:n, y = prob_rounded, labels = prob_rounded, pos = 3)

b)N=10, M=3, n=4 X=0,1,2

Para calcular la probabilidad de que a lo sumo fallen 2 misiles, necesitamos sumar las probabilidades de los siguientes casos:

Ningún misil falla, Exactamente 1 misil falla, Exactamente 2 misiles fallan.

Código

# parametros
n <- 4
M <- 3
N <- 10
# P(X <= 2)
# la parametrizacion de esta rutina es diferente a la presentada en la formula
sum(dhyper(x = 0:2, m = M, n = N-M, k = n))

[1] 0.9666667

Grafica de bastones

Código

# Parámetros
n <- 4
M <- 3
N <- 10
x <- 0:2
# Calcular las probabilidades de la distribución
prob <- dhyper(x = 0:n, m = M, n = N - M, k = n)
# f.m.p.
fx <- dhyper(x = x, m = M, n = N-M, k = n)
# f.d.a.
Fx <- phyper(q = x, m = M, n = N-M, k = n)
fxx<-c(1/6,1/2,3/10)
# f.m.p
plot(x = x, y = fxx, xlab = "Numero de misiles", ylab = "P(x)", pch = 15, col = "blue", xaxp = c(0, 2, 2),  yaxp = c(-70, 70, 5), main = "A lo sumo fallen 2 misiles (P<=2)")
axis(2, at = c(1/6, 1/2,3/10), labels = c("1/6", "1/2","3/10"), las = 1)
segments(x0 = x, y0 = 0, x1 = x, y1 = fx, lwd = 2, col = "blue")

Grafica de densidad

Código

# Parámetros
n <- 4
M <- 3
N <- 10
x <- 0:2
# Calcular las probabilidades de la distribución
prob <- dhyper(x = 0:n, m = M, n = N - M, k = n)
# f.m.p.
fx <- dhyper(x = x, m = M, n = N-M, k = n)
# f.d.a.
Fx <- phyper(q = x, m = M, n = N-M, k = n)
# f.d.a.
plot(x = c(0, x), y = c(0, Fx), type = "s", xlab = "Numero de misiles", ylab = "F(x)", col = "blue", lwd = 2, xaxp = c(0, 2, 2),main = "A lo sumo fallen 2 misiles (P<=2)")
points(x, Fx, col = "blue", pch = 15)

5.34 ¿Cuál es la probabilidad de que una camarera se rehúse a servir bebidas alcohólicas a sólo 2 menores si verifi ca al azar 5 identifi caciones de 9 estudiantes, de los cuales 4 son menores de edad?

Solución:

\(N=9, n=5,k=4,X=2\)

\(P(X=2)=\frac{\binom{4}{2}\binom{9-4}{5-2}}{\binom{9}{5}}\)

\(P(X=2)=\frac{\binom{4}{2}\binom{5}{3}}{\binom{9}{5}}\)

\(P(X=2)=\frac{6*10}{126}\)

\(P(X=2)=\frac{10}{21}\)

La probabilidad de que una camarera se rehúse a servir bebidas alcohólicas a sólo 2 menores es de 47.62%

Para realizarlo en R studio:

Código

# Parámetros
N <- 9
M <- 4
n <- 5
x <- 2

# Calcular la probabilidad utilizando la distribución hipergeométrica
prob <- dhyper(x = x, m = M, n = N - M, k = n)

# Imprimir el resultado
prob

[1] 0.4761905

Grafica en Ggplot2

Código

library(ggplot2)

# Parámetros
N <- 9
M <- 4
n <- 5
x <- 0:n

# Calcular las probabilidades de la distribución hipergeométrica
prob <- dhyper(x = x, m = M, n = N - M, k = n)

# Redondear las probabilidades a 4 decimales
prob_rounded <- round(prob, 4)

# Crear un data frame con los valores de x y las probabilidades redondeadas
df <- data.frame(x = x, prob = prob_rounded)

# Crear una variable para resaltar x = 2
df$resaltar <- ifelse(df$x == 2, TRUE, FALSE)

# Graficar la probabilidad
ggplot(df, aes(x = x, y = prob)) +
  geom_segment(aes(x = x, y = 0, xend = x, yend = prob, color = ifelse(df$x == 2, "red", "blue")), size = 1) +
  geom_text(aes(label = prob_rounded), vjust = -0.5, color = "blue", size = 4) +
  geom_point(aes(x = x, y = prob, color = ifelse(df$x == 2, "red", "blue")), size = 3) +
  scale_color_manual(values = c("blue", "red")) +
  scale_y_continuous(limits = c(0, max(prob) * 1.2), expand = c(0, 0.1)) +
  labs(title = "Probabilidad de que la camarera se rehúse a servir a 2 menores", x = "Número de menores de edad", y = "Probabilidad") +
  theme(legend.position = "none")

Grafica en R base

Código

# Parámetros
n <- 5
M <- 4
N <- 9
x <- 0:n

# Calcular las probabilidades de la distribución
prob <- dhyper(x = 0:n, m = M, n = N - M, k = n)
# f.m.p.
fx <- dhyper(x = x, m = M, n = N - M, k = n)

# Redondear las probabilidades a 4 decimales
prob <- round(prob, 4)

# f.m.p
plot(x = x, y = fx, xlab = "Número de menores de edad", ylab = "Probabilidad", pch = 15, col = "blue", main = "Probabilidad que la camarera se rehúse servir a 2 menores", ylim = c(0, 0.7))
segments(x0 = x, y0 = 0, x1 = x, y1 = fx, lwd = 2, col = "blue")

# Línea en x=2 de color rojo
segments(x0 = 2, y0 = 0, x1 = 2, y1 = fx[3], lwd = 2, col = "red")

# Punto en x=2 de color rojo
points(2, fx[3], pch = 15, col = "red")

# Agregar etiquetas de probabilidad
text(x = x, y = fx, labels = prob, pos = 3)

2.5 Distribución Geométrica:

Para denotar que una variable aleatoria 𝑋 sigue una distribución geométrica, se usa \(𝑋 ~ 𝐺𝑒𝑜𝑚(𝑝), 𝑖 = 1, ⋯ , 𝑁\). La fórmula para hallar la distribución geométrica es:

\(P(X = x)=p(1-p)^{x-1}\)

Donde: 𝑥 = 0,1,2, …, es el número de fallas en una secuencia antes de que ocurra el primer éxito. 𝑝 es la probabilidad de éxito en cada prueba.

X es el número de ensayos hasta obtener el primer éxito. Pruebas independientes.

Parámetros:

• Media: \(\mu=\frac{1}{p}\)

• Varianza: \(V(X)=\frac{(1-p)}{p^2}\)

Ejemplo:

Suponga que cada una de sus llamadas a una estación de radio popular tienen una probabilidad de 0.02 de ser respondida. Asumiendo que las llamadas son independientes

a) ¿Cuál es la probabilidad que le respondan a la decima llamada?

X= Número de llamadas a la estación hasta ser atendido

x=10

Éxito: llamada respondida

Fracaso: llamada no respondida

\(p=0.02 \rightarrow \: éxito\)

\(1- p=0.98 \rightarrow \: fracaso\)

\(f(x)=p(1-p)^{x-1}\)

\(f(10)=0.02(0.98)^{10-1} =0.0167=1.67%\)

La probabilidad que respondan a la décima llamada es del1.67%

b) ¿Cuál es el número medio de llamadas para conectar?

\(E\left [ X \right ]=\frac{1}{p}\)

\(E\left [ X \right ]=\frac{1}{p}=\frac{1}{0.02}=50\)

El numero medio de llamadas para conectar es 50

Para resolver por R:

a)

Grafica por R base

Código

# Cálculo de la probabilidad de la distribución geométrica
# Definir el parámetro p
p <- 0.02
# Definir el valor específico x
x <- 10
# Calcular la probabilidad de obtener el valor x
probabilidad <- p * (1 - p)^(x - 1)
# Redondear la probabilidad a 4 decimales
probabilidad_redondeada <- round(probabilidad, 4)
# Imprimir el resultado
cat("La probabilidad de obtener el valor", x, "en una distribución geométrica con p =", p, "es:", probabilidad_redondeada)

La probabilidad de obtener el valor 10 en una distribución geométrica con p = 0.02 es: 0.0167

Código

# Gráfico de la distribución geométrica
x_values <- 1:20  # Valores de x para graficar
probabilidades <- p * (1 - p)^(x_values - 1)  # Probabilidades correspondientes a cada valor de x

plot(x_values, probabilidades, type = "h", lwd = 2, col = "blue",
     xlab = "Número de llamadas", ylab = "Probabilidad",
     main = "Distribución Geométrica")

# Puntos en todas las líneas del gráfico
points(x_values, probabilidades, pch = 19, col = "blue")

# Punto en el valor x de color rojo
points(x, probabilidad, pch = 19, col = "red")
# Etiqueta de la probabilidad redondeada en el punto x
text(x, probabilidad, labels = probabilidad_redondeada, pos = 3)

# Línea en el valor x de color rojo
segments(x, 0, x, probabilidad, lwd = 2, col = "red")

Grafica por ggplot2

Código

# Cargar la librería necesaria
library(ggplot2)

# Definir el parámetro p
p <- 0.02

# Crear un vector de valores de x
x <- 0:60

# Calcular la probabilidad para cada valor de x
probabilidad <- p * (1 - p)^(x)

# Crear un data frame con los valores de x y sus probabilidades
data <- data.frame(x, probabilidad)

# Graficar la distribución geométrica
ggplot(data, aes(x = x, y = probabilidad)) +
  geom_bar(stat = "identity", fill = ifelse(x == 10, "red", "steelblue")) +
  geom_text(data = data[data$x == 10, ], aes(x = x, y = probabilidad, label = round(probabilidad, 4)), vjust = -0.5, color = "red") +
  xlab("Valor") +
  ylab("Probabilidad") +
  ggtitle("Distribución Geométrica (p = 0.02)")

b)

Código

# Cálculo de la esperanza matemática de una distribución geométrica
# Definir el parámetro p
p <- 0.02
# Calcular la esperanza matemática
esperanza <- 1 / p
# Imprimir el resultado
cat("La esperanza matemática de la distribución geométrica con p =", p, "es:", esperanza)

La esperanza matemática de la distribución geométrica con p = 0.02 es: 50

2.6 Distribución Binomial negativa:

La distribución binomial negativa es una distribución discreta, que simula el número de ensayos necesarios para producir un número específico de eventos. Cada prueba tiene dos posibles resultados. La distribución binomial negativa también puede modelar el número de no eventos que deben ocurrir para observar un número específico de resultados. Para representar que una variable aleatoria 𝑋 sigue una distribución binomial negativa, se usa \(X\sim BN(k,p)\)

Donde 𝑋 es el número de ensayos independientes hasta obtener el 𝒌-ésimo éxito y p, la probabilidad de éxito en los ensayos.

\(P(X = x) = {x-1 \choose k-1}p^k (1-p)^{x-k}\)

\(0 < p <1\)

Parámetros:

Media:

\(\mu=\frac{k(1-p)}{p}\)si se piensa en el número de fracasos.

\(\mu=\frac{k}{p}\) si se cuenta también los \(k−1\) éxitos.

Varianza: \(V(X)=\frac{k(1-p)}{p^2}\)

ejemplo:

5.49 La probabilidad de que una persona que vive en cierta ciudad tenga un perro es de 0.3. Calcule la probabilidad de que la décima persona entrevistada al azar en esa ciudad sea la quinta que tiene un perro.

Solución:

\(k\)=5, \(X\)=10, \(P\)=0.3

\(P(X = 5) = {10-1 \choose 5-1}0.3^5(1-0.3)^{5}\)

\(P(X = 5) = (126)0.3^5(0.7)^{5}\) \(P(X = 5) =0.0515\)

La probabilidad que la decima persona entrevistada sea la quinta que tiene un perro es de 5.15%

Para realizarlo en R studio:

Resolución por R

Código

x<-10
k<-5
p<-0.3
pbn<-dnbinom(x-k,size = k, prob = p)
paste0("La probabilidad que la decima persona entrevistada sea la quinta es ",format(100*pbn,digits = 3),"%.")

[1] "La probabilidad que la decima persona entrevistada sea la quinta es 5.15%."

Gráfica

Código

library(ggplot2)

p <- 0.3
k <- 5

x <- 0:20  # Valores de x para graficar
pbn <- dnbinom(x - k, size = k, prob = p)  # Probabilidades correspondientes a cada valor de x

data <- data.frame(x, pbn)

ggplot(data, aes(x = as.factor(x), y = pbn, fill = as.factor(x == 10))) +
  geom_col(color = "black", width = 0.5) +
  scale_fill_manual(values = c("FALSE" = "steelblue", "TRUE" = "red")) +
  geom_text(aes(label = ifelse(x == 10, round(pbn, 4), "")),
            vjust = -0.5, color = "black", size = 3) +
  xlab("Número de personas que tienen perro") +
  ylab("Probabilidad") +
  ggtitle("Distribución Binomial Negativa") +
  theme_minimal() +
  theme(legend.position = "none")

5.52 Un científico inocula a varios ratones, uno a la vez, el virus que produce una enfermedad, hasta que encuentra a 2 que contraen la enfermedad. Si la probabilidad de contraer la enfermedad es de 1/6, ¿cuál es la probabilidad de que tenga que inocular a 8 ratones?

\(k=8, X=2, P=0.16\)

\(P(X = 8) = {8-1 \choose 2-1}\frac{1}{6}^2(1-\frac{1}{6})^{6}\)

\(P(X = 8) = {7 \choose 1}\frac{1}{6}^2(\frac{1}{6})^{6}\)

\(P(X = 8) =0.0651\)

La probabilidad que tenga que inocular a 8 ratones es 6.51%

Para realizarlo en R studio

Resolución en R

Código

x<-8
k<-2
p<-1/6
pbn<-dnbinom(x-k,size = k, prob = p)
paste0("La probabilidad que tenga que inocular a 8 ratones es ",format(100*pbn,digits = 3),"%.")

[1] "La probabilidad que tenga que inocular a 8 ratones es 6.51%."

Grafica

Código

library(ggplot2)

x <- 0:15
k <- 2
p <- 1/6
pbn <- dnbinom(x - k, size = k, prob = p)

# Crear un data frame con los valores
data <- data.frame(x = x, pbn = pbn)

# Filtrar el data frame para mantener solo x=8
data_filtered <- data[data$x == 8, ]

# Graficar la distribución binomial negativa
ggplot(data, aes(x = as.factor(x), y = pbn, fill = as.factor(x == 8))) +
  geom_col(color = "black", width = 0.5) +
  scale_fill_manual(values = c("FALSE" = "steelblue", "TRUE" = "red")) +
  geom_text(data = data_filtered,
            aes(label = paste0("Probabilidad: ", format(100 * pbn, digits = 3), "%")),
            vjust = -0.5, color = "black", size = 3) +
  xlab("Número de inyecciones requeridas") +
  ylab("Probabilidad") +
  ggtitle("Distribución Binomial Negativa") +
  theme_minimal() +
  theme(legend.position = "none")

3 Modelos probabilísticos de variable continua

3.1 Distribución normal:

La distribución normal se caracteriza por tener una forma de campana simétrica alrededor de su media. Su gráfica es una curva suave y continua, que se extiende desde menos infinito hasta más infinito. La curva alcanza su punto máximo en la media y tiene una desviación estándar que determina qué tan “ancha” o “angosta” es la distribución.

La función de densidad de probabilidad de la distribución normal es conocida como la curva de campana o curva gaussiana, y está determinada por dos parámetros: la media (μ) y la desviación estándar (σ). La media representa el centro de la distribución, mientras que la desviación estándar indica la dispersión de los datos alrededor de la media.

\[f(x)\;= \; \frac{1}{\sqrt{2 \pi \sigma^2}} \, e^{-\frac{(x-\mu)^2} {2 \sigma^2}}, \qquad \text{para todo $x$ real}\]

Ejemplo: La tasa de remuneración media por hora para administrativos financieros en una determinada región es $32.62 y la desviación estándar es $2.32 (Bureau of Labor Statistics, septiembre de 2005). Suponga que estas tasas de remuneración están distribuidas normalmente.

a. ¿Cuál es la probabilidad de que un directivo financiero tenga una remuneración entre $30 y $35 por hora?

Solución:

Para resolver este problema, necesitamos estandarizar los valores de la remuneración utilizando la media y la desviación estándar. Definamos:

- \(\mu\): media de la remuneración \(\mu = 32.62\)

- \(\sigma\): desviación estándar de la remuneración \(\sigma = 2.32\)

- \(X\): variable aleatoria que indica la remuneración de un director financiero

Luego, podemos estandarizar los valores de $X$ utilizando la fórmula de la distribución normal estándar:

\(Z = \frac{X - \mu}{\sigma}\)

De esta forma, podemos calcular la probabilidad de que la remuneración de un director financiero esté entre $30$ y $35$ dólares por hora como la probabilidad de que $Z$ esté entre ciertos límites:

\(P(30 < X < 35) = P\left(\frac{30 - \mu}{\sigma} < Z < \frac{35 - \mu}{\sigma}\right)\)

Sustituyendo los valores correspondientes, obtenemos:

\(P(30 < X < 35) = P\left(\frac{30 - 32.62}{2.32} < Z < \frac{35 - 32.62}{2.32}\right) = P(-1.13 < Z < 1.01)\)

Ahora, necesitamos utilizar una tabla de la distribución normal estándar (o un software estadístico) para obtener la probabilidad correspondiente. La probabilidad de que \(Z\) esté entre \(-1.13\) y \(1.01\) es aproximadamente \(0.7734\). Por lo tanto, podemos concluir que la probabilidad de que un director financiero tenga una remuneración entre \(30\) y \(35\) dólares por hora es del \(77.34\%\)

Solución en R studio:

Solución

Código

pnorm(35,32.62,2.32)-pnorm(30,32.62,2.32)

[1] 0.7181382

Grafica

Código

library(mosaic)

# Parámetros de la distribución normal
media <- 32.62
desviacion <- 2.32

# Cálculo de la probabilidad acumulada
probabilidad <- pnorm(35, media, desviacion) - pnorm(30, media, desviacion)

# Rango de valores para graficar
x <- seq(media - 4 * desviacion, media + 4 * desviacion, 0.01)

# Valores de la densidad de la distribución normal
y <- dnorm(x, media, desviacion)

# Crear un data frame con los valores de x y y
data <- data.frame(x, y)

# Graficar la distribución normal utilizando curve
curve(dnorm(x, media, desviacion), from = media - 4 * desviacion, to = media + 4 * desviacion,
      xlab = "Valores", ylab = "Densidad", main = "Distribución Normal",
      col = "steelblue", lwd = 2)

# Resaltar el área bajo la curva correspondiente a la probabilidad acumulada
polygon(c(30, seq(30, 35, 0.01), 35), c(0, dnorm(seq(30, 35, 0.01), media, desviacion), 0),
        col = "steelblue", border = NA)

# Agregar texto para mostrar la probabilidad acumulada
text(32.5, 0.08, paste("Probabilidad:", round(probabilidad, 4)), adj = c(0, 1), col = "red")