¡Descarga practica estadistica y más Ejercicios en PDF de Estadística solo en Docsity!
Estad´ıstica
Pr´actica 5: Simulaci´on de variables aleatorias
J. Gibergans y M. Ortego - DMA III/EUETIB
Curso 2010-2011/T
Resumen En esta pr´actica se introduce la simulaci´on de variables aleatorias como una poderosa t´ecnica para la resoluci´on de problemas. Las primeras referencias sobre simulaci´on se en- cuentran hacia el a˜no 1940 cuando Von Neumann y Ullman trabajaban sobre la simulaci´on del flujo de neutrones para la construcci´on de la bomba at´omica en el proyecto secreto MonteCarlo. En aquel entonces, las t´ecnicas de simulaci´on se conoc´ıan como procesos de MonteCarlo, pero en la actualidad se diferencian ambas cosas, siendo el segundo un tipo particular de simulaci´on.
1. Objetivos, conocimientos previos y evaluaci´on
1.1. Objetivos
Al finalizar esta pr´actica, el alumno ha de ser capaz de:
Saber generar valores de variables aleatorias uniformes y utilizarlos para la resoluci´on de problemas de probabilidad.
Aplicar el m´etodo de la transformada inversa para generar variables aleatorias discretas y continuas.
Simular sistemas de componentes y resolver problemas probabil´ısticos relacionados con la fiabilidad.
Simular sistemas f´ısicos y resolver problemas de la mec´anica desde el punto de vista probabil´ıstico.
Implementar los puntos anteriores con el software R.
1.2. Conocimientos previos
Para realizar esta pr´actica, es necesario que el alumno haya estudiado y realizado ejercicios del tema de variables aleatorias. As´ı mismo, debe haber le´ıdo el siguiente resumen te´orico y haber trabajado los tutoriales resolviendo los ejercicios propuestos.
2. Resumen te´orico
2.1. Generaci´on de n´umeros aleatorios
El paso inicial de los experimentos de simulaci´on es la selecci´on de una fuente de n´umeros aleatorios y la obtenci´on de los mismos. As´ı pues, disponer de un buen generador de n´umeros aleatorios es fundamental en Simulaci´on.
En un principio, los n´umeros aleatorios eran generados de manera manual o mec´anica utilizando ruedas giratorias, lanzamientos de dados o barajas de cartas. En la actualidad se hacen servir ordenadores para llevar a cabo este trabajo.
Estos n´umeros constituyen una sucesi´on de valores que, aunque son obtenidos de forma de- terminista, tienen toda la apariencia de variables aleatorias uniformes U (0, 1) e independientes reciben el nombre de n´umeros pseudoaleatorios aunque en muchas ocasiones simplemente se denominan n´umeros aleatorios.
Uno de los m´etodos m´as comunes para obtener estos n´umeros aleatorios es el m´etodo conguencial. Comienza con un valor inicial x 0 y despu´es va calculando de forma recursiva los valores sucesivos xn, n ≥ 1, haciendo
xn = axn− 1 modulo m
donde a y m son enteros positivos. La cantidad xn/m se considera como una aproximaci´on del valor de una variable aleatoria U (0, 1).
Ejemplo 1. En las subroutinas IMSL de FORTRAN el m´etodo congruencial para el c´alculo de n´umeros pseudoaleatorios utiliza a = 39720494, m = 2^31 − 1.
No es objetivo de estas pr´acticas el c´alculo de n´umeros aleatorios de manera que, para encontrar estos valores, haremos servir el programa R.
2.2. Generaci´on de variables aleatorias discretas
El objetivo de este apartado es el de generar una variable aleatoria discreta X, que tome valores xi con probabilidades pi = P [X = xi]. El procedimiento consiste en dividir el intervalo (0, 1) en tantas partes como valores tome la variable aleatoria X, de manera que podemos establecer una correspondencia entre la longitud del subintervalo Ii y la probabilidad asociada a ´el, pi. Posteriormente, se asigna a la variable aleatoria X el valor xi si el valor U simulado de la variable aleatoria uniforme U (0, 1) es tal que:
∑^ i−^1
k=
pk ≤ U <
∑^ i
k=
pk
Si los xi estan ordenados de manera que x 0 < x 1 <... y si F es la funci´on de distribuci´on de X, entonces F (xk) =
∑k i=0 pi^ y tenemos que
F (xj− 1 ) ≤ U < F (xj ) ⇒ X = xj Es decir, generamos el valor aleatorio U y luego determinamos el valor de X encontrando el intervalo (F (xj− 1 ), F (xj )) en el que se encuentra U ; o lo que es lo mismo, encontrando la inversa de F (U ). Esta forma de actuar justifica el nombre de la transformada inversa que recibe el m´etodo. En la figura 1 se muestra su interpretaci´on gr´afica.
2.3. Generaci´on de variables aleatorias continuas
Un m´etodo an´alogo es aplicable a variables aleatorias continuas con funci´on de distribuci´on invertible, por lo que tambi´en se denomina m´etodo de la transformada inversa.
As´ı pues, para generar una variable aleatoria X a partir de la funci´on de distribuci´on continua F , generaremos un n´umero aleatorio U y calculamos:
X = F −^1 (U )
tal y como se muestra en la figura 2.
sample(c(a,b),n,rep=T,prob=c(pa,pb)): remuestrea al azar y con reemplazo n ele- mentos entre los dos sucesos a y b, donde la probabilidad del suceso a es pa y la del suceso b es pb.
sapply: instrucci´on utilizada para aplicar funciones a las filas o columnas de una matriz o dataframe.
sum(x): proporciona la suma de los elementos de x.
table(x): devuelve una tabla con los distintos valores de los elementos de x y sus corre- spondientes frecuencias absolutas.
Tambi´en se utilizar´an algunos operadores comparativos:
< menor que
mayor que
<= menor o igual que
= mayor o igual que
== igual
!= diferente de
y operadores l´ogicos:
x & y Y l´ogico: x Y y
x | y O l´ogico: x O y
- Tutoriales
4.1. N´umeros aleatorios, dados y urnas
La funci´on sample es de gran utilidad para efectuar simulaciones
sample(1:10)
extrae del vector 1 : 10 un elemento al azar.
Problema 4.1.1. ¿Qu´e resultados proporcionan las instrucciones:
a) sample(1:6,3)? Rep´ıtela varias veces para observar qu´e ocurre.
b) sample(1:6,3,rep=T)? Rep´ıtela varias veces para observar qu´e ocurre. ¿Qu´e diferencia hay con la instrucci´on anterior?
Para hallar la frecuencia de caras en el lanzamiento de un dado, asignaremos cara = 1 y cruz = 0, y podremos hacer la simulaci´on del lanzamiento de 10000 dados mediante la instrucci´on:
sum(sample(0:1,10000,rep=T))/
Si la variable aleatoria puede tomar m´as de dos resultados, entonces la estrategia anterior no nos sirve. La funci´on table proporciona una tabla de frecuencias a partir de un vector. Por ejemplo, si queremos obtener la tabla de frecuencias absolutas correspondiente a la simulaci´on de 100 lanzamientos de un dado hacemos:
t<-table(sample(1:6,100,rep=T))
Dado que la puntuaci´on obtenida en el lanzamiento de un dado es una variable aleatoria discreta, podemos obtener su diagrama de barras haciendo
barplot(t)
Problema 4.1.2. Interpreta el diagrama de barras obtenido. ¿Qu´e distribuci´on sigue la variable aleatoria correspondiente a la puntuaci´on obtenida en el lanzamiento de un dado?
Si ahora consideramos la variable aleatoria X =suma de las puntuaciones obtenidas en el lanzamiento de tres dados, podemos simular valores de dicha variable mediante la instrucci´on:
y<-sapply(1:10000,function(x){sum(sample(1:6,3,rep=T))})
donde la funci´on sapply aplica a un vector de tama˜no 10000 una funci´on sin nombre, generando a su vez otro vector de tama˜no 10000. La funci´on obtiene muestras del lanzamiento del dado, de tama˜no 3 y seguidamente, suma las puntuaciones obtenidas en la muestra.
Problema 4.1.3. A partir de la muestra obtenida, haz la representaci´on gr´afica correspondi- ente. Utilizando las instrucciones de R que ya conoces de pr´acticas anteriores, calcula la media y la varianza de X.
Problema 4.1.4. Considerando de nuevo la variable aleatoria X =suma de las puntuaciones obtenidas en el lanzamiento de tres dados, calcula E(X) y V ar(X) de forma te´orica. Compara los resultados obtenidos con los obtenidos por simulaci´on en el ejercicio anterior.
Hasta ahora hemos estado considerando sucesos equiprobables. Veamos qu´e podemos hacer cuando no ocurra ´esto. Supongamos que tenemos una urna con 4 bolas blancas (que repre- sentaremos mediante 0) y 6 bolas negras (que representaremos mediante 1). Para simular la extracci´on de una bola hacemos:
sample(c(1,0),1,prob=c(4/10,6/10))
Si queremos simular la extracci´on de 6 bolas con reposici´on:
sample(c(1,0),6,rep=T,prob=c(4/10,6/10))
4.2. El m´etodo de la transformada inversa
En este tutorial vamos a ilustrar c´omo funciona el m´etodo de la transformada inversa.
Problema 4.2.1. Variables aleatorias discretas: Sea una variable aleatoria que toma los valores X = 0, 1 , 2 , 3 , 4 con probabilidades P [X = 0] = 0, 1 ; P [X = 1] = 0, 2 ; P [X = 2] = 0, 25 ; P [X = 3] = 0, 2 ; P [X = 4] = 0, 25. Simular 2000 valores de X, hacer la tabla de frecuencias y el diagrama de barras correspondiente. Justifica una a una las siguientes expresiones escritas con R.
Deberemos ejecutar las siguientes ´ordenes:
n<- x<-runif(n,0,1) y<-numeric(n) y[x<0.10]<- y[0.10 <=x & x<0.30]<- y[0.30 <=x & x<0.55]<- y[0.55 <=x & x<0.75]<- y[x>=0.75]<- table(y)/n barplot(table(y))
Problema 4.2.2. Variables aleatorias continuas: Sea una variable aleatoria X con funci´on de densidad de probabilidad:
f (x) =
1 + x, − 1 ≤ x < 0 1 − x, 0 ≤ x < 1
Genera 100000 valores aleatorios de X y haz la representaci´on gr´afica correspondiente.
Figura 3: Aproximando el valor de π.
n<- x<-runif(n,-1,1) y<-runif(n,-1,1) t<-numeric(n) t[x^2+y^2<=1]<- 4*sum(t)/n
Observa que aparece la instrucci´on runif que nos permite obtener valores de una variable aleatoria uniforme en un intervalo determinado.
Problema 4.3.1. Interpreta las instrucciones anteriores. ¿Qu´e cambios realizar´ıas para hallar el valor de π con m´as precisi´on?
4.4. El Teorema de Bayes
Una f´abrica produce piezas por medio de dos m´aquinas A y B. La probabilidad que una pieza producida por la m´aquina sea defectuosa es del 5 %, mientras que si la ha producido la m´aquina B esta probabilidad es del 7 %. La producci´on de la m´aquina A es 4 veces la de la m´aquina B. Si se toma una pieza de la producci´on total y resulta que no es defectuosa, ¿cu´al es la probabilidad que haya sido producida por la m´aquina A?
Problema 4.4.1. Resuelve este problema de forma exacta aplicando el Teorema de Bayes.
Tambi´en puede ser resuelto de forma aproximada por simulaci´on:
n<- nA<-n0. nB<-n0. produccion<-c(rep("A",nA),rep("B",nB)) piezas<-c("D","ND") fA<-sample(piezas,nA,rep=TRUE,prob=c(0.05,0.95)) fB<-sample(piezas,nB,rep=TRUE,prob=c(0.07,0.93)) calidad<-c(fA,fB) ADn<-numeric(nA+nB) ADn[produccion=="A" & calidad=="ND"]<- Dn<-numeric(nA+nB) Dn[calidad=="ND"]<- sum(ADn)/sum(Dn)
Problema 4.4.2. Interpreta estas instrucciones y compara los resultados obtenidos mediante ambos m´etodos.
4.5. Fiabilidad (I)
Considera el sistema de la figura 4, formado por componentes en serie y en paralelo, en el que cada componente puede funcionar con las siguientes probabilidades: A: 0.9, B: 0.8, C:0.95, D:0.9, E: 0.9, F:0.5. ¿Cu´al es la probabilidad que el sistema funcione?
Figura 4: Sistema con componentes en serie y en paralelo.
Podemos simular un gran n´umero de circuitos y ver la proporci´on de ellos que funcionan. Esto nos dar´a una probabilidad aproximada de que el sistema funcione:
n<- A<-runif(n,0,1) B<-runif(n,0,1) C<-runif(n,0,1) D<-runif(n,0,1) E<-runif(n,0,1) F<-runif(n,0,1) CIRC<-numeric(n) CIRC[(A<0.9|B<0.8)&(C<0.95)&(D<0.9|E<0.9|F<0.5)]<- sum(CIRC)/n
Problema 4.5.1. Resuelve el problema de forma anal´ıtica y compara el resultado con el obtenido mediante simulaci´on.
