EJERCICIOS PROPUESTOS
- Cuando se produce una sobrecarga eléctrica en una cadena productiva en número de componentes averiadas está descrito mediante una variable aleatoria, \(X\), que puede tomar los valores 1,2,3,4 y 5, con probabilidades \[
p_1=0.11, \space p_2=0.22, \space p_3=0.45, \space p_4=0.13, \space p_5=0.09.
\] Para simular dicha variable se usa el algoritmo propuesto. ¿Cuál es el número medio de operaciones para simular un valor de \(X\)? Construir un algoritmo basado en un árbol binario óptimo. ¿Cuál es ahora el número medio de comparaciones? Dar un algoritmo basado en el método de la tabla guía. Calcular, también en este caso, el número medio de comparaciones.
#Algoritmo propuesto
ncomp <<- 0
sim1N <- function(n=1000){
X <- numeric(n)
for(isim in 1:n)
{
U <- runif(1,0,1)
if(U<=0.11){X[isim]<-1
ncomp<<-ncomp+1}
else{if(U<=0.33){X[isim]<-2
ncomp<<-ncomp+2}
else{if(U<=0.78){X[isim]<-3
ncomp<<-ncomp+3}
else{if(U<=0.91){X[isim]<-4
ncomp<<-ncomp+4}
else{X[isim]<-5
ncomp<<-ncomp+4}}}}
}
return(X)
}
set.seed(500)
nsim <- 10^4
x <- sim1N(nsim)
res <- as.data.frame(table(x)/nsim)
names(res) <- c("x", "psim")
plot(as.matrix(res), col='red',lwd=3, type="h")
points(c(1:5),c(0.11,0.22,0.45,0.13,0.09),pch=22)

cat("Nº medio de comparaciones = ", ncomp/nsim)
## Nº medio de comparaciones = 2.773
#Algoritmo basado en un arbol binario optimo
ncomp <<- 0
sim2N <- function(n=1000){
X <- numeric(n)
for (isim in 1:n) {
U <- runif(1,0,1)
if (U>0.33){if(U>0.78){X[isim]<-2
ncomp <<- ncomp+2}
else{X[isim]<-3
ncomp <<- ncomp+2}}
else{if(U>0.2){X[isim]<-4
ncomp <<- ncomp+2}
else{if(U>0.11){X[isim]<-5
ncomp <<- ncomp+3}
else{X[isim]<-1
ncomp <<- ncomp+3}}}
}
return(X)
return(ncomp)
}
set.seed(500)
nsim <- 10^4
u <- sim2N(nsim)
res <- as.data.frame(table(u)/nsim)
names(res) <- c("u", "psim")
plot(as.matrix(res), col='red',lwd=3,type="h")
points(c(1:5),c(0.11,0.22,0.45,0.13,0.09),pch=18)

cat("Nº medio de comparaciones = ", ncomp/nsim)
## Nº medio de comparaciones = 2.2038
#datos
x <- c(1:5)
p <- c(0.11,0.22,0.45,0.13,0.09)
#Algoritmo basado en Tabla Guia
ncomp <<- 0
rfmp.tabla <- function(x, prob = 1/length(x), m, nsim = 1000) {
# Simulación v.a. discreta a partir de función de masa de probabilidad
# por tabla guia de tamaño m
# Inicializar tabla y FD
Fx <- cumsum(prob)
g <- rep(1,m)
i <- 1
for(j in 2:m) {
while (Fx[i] < (j-1)/m) i <- i+1
g[j] <- i
}
# Generar valores
X <- numeric(nsim)
U <- runif(nsim)
for(j in 1:nsim) {
i <- g[floor(U[j]*m)+1]
while (Fx[i] < U[j]) i <- i + 1
X[j] <- x[i]
ncomp <<- ncomp+1
}
return(X)
}
set.seed(123)
nsim <- 10^4
v <- rfmp.tabla(x,p,4,nsim)
res <- as.data.frame(table(v)/nsim)
names(res) <- c("v", "psim")
plot(as.matrix(res), col='red',lwd=3, type="h")
points(c(1:5),c(0.11,0.22,0.45,0.13,0.09),pch=5)

cat("Nº medio de comparaciones = ", ncomp/nsim)
## Nº medio de comparaciones = 1
#Notemos que el algoritmo mejora considerablemente al usar la tabla guia
- El número de procesos de usuarios activos en una estación de trabajo es una variable, \(X\), con distribución dada por: \[
P(X=0)=0.33, \space P(X=1)=0.22, \space P(X=2)=0.15,\\ \space P(X=3)=0.11, \space P(X=4)=0.09, \space P(X=5)=0.06,\\ \space P(X=6)=0.04.
\] Dar un algoritmo para simular esta variable mediante búsqueda en forma de árbol de Huffman. Calcular el número medio de comparaciones para simular cada valor de la variable \(X\). ¿Cuál es una cota para dicho número medio de comparaciones en caso de usar el método de la tabla guía?. Para este último método, calcular el valor exacto del número medio de comparaciones.
#datos
s <- c(0:6)
p <- c(0.33,0.22,0.15,0.11,0.09,0.06,0.04)
ncomp0 <<- 0
ncomp1 <<- 0
ncomp2 <<- 0
ncomp3 <<- 0
ncomp4 <<- 0
ncomp5 <<- 0
ncomp6 <<- 0
sim21N <- function(n=10000){
X <- numeric(n)
for(i in 1:n){
U <- runif(1)
if(U>0.45){if(U>0.78){X[i]<-1
ncomp1<<-ncomp1+2}
else{X[i]<-0
ncomp0<<-ncomp0+2}}
else{if(U>0.3){X[i]<-2
ncomp2<<-ncomp2+3}
else{if(U>0.19){X[i]<-3
ncomp3<<-ncomp3+4}
else{if(U>0.1){X[i]<-4
ncomp4<<-ncomp4+5}
else{if(U>0.04){X[i]<-5
ncomp5<<-ncomp5+6}
else{X[i]<-6
ncomp6<<-ncomp6+6}}}}}
}
return(X)
}
set.seed(500)
nsim <- 10^4
u <- sim21N(nsim)
res <- as.data.frame(table(u)/nsim)
names(res) <- c("u", "psim")
plot(as.matrix(res), col='red',lwd=3,type="h")
points(s,p,pch=18)

cat("Número medio de comparaciones para 0 =", ncomp0/(nsim*p[1]))
## Número medio de comparaciones para 0 = 1.970303
cat("\n Número medio de comparaciones para 1 =", ncomp1/(nsim*p[2]))
##
## Número medio de comparaciones para 1 = 2.02
cat("\n Número medio de comparaciones para 2 =", ncomp2/(nsim*p[3]))
##
## Número medio de comparaciones para 2 = 2.976
cat("\n Número medio de comparaciones para 3 =", ncomp3/(nsim*p[4]))
##
## Número medio de comparaciones para 3 = 4.010909
cat("\n Número medio de comparaciones para 4 =", ncomp4/(nsim*p[5]))
##
## Número medio de comparaciones para 4 = 4.888889
cat("\n Número medio de comparaciones para 5 =", ncomp5/(nsim*p[6]))
##
## Número medio de comparaciones para 5 = 6.48
cat("\n Número medio de comparaciones para 6 =", ncomp6/(nsim*p[7]))
##
## Número medio de comparaciones para 6 = 6.12
#Ahora, usando el método de la tabla guia
ncomp <<- 0
rfmp.tabla <- function(x, prob = 1/length(x), m, nsim = 1000) {
# Simulación v.a. discreta a partir de función de masa de probabilidad
# por tabla guía de tamaño m
# Inicializar tabla y FD
Fx <- cumsum(prob)
g <- rep(1,m)
i <- 1
for(j in 2:m) {
while (Fx[i] < (j-1)/m) i <- i+1
g[j] <- i
}
# Generar valores
X <- numeric(nsim)
U <- runif(nsim)
for(j in 1:nsim) {
i <- g[floor(U[j]*m)+1]
while (Fx[i] < U[j]) i <- i + 1
X[j] <- x[i]
ncomp <<- ncomp+1
}
return(X)
}
set.seed(500)
nsim <- 10^4
v <- rfmp.tabla(s,p,6,nsim)
res <- as.data.frame(table(v)/nsim)
names(res) <- c("v", "psim")
plot(as.matrix(res), col='red',lwd=3,type="h")
points(s,p,pch=17)

cat("Nº medio de comparaciones = ", ncomp/nsim)
## Nº medio de comparaciones = 1
#Como se esperaba el numero de comparaciones es 1.
- El número deveces que se produce una caida de un sitema informático durante un mes es una variable aleatoria, \(X\), con distribución dada por \(P(X=0)=0.09, \space P(X=1)=0.21, \space P(X=2)=0.39, \space P(X=3)=0.19, \space P(X=4)=0.08, \space P(X=5)=0.04.\) Dar un algoritmo para simular esta variable mediante búsqueda en forma de árbol de Huffman. Calcular el número medio de comparaciones para simular cada valor de la variable \(X\). ¿Qué ventajas se obtienen para este método de simulación en comparación con la implementación mediante el método secuencial (tanto en su versión directa como en la del etiquetado más eficiente)?
#Datos
s <- c(0:5)
p <- c(0.09,0.21,0.39,0.19,0.08,0.04)
ncomp0 <<- 0
ncomp1 <<- 0
ncomp2 <<- 0
ncomp3 <<- 0
ncomp4 <<- 0
ncomp5 <<- 0
sim31N <- function(n=1000){
X <- numeric(n)
for (i in 1:n) {
U <- runif(1)
if(U>0.4){if(U>0.79){X[i]<-1
ncomp1 <<- ncomp1+2}
else{X[i]<-2
ncomp2 <<- ncomp2+2}}
else{if(U>0.21){X[i]<-3
ncomp3 <<- ncomp3+3}
else{if(U>0.12){X[i]<-0
ncomp0 <<- ncomp0+4}
else{if(U>0.08){X[i]<-5
ncomp5 <<- ncomp5+5}
else{X[i]<-4
ncomp4 <<- ncomp4+5}}}
}
}
return(X)
}
set.seed(500)
nsim <- 10^4
system.time(u <- sim31N(nsim))
## user system elapsed
## 0.17 0.00 0.31
res <- as.data.frame(table(u)/nsim)
names(res) <- c("u", "psim")
plot(as.matrix(res), col='red',lwd=3,type="h")
points(s,p,pch=18)

cat("Numero de comparaciones para 0 =", ncomp0/(nsim*p[1]))
## Numero de comparaciones para 0 = 3.973333
cat("\n Numero de comparaciones para 1 =", ncomp1/(nsim*p[2]))
##
## Numero de comparaciones para 1 = 2.032381
cat("\n Numero de comparaciones para 2 =", ncomp2/(nsim*p[3]))
##
## Numero de comparaciones para 2 = 1.962051
cat("\n Numero de comparaciones para 3 =", ncomp3/(nsim*p[4]))
##
## Numero de comparaciones para 3 = 3.004737
cat("\n Numero de comparaciones para 4 =", ncomp4/(nsim*p[5]))
##
## Numero de comparaciones para 4 = 5.0875
cat("\n Numero de comparaciones para 5 =", ncomp5/(nsim*p[6]))
##
## Numero de comparaciones para 5 = 5.3625
#Ahora, usando el método de busqueda secuencial
rfmp <- function(x, prob = 1/length(x), nsim = 1000) {
# Simulación nsim v.a. discreta a partir de fmp
# por inversión generalizada (transformación cuantil)
# Inicializar FD
Fx <- cumsum(prob)
# Simular
X <- numeric(nsim)
U <- runif(nsim)
for(j in 1:nsim) {
i <- 1
while (Fx[i] < U[j]) i <- i + 1
X[j] <- x[i]
ncomp <<- ncomp + i
}
return(X)
}
ncomp <<- 0
system.time(v <- rfmp(s,p,nsim))
## user system elapsed
## 0.06 0.00 0.06
res <- as.data.frame(table(u)/nsim)
names(res) <- c("u", "psim")
plot(as.matrix(res), col='red',lwd=3,type="h")
points(s,p,pch=18)

cat("Nº medio de comparaciones = ", ncomp/nsim)
## Nº medio de comparaciones = 3.0997
#Método secuencial
so <- c(2,1,3,0,4,5)
po <- c(0.39,0.21,0.19,0.09,0.08,0.04)
ncomp <<- 0
system.time(v <- rfmp(so,po,nsim))
## user system elapsed
## 0.01 0.00 0.01
res <- as.data.frame(table(u)/nsim)
names(res) <- c("u", "psim")
plot(as.matrix(res), col='red',lwd=3,type="h")
points(s,p,pch=18)

cat("Nº medio de comparaciones = ", ncomp/nsim)
## Nº medio de comparaciones = 2.3544
#Notemos que el método de arbol binario es muy lento y por tanto los otros métodos son mas eficientes
- El código de error que puede siministrar una rutina implementada en un lenguaje de alto nivel puede describirse mediante una varibale aleatoria, \(C\), que puede tomar los valores 0,1,2,4 y 7, con probabilidades: \[
p_0=0.28, \space p_1=0.35, \space p_3=0.16,\space p_4=0.12, \space p_7=0.09.
\] Dentro de un complejo modelo de simulación se quiere generar dicha variable usando, para cada valor generado, un único número pseudoaleatorio \(U(0,1)\). Encontrar algoritmos, para tal fin, basados en la construcción de un árbol Huffman y en el método de las tablas guías. Estudiar la eficiencia de ambos algoritmos y compararla.
#Datos
s <- c(0,1,3,4,7)
p <- c(0.28,0.35,0.16,0.12,0.09)
sim41N <- function(n=10000){
X<-numeric(n)
for (i in 1:n) {
U <- runif(1)
if(U>0.37){if(U>0.72){X[i]<-0}
else{X[i]<-1}}
else{if(U>0.21){X[i]<-3}
else{if(U>0.12){X[i]<-7}
else{X[i]<-4}}}
}
return(X)
}
set.seed(500)
nsim <- 10^4
system.time(u <- sim41N(nsim))
## user system elapsed
## 0.14 0.00 0.14
res <- as.data.frame(table(u)/nsim)
names(res) <- c("u", "psim")
plot(as.matrix(res), col='red',lwd=3,type="h")
points(s,p,pch=18)

ncomp <<- 0
rfmp.tabla <- function(x, prob = 1/length(x), m, nsim = 1000) {
# Simulación v.a. discreta a partir de función de masa de probabilidad
# por tabla guía de tamaño m
# Inicializar tabla y FD
Fx <- cumsum(prob)
g <- rep(1,m)
i <- 1
for(j in 2:m) {
while (Fx[i] < (j-1)/m) i <- i+1
g[j] <- i
}
# Generar valores
X <- numeric(nsim)
U <- runif(nsim)
for(j in 1:nsim) {
i <- g[floor(U[j]*m)+1]
while (Fx[i] < U[j]) i <- i + 1
X[j] <- x[i]
ncomp <<- ncomp+1
}
return(X)
}
set.seed(123)
nsim <- 10^4
system.time(v <- rfmp.tabla(s,p,4,nsim))
## user system elapsed
## 0.08 0.00 0.09
res <- as.data.frame(table(v)/nsim)
names(res) <- c("v", "psim")
plot(as.matrix(res), col='red',lwd=3,type="h")
points(s,p,pch=17)

cat("Nº medio de comparaciones = ", ncomp/nsim)
## Nº medio de comparaciones = 1
#Ambos procesos son muy similares
- Considérese la variable aleatoria discreta, \(X\), número de fallos mensuales en el suministro eléctrico de un equipo informático, que tiene la siguiente masa de probablilidad \[
P_0=0.35, \space P_1=0.20, \space P_2=0.16, \space P_3=0.12, \space P_4=0.08, \space P_5=0.05, \space P_6=0.03, \space P_7=0.01.
\] Construir una tabla guía y, basándose en ella, dar un algoritmo para simular esta variable. ¿Cuál es el número medio exacto de comparaciones para generar un valor mediante ese algoritmo?
#Datos
s <- c(0:7)
p <- c(0.35,0.20,0.16,0.12,0.08,0.05,0.03,0.01)
ncomp <<- 0
rfmp.tabla <- function(x, prob = 1/length(x), m, nsim = 1000) {
# Simulación v.a. discreta a partir de función de masa de probabilidad
# por tabla guia de tamaño m
# Inicializar tabla y FD
Fx <- cumsum(prob)
g <- rep(1,m)
i <- 1
for(j in 2:m) {
while (Fx[i] < (j-1)/m) i <- i+1
g[j] <- i
}
# Generar valores
X <- numeric(nsim)
U <- runif(nsim)
for(j in 1:nsim) {
i <- g[floor(U[j]*m)+1]
while (Fx[i] < U[j]) i <- i + 1
X[j] <- x[i]
ncomp <<- ncomp+1
}
return(X)
}
set.seed(500)
nsim <- 10^4
system.time(v <- rfmp.tabla(s,p,7,nsim))
## user system elapsed
## 0.07 0.00 0.08
res <- as.data.frame(table(v)/nsim)
names(res) <- c("v", "psim")
plot(as.matrix(res), col='red',lwd=3,type="h")
points(s,p,pch=17)

cat("Nº medio de comparaciones = ", ncomp/nsim)
## Nº medio de comparaciones = 1
- Dando por supuesto un algoritmo que permita simular la variable aleatoria del ejercicio anterior, deducir un método para aproximar por simulación la probabilidad de que en un año haya más de 30 fallos en el suministro eléctrico. Plantear el cálculo analítico de la cantidad anterior. ¿Sería sencillo obtener teóricamente la masa de probabilidad de la variable aleatoria número de fallos anuales?
#Datos
s <- c(0:7)
p <- c(0.35,0.20,0.16,0.12,0.08,0.05,0.03,0.01)
ncomp <<- 0
rfmp.tabla <- function(x, prob = 1/length(x), m, nsim = 1000) {
# Simulación v.a. discreta a partir de función de masa de probabilidad
# por tabla guia de tamaño m
# Inicializar tabla y FD
Fx <- cumsum(prob)
g <- rep(1,m)
i <- 1
for(j in 2:m) {
while (Fx[i] < (j-1)/m) i <- i+1
g[j] <- i
}
# Generar valores
X <- numeric(nsim)
U <- runif(nsim)
for(j in 1:nsim) {
i <- g[floor(U[j]*m)+1]
while (Fx[i] < U[j]) i <- i + 1
X[j] <- x[i]
ncomp <<- ncomp+1
}
return(X)
}
set.seed(123)
n <- 12
nsim <- 1000
X<<-numeric(nsim)
for (i in 1:n) {
u <- rfmp.tabla(s,p,7,nsim)
X <<- X+u
}
probabilidad <- function(X,nsim){
aux1 <-0
for (i in 1:nsim) {
if(X[i]>=30){aux1 <- aux1+1}
}
return(aux1/nsim)
}
probabilidad(X,nsim)
## [1] 0.07