Algoritmos Probabilísticos
(Parte III)
Algoritmos Probabilísticos Numéricos
„
A agulha de Buffon
No século XVIII Georges
Louis Leclerc, conde de
Buffon, demonstrou que se
lançarmos
aleatoriamente,
uma agulha de comprimento
L sobre uma superfície com
linhas paralelas separadas a
uma distância D, igual ou
maior que comprimento da
agulha, a probabilidade de
que a agulha cruze uma das
linhas é:
2L
πD
Algoritmos Probabilísticos
2
Algoritmos Probabilísticos Numéricos
A agulha de Buffon
„
Caso de L=D
rSenAng<-function(){
repeat{
x <-runif(1,-1,1)
y <-runif(1)
h<-sqrt(x^2+y^2)
if (h<=1) break
}
if(x==0) SinTheta<-1
else SinTheta<-y/h
return(SinTheta)
}
Algoritmos Probabilísticos
AgulhaBuffon<-function(n,d){
cruz<-0
for (i in 1:n){
if(runif(1,0,d/2)<d/2*rSenAng())
cruz<-cruz+1
}
return(2*n/cruz)
}
3
Algoritmos Probabilísticos Numéricos (cont.)
„
Integração Numérica
… Suponha que pretende calcular
b
I = ∫ f ( x)dx
a
Considere um rectângulo de largura (b-a) e altura I/(b-a)
… A área deste rectângulo é I
… Uma vez que o rectângulo e a superfície abaixo da curva f(x)
têm a mesma área e a mesma largura então devem ter a
mesma altura média
… Concluímos assim que a altura média da curva entre a e b é
I/(b-a)
…
Algoritmos Probabilísticos
4
Algoritmos Probabilísticos Numéricos
Integração Numérica
MCint<-function(f,n,a,b){
soma<-0
for(i in 1:n){
x<-runif(1,a,b)
soma<-soma+eval(f)
}
areaInt<-(b-a)*(soma/n)
return(areaInt)
}
Obs.: Verifica-se que a variância da estimativa calculada por este
algoritmo é proporcional ao número de pontos gerados e que a
distribuição da estimativa é aproximadamente Normal para valores
de n elevados.
Algoritmos Probabilísticos
5
Algoritmos Probabilísticos Numéricos
Integração Numérica (cont.)
„
Problema: Uma outra forma de estimar o valor de π
é usar a integração de Monte Carlo recorrendo-se à
seguinte relação:
π =∫
2
0
Algoritmos Probabilísticos
4 − x dx
2
6
Algoritmos de Monte Carlo
„
Exemplo 1: Verificar a Multiplicação de
Matrizes
Input: Três matrizes A, B e C (nxn)
Ouput: Verdade se C=AB ou Falso, caso contrário
Algoritmos Probabilísticos
7
Algoritmos de Monte Carlo: Exemplo 1 (cont.)
verificarMult<-function(A,B,C,n){
X<-integer(n)
for (j in 1:n) X[j]<-uniformInt(0,1)
if (identical((X%*%A)%*%B,X%*%C)))
return(TRUE)
else
return(FALSE)
}
Nota: Este algoritmo é ½-correct
Algoritmos Probabilísticos
8
Algoritmos de Monte Carlo: Exemplo 1 (cont.)
RepVerificarMult<-function(A,B,C,n,k){
for (i in 1:k)
if (verificarMult(A,B,C,n)==FALSE)
return(FALSE)
return(TRUE)
}
Questão: Neste caso, qual o valor de p?
Algoritmos Probabilísticos
9
Algoritmos de Monte Carlo: Exemplo 1 (cont.)
RepVerificarMultErro<-function(A,B,C,n,e){
k<-ceiling(log(1/e,base=2))
return(RepVerificarMult(A,B,C,n,k))
}
Algoritmos Probabilísticos
10
Algoritmos de Monte Carlo
„
Exemplo 2: Verificar se um número é primo
Input: Um número inteiro n
Ouput: Verdade se n é primo ou Falso se o número é
composto
Algoritmos Probabilísticos
11
Algoritmos de Monte Carlo: Exemplo 2 (cont.)
Teorema (Fermat’s little theorem)
Seja n um número primo. Então
an-1 mod n = 1
para qualquer inteiro a tal que 1≤ a ≤ n-1
Algoritmos Probabilísticos
12
Algoritmos de Monte Carlo: Exemplo 2 (cont.)
Função para calcular an mod z
expomod<-function(a,n,z){
i<-n; r<-1; x<-a%%z
while(i>0){
if (i%%2==1) r<-(r*x)%%z
x<-x^2%%z
i<-i%/%2
}
return(r)
}
Algoritmos Probabilísticos
13
Algoritmos de Monte Carlo: Exemplo 2 (cont.)
Seja n um número maior ou igual a 2.
fermat<-function(n){
a<-uniformInt(1,n-1)
if (expomod(a,n-1,n)==1)
return(TRUE)
else
return(FALSE)
}
Algoritmos Probabilísticos
14
Algoritmos de Monte Carlo: Exemplo 2 (cont.)
Definição de conjunto B(n):
Seja n um inteiro impar maior que 4 e sejam s e t inteiros
tais que n-1=2st, com t impar. Seja B(n) o conjunto dos
inteiros definidos da seguinte forma:
a∈B(n) se e só se 2≤a≤n-2 e
„ at mod n=1
ou
„ existe um inteiro i, 0≤i<s, tal que
a
Algoritmos Probabilísticos
2i t
mod n = n − 1
15
Algoritmos de Monte Carlo: Exemplo 2 (cont.)
„
Uma extensão ao teorema de Fermat mostra
que
… a∈B(n)
„
para todo 2≤a≤n-2 quando é primo
Diz-se que n é um pseudoprimo forte (strong
pseudoprime) em relação à base a e que a é um
strong false witness de primalidade para n,
sempre que n>4 for um inteiro composto impar e
a∈B(n)
Algoritmos Probabilísticos
16
Algoritmos de Monte Carlo: Exemplo 2 (cont.)
Btest<-function(a,n){
# n tem de ser um número impar maior que 4
s<-0; t1<-n-1
repeat{
s<-s+1; t1<-t1%/%2
if (t1%%2==1) break
}
x<-expomod(a,t1,n)
if (x==1 || x==n-1) return(TRUE)
for(i in 1:(s-1)){
x<-x^2%%n
if (x==n-1) return(TRUE)
}
return(FALSE)
}
Algoritmos Probabilísticos
17
Algoritmos de Monte Carlo: Exemplo 2 (cont.)
„
A função Btest(a,n) retorna verdade quando n é um
número primo maior que 4 e 2≤a≤n-2, enquanto retorna
falso com probabilidade superior a ¾ quando n é um
número composto impar maior que 4 e a é escolhido
aleatoriamente entre 2 e n-2
„
Teste de Miller-Rabin
MillRab<-function(n){
a<-uniformInt(2,n-2)
return(Btest(a,n))
}
Algoritmos Probabilísticos
18
Algoritmos de Monte Carlo: Exemplo 2 (cont.)
„
Porque a resposta falsa é garantidamente
correcta, podemos reduzir a probabilidade do
erro do algoritmo repetindo o mesmo diversas
vezes
RepMillRab<-function(n,k){
for (i in 1:k)
if(MillRab(n)==FALSE) return(FALSE)
return(TRUE)
}
Algoritmos Probabilísticos
19