Mostrando postagens com marcador combinatória. Mostrar todas as postagens
Mostrando postagens com marcador combinatória. Mostrar todas as postagens

sexta-feira, 16 de abril de 2021

As Figurinhas dos Vingadores são Honestas? (parte 2)

No último post do blog, eu comprei vários pacotinhos de figurinhas dos Vingadores, para tentar responder à pergunta: as figurinhas são honestas? Após fazer dois testes de hipótese, eu concluí que as figurinhas possuem distribuição uniforme, então parece que sim. Mas tem uma pegadinha: para serem honestas, garantir que são uniformes não é o suficiente, ainda tem outro teste que precisamos fazer!


A Distribuição dos Pacotinhos


Para fazer a análise das figurinhas, eu comprei 106 pacotinhos na banca. Mas enquanto eu abria os pacotinhos, notei uma propriedade curiosa: todos os pacotinhos tinham sempre 4 figurinhas diferentes. Mas olha só, se os pacotinhos não possuem repetidas entre si, então as figurinhas podem até ser uniformes, mas não são independentes! Cada pacotinho tem quatro figurinhas, e se eu tirei a figurinha 42 em um pacotinho, então eu sei que as outras três não são a figurinha 42.

Mas pode ter sido só o acaso né? Vai ver eu dei sorte de comprar exatamente 106 pacotinhos que não tinham repetidas. Mas felizmente nós sabemos como testar se foi sorte ou não, basta fazer o teste do ##\chi^2##, que vimos no último post!

Acho que o jeito mais fácil de testar é contar quantas figurinhas diferentes tem dentro do pacotinho. Quantos pacotinhos possíveis existem com 4 figurinhas diferentes? Esse é fácil, o álbum tem 180 figurinhas, então são ##180\choose 4## pacotinhos. Com 3 diferentes precisa de mais cuidado, existem ##180\choose 3## combinações de figurinhas, mas você precisa contar qual dessas três diferentes está repetida (os casos são AABC, ABBC, ABCC). Então para 3 diferentes são ##3 {180\choose 4}## pacotinhos.

Com 2 diferentes também tem três casos: AAAB, AABB, ABBB, logo temos ##3{180\choose 2}## pacotinhos. Para uma única figurinha diferente, são apenas ##180\choose 1## pacotinhos. Agora podemos montar a tabela com quantos pacotinhos eu deveria ter visto, e quantos pacotinhos eu vi:

1 diferente2 diferentes3 diferentes4 diferentes
Total de pacotinhos##180\choose 1####3{180\choose 2}####3{180\choose 3}####180\choose 4##
Total (numérico)18048330286758042296805
Total (normalizado)3.98e-61.07e-30.6340.936
Esperado (total*106)422e-60.1136.7299.163
Observado000106

Opa, nas duas primeiras categorias ##np\lt 5##, então eu não tenho amostras suficientes para usar o ##\chi^2##! Um jeito de consertar isso sem comprar mais pacotinhos é juntando as três primeiras categorias em uma só:

Menos de 4 diferentes4 diferentes
Esperado6.83799.163
Observado0106

Com isso conseguimos calcular o ##\chi^2## pela definição, chegando no valor de ##\chi^2=7.308##. Qual é o limite para rejeitar a hipótese? Se estivermos trabalhando com ##\alpha=5\%##, e um grau de liberdade, então o limite é:
InverseCDF[ChiSquareDistribution[1], 1 - 0.05]
>> 3.84146
Como ##7.308\gt 3.84146##, então temos que rejeitar a hipótese, e concluir que as figurinhas dos Vingadores não são independentes, e, como consequência, não são honestas, apesar de serem uniformes!


Felizmente para o colecionador, elas são desonestas mas estão jogando do seu lado. Se os pacotinhos diminuem a chance de ter repetidas, então o total esperado de pacotinhos que você precisa comprar é menor que no caso puramente honesto.

Dá para calcular o quanto você economiza com as figurinhas desonestas? Sure, vou mostrar dois métodos para fazer isso:

Cadeias de Markov


Esse primeiro método é mais complicado, mas vou começar por ele porque é mais genérico, e se aplica a mais casos. Imagine que eu já tenho 50 figurinhas no meu álbum, e eu abro mais um pacotinho. Se eu der azar, vão sair quatro figurinhas que eu já tenho. Se eu der sorte, pode sair até quatro que eu ainda não tinha. 

O método mais natural para modelar esse problema é com uma Cadeia de Markov: você considera o estado como sendo quantas figurinhas você tem no álbum, e as transições possíveis são quantas figurinhas você ainda não tinha no pacotinho.

Para montar a matriz de transição, nós precisamos da probabilidade de ir de um estado ##n## até um estado ##m##. Como cada pacotinho tem 4 figurinhas diferentes, então o ##m## precisa ser necessariamente ##n##, ##n+1##, ..., ##n+4##, e outras transições não são possíveis. 

Qual a chance de ir do estado ##n## para ##n+k##? Nós sabemos que as quatro figurinhas do pacotinho são diferentes entre si, então temos ##k## figurinhas novas e ##4-k## figurinhas que já tínhamos. A quantidade de pacotinhos que fazem essa transição é igual a ##{180-n\choose k}{n\choose 4-k}##, já que precisamos escolher ##k## figurinhas dentre que as não temos, e ##4-k## dentre as que temos. O total de pacotinhos é ##180\choose 4##, então a probabilidade total do estado ##n## para ##n+k## é ##{180-n\choose k}{n\choose 4-k}/{180\choose 4}## . Com isso podemos montar a matriz de transição ##M##:


$$ M= \begin{bmatrix} 0 & 0 & 0 & 0 & 1 & 0 & 0 & 0 & \dots & 0 &  0 & 0 \\ 0 & 0 & 0 & 0 & 0.02 & 0.97 & 0 & 0 & \dots & 0 &  0 & 0 \\ 0 & 0 & 0 & 0 & \approx0 & 0.04 & 0.95 & 0 & \dots & 0 &  0 & 0 \\ 0 & 0 & 0 & 0 & \approx0 & \approx0 & 0.06 & 0.93 & \dots & 0 &  0 & 0 \\ \vdots &\vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & & \vdots & \vdots & \vdots  \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & \dots & 0.95 &  0.04 & \approx0 \\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & \dots & 0 &  0.97 & 0.02 \\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & \dots & 0 &  0 & 1 \\ \end{bmatrix}$$


Pela matriz podemos conferir a intuição: se você não tem nenhuma figurinha, o primeiro pacotinho vai garantidamente te dar 4 figurinhas novas, então o estado 0->4 na matriz vale 1, e todos os outros nessa linha são zero. Por outro lado, se está faltando só uma figurinha para completar, comprar um pacotinho novo dá um chance de 2% de conseguir a faltante. E, claro, se você tem todas as figurinhas, então a chance de você não conseguir mais nenhuma é 100%.

O que queremos calcular é quanto economizamos devido às figurinhas desonestas. Para isso, vamos calcular qual o valor esperado do número de pacotinhos que precisamos comprar para completar o álbum. Tendo a matriz de transição em mãos, isso é equivalente ao número médio de transições para ir do estado 0 (álbum vazio) ao estado 180 (álbum completo). E como calcular isso? Usando um algoritmo secreto!

Really, um algoritmo secreto. Não sei o motivo, mas você não acha esse algoritmo em lugar nenhum; não tem na wikipedia e nem nos livros de processos estocásticos que eu tenho aqui em casa. Então aproveitem para anotar porque é bem útil:

  1. Dada a matriz de transição ##M##, calcule a matriz reduzida ##R##, retirando da matriz original todas as linhas e colunas associadas aos estados finais. (No nosso caso só a última linha e coluna são finais).
  2. Calcule a matriz fundamental ##N## do processo, através da fórmula ##N=(I-R)^{-1}##.
  3. O valor esperado do número de transições é a soma dos elementos da primeira linha da matriz ##N##.
Implementando o algoritmo acima com um scriptinho rápido no Mathematica, concluímos que o valor esperado do número de pacotinhos que você precisa comprar é ##257.98##. Como cada pacotinho custa R$ 1.20, então você precisa de R$309.58 em média para completar o álbum.

E quanto você gastaria com figurinhas 100% honestas? Esse cálculo eu já fiz aqui no blog anteriormente, se as figurinhas são uniformes e independentes então esse é o Coupon Collector Problem: para ##n## figurinhas a solução é ##n H(n)##, onde ##H(n)## é o número harmônico de ##n##. Essa fórmula dá o total de figurinhas, como elas são independentes, para o total de pacotinhos é só dividir por 4, chegando no total de ##259.78## pacotinhos, ou R$ 311.74, uma economia de 2 reais! (Eu falei que economizava, só não falei que era muito).

Embora o algoritmo secreto seja bacana e funcione em vários casos, ele é meio lento. Você precisa inverter uma matriz 181x181, o que não é muito fácil. O fato dela ser triangular e esparsa ajuda, mas tem um método melhor, que calcula o valor em O(n)!

Linearidade da Expectativa


Qual é o teorema mais surpreendente que você conhece? Só em probabilidade tem vários que são difíceis de acreditar à primeira vista, como o Monty Hall Problem e o Paradoxo do Aniversário. Mas um que sempre me espanta é a linearidade da expectativa (o valor esperado da soma é a soma dos valores esperados):

$$E[a+b]=E[a]+E[b]$$

Para mim a parte curiosa é que ele funciona mesmo se as variáveis ##a## e ##b## não forem independentes! Para o nosso problema isso ajuda bastante: o valor esperado para conseguir 180 figurinhas claramente é correlacionado com o valor esperado para conseguir 179 figurinhas, mas eu posso usar relações lineares entre eles sem medo.

Para essa análise eu acho mais fácil pensar invertido. Se estiver faltando 0 figurinhas para completar o album, quantos pacotinhos eu preciso? Nenhum né, ele já está completo. E se estiver faltando só uma figurinha? Aí depende do que sair no pacotinho, se sair a figurinha faltante, um pacotinho basta, se não sair, então você adiciona um pacotinho no seu valor esperado e repete o procedimento:

$$\begin{align*}\\E[1]&=p_0(1+E[0]) + p_1(1+E[1])\\ (1-p_1)E[1] &= p_1 + p_0(1+E[0]) \\ E[1] &=\frac{1}{1-p_1}\left(p_1+p_0(1+E[0])\right)\end{align*}$$

Olha só, eu sei que ##E[0]=0##; para calcular ##E[1]## eu preciso só do ##E[0]##; para o ##E[2]## eu preciso do ##E[1]## e ##E[0]##; e assim por diante. Como cada pacotinho tem no máximo 4 figurinhas que eu não tenho, então eu só preciso olhar quatro valores para trás na recorrência. Portanto, esse algoritmo é linear!

Preenchendo as probabilidades corretas, o algoritmo completo fica assim:

$$\begin{align*} \\ E[0] &= 0 \\ E[n] &= \frac{1}{1-{{180-n}\choose 4}/{180\choose 4}}\left(\frac{{180-n}\choose 4}{180\choose 4}+\sum_{1\le k \le 4}\frac{{n\choose k}{{180-n}\choose{4-k}}}{180\choose 4}E[n-k]\right) \\ \end{align*} $$

Essa implementação realmente é bem mais rápida que o anterior. A versão com matrizes no Mathematica leva 1s para rodar, enquanto que um scriptinho python com o algoritmo linear leva só 55ms!

Agora que temos os pacotinhos corretamente modelados, finalmente vamos conseguir bolar uma estratégia para minimizar o custo de completar o álbum. Mas vai ficar para o próximo post porque esse já ficou comprido demais :)


segunda-feira, 16 de julho de 2018

As Figurinhas dos Vingadores são Honestas? (parte 1)

O Brasil é um país muito estranho. Uma das inúmeras bizarrices brasileiras é que temos tabus temporários. Por exemplo, homem se vestir de mulher é tabu (vai acabar com a família, fim dos tempos, etc), mas no carnaval pode! Da mesma maneira, adulto comprando álbum de figurinhas é tabu (precisa crescer, largar os brinquedos, etc.), mas em época de Copa pode!


Sempre que chega época de Copa eu fico espantado com a dedicação que meus amigos tem para completar o álbum com os jogadores. E pensando um pouco no fenômeno, me veio a dúvida: quanto dinheiro você precisa para terminar um álbum desses?

Resolvi fazer um experimento para tentar calcular esse valor. Mas, ao invés de usar o álbum da Copa, escolhi completar o álbum dos Vingadores, já que é um universo que eu tenho mais afinidade. Os valores originais são os seguintes:

  • O álbum em si custa R$ 6.90, e você precisa de 180 figurinhas para completar. 
  • Cada pacotinho com 4 figurinhas custa R$ 1.20.
  • Existe um kit com 15 pacotinhos que você pode comprar com desconto por R$ 13.90. Mas é difícil de achar: eu só encontrei na loja Geek na av. Paulista.
  • Por fim, você sempre pode pedir figurinhas avulsas no site da editora (nesse caso você pode comprar até 40 figurinhas com preço individual de R$ 0.30, mais o frete de R$ 10.00).

Se você for uma pessoa sem imaginação, de cara já tem uma estratégia fácil para completar o álbum: é só pedir tudo pelo correio. Nesse caso, basta comprar o álbum e mais 180 figurinhas em 5 fretes, num total de R$110.90.

Tem como gastar menos que isso? Para descobrir, temos que analisar as figurinhas!

A Distribuição das Figurinhas


Para descobrir a distribuição das figurinhas, eu comprei ##t=106## pacotinhos (ou seja, ##n=106\times 4=424## figurinhas), e anotei o conteúdo de cada um deles (os dados originais estão no google docs). A primeira coisa que fiz foi traçar o histograma para ver o jeitão dos dados:

Pelo histograma dá para ver que faltam nove figurinhas para completar meu álbum (são os nove buracos no histograma). Por outro lado tem três figurinhas que repetiram 6 vezes, o que indica que possivelmente esse álbum não é honesto, já que algumas figurinhas são mais fáceis que outras.

Ou... talvez não? A intuição humana é notoriamente falha em intuir probabilidades, de repente essa impressão de que algumas figurinhas são mais fáceis é falsa. O jeito correto de decidir isso é fazendo um teste de hipótese. O teste em si é simples, mas eu vou fazer em câmera lenta caso alguém nunca tenha visto.

Para começar, precisamos de uma hipótese nula ##H_0## e um nível de significância ##\alpha##.

Nossa hipótese nula ##H_0## é que as figurinhas tem distribuição uniforme e são independentes. O teste vai nos falar se esse hipótese precisa ser rejeitada. O nível de significância nós podemos escolher, mas escolher o ##\alpha## é uma arte: você precisa levar em conta o número de amostras, e as taxas aceitáveis de falsos positivos e falsos negativos. Usualmente, ##\alpha##=5% funciona bem.

Qual o significado desses parâmetros? Fica simples de entender com um exemplo: suponha que todas as 424 figurinhas que eu comprei vieram repetidas. Em uma distribuição é uniforme, isso não é impossível, mas é extremamente improvável: a chance é de 1 em ##180^{424}\sim 1.7\times 10^{956}##. Esse valor é maior que o número de átomos no universo (muito maior!), então a chance de ter acontecido ao acaso é quase nula.

Se eu tirasse 424 figurinhas repetidas, eu teria bastante convicção de que as figurinhas não são uniformes. Quando fazemos um teste de hipótese, é essa ideia que estamos tentando capturar: os dados que coletamos são tão raros que não poderiam ter acontecido na prática? A chance de terem acontecido ao acaso é menor que ##\alpha##=5%? Se for, então talvez a hipótese ##H_0## seja falsa.

Existem vários testes de hipótese possíveis. Como queremos checar se nossas amostras seguem uma distribuição dada (uniforme), então o teste mais natural a se fazer é o teste do ##\chi^2## (pronuncia-se quiquadrado), que serve justamente para verificar se duas distribuições são iguais.

Você começa o teste do ##\chi^2## separando suas amostras em categorias. No nosso caso, o mais natural seria fazer cada figurinha ser uma categoria, mas isso não vai funcionar. O motivo é que o teste do ##\chi^2##, estritamente falando, só funciona com infinitas amostras. Como é impossível coletar infinitas amostras, nós assumimos que o resultado não muda muito de "infinitas amostras" para "muitas amostras".

Mas quantas amostras são "muitas amostras"? Em geral, a regra usada é que ##np>5##, onde ##n## é o número de amostras e ##p## é a menor probabilidade de uma amostra estar em uma categoria. Vamos fazer as contas: se cada figurinha é uma amostra, e as 180 figurinhas são uniformes, então ##p=1/180##. Como tenho ##n=424## figurinhas, então ##np=424/180\sim 2.35##. Oops, 2.35 é menor que 5, então 180 categorias é muito.

Vamos reduzir então. Eu vou escolher ##c=5## categorias, porque 180 divide certinho por 5. Cada categoria vai ter ##180/5=36## figurinhas, logo ##p=1/36## e ##np=424/36\sim 11.77##, agora sim! O histograma com cinco categorias é assim:


Agora nós vamos tabelar nossas amostras. A linha Observado mostra quantas amostras nós coletamos em cada categoria, a linha Esperado mostra quantas amostras uma distribuição ideal teria (no caso uniforme, todas são iguais a ##n/c=424/5=84.8##):

Categoria 1Categoria 2Categoria 3Categoria 4Categoria 5
Observado7590809584
Esperado84.884.884.884.884.8

Podemos calcular o ##\chi^2## agora. A intuição do teste é simples, vamos calcular qual o quadrado da diferença entre o observado e o esperado, e normalizar. Se esse número for muito grande, então tem algo errado. Calculando pela fórmula, temos:

$$\chi^2=\sum\frac{(Y_i-E_i)^2}{E_i}\sim 2.957$$
Esse valor, 2.957 é grande ou é pequeno? Para decidir isso precisamos primeiro calcular quantos graus de liberdade tem esse modelo. O número de graus de liberdade é o número de parâmetros que estamos estimando. Como nossa tabela tem 5 categorias, então o número de graus de liberdade é 4.

Opa! Se são 5 categorias, não tinha que ser 5 graus de liberdade? É verdade que estamos tentando estimar ##Y_1##, ##Y_2##, ##Y_3##, ##Y_4## e ##Y_5##. Mas nós sabemos que a soma de todos eles é ##n##, então se tivermos os quatro primeiros, o quinto é só calcular!

Em tempos idos esse era o momento em que você precisava sacar do bolso uma tabela de ##\chi^2## (como curiosidade, coloquei online a tabela que eu usava na Poli). Mas hoje em dia você pode ir direto no Mathematica para descobrir o valor de ##\chi^2## com ##\alpha##=5% e quatro graus de liberdade:
InverseCDF[ChiSquareDistribution[4], 1 - 0.05]
>> 9.48773

Como 2.957 é menor que 9.488, então não podemos rejeitar a hipótese de que as figurinhas são uniformes e independentes, e concluímos que sim, as figurinhas provavelmente são honestas.

Alguns de vocês irão chiar, eu sei. "Mas veja, algumas figurinhas saíram repetidas seis vezes, outras não apareceram nenhuma vez! Não tem como isso ser uniforme, claramente tem figurinhas que são mais fáceis!".

Bem, se você ainda estiver na dúvida, você pode comprar mais pacotinhos para deixar o teste mais robusto. Ao invés disso, eu vou usar só as figurinhas que eu tenho para ir direto na raiz da dúvida: afinal de contas, é normal tirar muitas figurinhas repetidas?

A Análise das Repetidas


A intuição nos diz que, em uma distribuição uniforme, todas as figurinhas tem que aparecer mais ou menos com a mesma frequência. Só que esse é mais um caso onde a intuição é falha! Quem lê o blog há mais tempo deve lembrar de um post sobre o paradoxo do aniversário, que diz exatamente isso: mesmo em distribuições uniformes, a chance de colisões é muito alta.

Podemos usar esse fato para tentar um novo teste de hipótese. Dessa vez, vamos contar quantas repetidas nós temos para cada figurinha, e separá-las por número de repetições. O histograma da quantidade de figurinhas repetidas fica assim:


O histograma foi traçado da quantidade de figurinhas repetidas:

Quantidade de repetições0123456
Figurinhas94150482453

O álbum tem 180 figurinhas, então a soma da segunda linha vale 180. Do total de figurinhas, nove delas não apareceram, três apareceram 6 vezes, e a maioria apareceu entre 1 e 3 vezes.

O formato do histograma parece familiar né? Se nós soubéssemos qual a é a fórmula exata para o número de figurinhas repetidas em uma distribuição uniforme, poderíamos usar o teste do ##\chi^2## para ver se a hipótese ainda se mantém. Só que agora não tem jeito, para descobrir a fórmula precisamos fazer as contas. Pule a caixa azul se você tem medo de combinatória:

O truque para usar combinatória analítica numerada (labelled) é sempre tentar achar um isomorfismo entre o seu problema e uma lista de inteiros. Para o nosso caso, imagine o seguinte experimento: temos ##n## bolinhas numeradas de ##1## a ##n##, e temos também ##m## caixinhas, cada caixa representando uma figurinha. A bolinha representa a ordem em que eu comprei a figurinha. Se eu atirar as bolinhas em caixas aleatórias, tenho um problema isomórfico a comprar figurinhas aleatórias.

Por exemplo: se a bolinha 5 cair na caixa 32, então significa que a quinta figurinha que eu comprei foi a de número 32. Se no fim do experimento a caixa 32 tiver seis bolinhas, então eu tenho a figurinha 32 repetida seis vezes.

No nosso experimento, as caixinhas são diferentes umas das outras, então vamos usar o operador ##\text{SEQ}_m(z)## para modelar as caixinhas. A ordem das bolinhas dentro da caixinha não importa (tanto faz se dissermos que a caixa 32 tem as bolinhas 3, 5, 9, 10 e 15; ou se tiver as bolinhas 9, 5, 3, 15, e 10, é a mesma coisa). Então dentro das caixinhas, usaremos o operador ##\text{SET}(z)##. O experimento completo é ##\text{SEQ}_m(\text{SET}(z))##.

Digamos que estamos interessados em saber quantas caixinhas vão terminar com ##s## bolinhas. Um jeito de medir isso é colocar um marcador nas caixinhas com ##s## elementos. Para isso, nós tiramos as das caixinhas os conjuntos de ##s## elementos, e colocamos eles de volta com o marcador, digamos que o marcador seja ##u##. Então o experimento para medir ##s## repetidas é ##\text{SEQ}_m(\text{SET}(z)-\text{SET}_s(z)+u \text{SET}_s(z))##.

Os experimentos vão mapear para funções geradoras exponenciais de duas variáveis. Das tabelas de combinatória analítica:
$$\begin{align*}
\text{SEQ}_m(z) &= z^m \\
\text{SET}(z) &= e^z \\
\text{SET}_s(z) &= \frac{z^s}{s!}\\
\text{SEQ}_m(\text{SET}(z)-\text{SET}_s(z)+u \text{SET}_s(z)) &=
\left(e^z-\frac{z^s}{s!}+u \frac{z^s}{s!}\right)^m \\
&= \left(e^z+ (u-1) \frac{z^s}{s!}\right)^m
\end{align*}$$
Se abrirmos essa função geradora em série, então o coeficiente do termo ##z^n u^k## indica quantas combinações possíveis de ##n## bolinhas vão ter ##k## caixinhas com ##s## bolinhas repetidas (como a função geradora é exponencial, você precisa multiplicar o coeficiente por ##n!## para ter o valor real). Se chamarmos esse coeficiente de ##R_s(n, k)##, então temos:
$$\sum_{n}\sum_{k}R_s(n,k) z^n u^k = \left(e^z+ (u-1) \frac{z^s}{s!}\right)^m $$
Vamos listar com cuidado o que temos. ##n! R_s(n,k)## é número de combinações possíveis com ##k## caixinhas contendo exatamente ##s## bolinhas. Logo, de todas as combinações possíveis de ##n## bolinhas, existem ##\sum k n! R_s(n,k)## caixinhas com ##s## bolinhas.

Por outro lado, o número total de combinações de ##n## bolinhas em ##m## caixinhas é ##m^n##. Como cada combinação dessas tem ##m## caixinhas, então o número total de caixinhas que podem ser preenchidas em todas as combinações é ##m^{n+1}##.

Juntando as duas contagens, a probabilidade de um caixinha tem ##s## bolinhas em ##n## arremessos é:
$$p_s(n)=\sum_{k>0}\frac{k \;n!}{m^{n+1}}R_s(n, k)$$
Mas olha só! Nós podemos conseguir essa somatória derivando a função geradora em relação a ##u##, e depois setando ##u=1##.
$$\left.\frac{\partial}{\partial u}\sum_{n}\sum_{k}R_s(n,k) z^n u^k \right|_{u=1}=\sum_n \sum_{k}k R_s(n,k) z^n$$
Agora é só fazer as contas:
$$\begin{align*}p_s(n)
&=\sum_k \frac{n!}{m^{n+1}}k R_s(n, k) \\
&=\frac{n!}{m^{n+1}} [z^n] \left.\frac{\partial}{\partial u}\sum_{n}\sum_{k}R_s(n,k) z^n u^k \right|_{u=1}\\
&=\frac{n!}{m^{n+1}} [z^n] \left.\frac{\partial}{\partial u}\left(e^z+ (u-1) \frac{z^s}{s!}\right)^m  \right|_{u=1}\\
&=\frac{n!}{m^{n+1}} [z^n] \left. m\left(e^z+ (u-1) \frac{z^s}{s!}\right)^{m-1} \frac{z^s}{s!} \right|_{u=1}\\
&=\frac{n!}{m^{n+1}} [z^n]  m\left(e^z\right)^{m-1} \frac{z^s}{s!} \\
&=\frac{n!}{s!\;m^n} [z^n]z^s e^{z(m-1)}\\
\end{align*}$$
Abrindo a exponencial em série:
$$\begin{align*}p_s(n)
&=\frac{n!}{s!\;m^{n}} [z^n]z^s\sum_{k\ge 0}\frac{(z(m-1))^k}{k!} \\ &=\frac{n!}{s!\;m^{n}} [z^n]\sum_{k\ge 0}\frac{z^{k+s}(m-1)^k}{k!} \\ &=\frac{n!}{s!\;m^{n}} [z^n]\sum_{p\ge s}\frac{z^p(m-1)^{p-s}}{(p-s)!} \\ &=\frac{n!}{s!\;m^{n}}\; \frac{(m-1)^{n-s}}{(n-s)!} \\ &=\frac{n!}{s!(n-s)!}\; \frac{(m-1)^{n-s}}{m^{n}} \\ &={{n}\choose{s}}\frac{(m-1)^{n-s}}{m^{n}} \\ \end{align*}$$
Essa é a solução exata, em forma fechada, que queríamos! Mas tem uma massagem que podemos fazer que vai ajudar a entender por que essa fórmula é assim. Em geral ##n## é bem maior que ##m##, então vamos reescrever ##m## de modo que ##n=\lambda m##, e calcular o limite quando ##n## é muito grande: $$\begin{align*}\lim_{n\to\infty}p_s(n) &=\lim_{n\to\infty}\frac{n!}{s!(n-s)!}\; \frac{(n/\lambda-1)^{n-s}}{(n/\lambda)^{n}} \\ &=\lim_{n\to\infty}\frac{n!}{s!(n-s)!}\; \left(\frac{n-\lambda}{\lambda}\right)^{n-s} \left(\frac{\lambda}{n}\right)^{n}\\ &=\lim_{n\to\infty}\frac{n!}{s!(n-s)!}\; \lambda^{s} \left(\frac{n-\lambda}{n}\right)^{n-s}\left(\frac{1}{n}\right)^{s} \\ &=\lim_{n\to\infty}\frac{n!\;\lambda^{s}}{s!(n-s)!\;n^{s}}\; \left(1-\frac{\lambda}{n}\right)^{n-s}\\ &=\lim_{n\to\infty}\left(\frac{\lambda^{s}}{s!}\right) \left(\frac{n!}{(n-s)!\;n^s}\right) \left(1-\frac{\lambda}{n}\right)^{n} \left(1-\frac{\lambda}{n}\right)^{-s} \\ \end{align*}$$ O limite do produto é o produto dos limites, desde que todos existam e sejam finitos, então vamos conferir um por um. O primeiro termo é constante, easy. No segundo, o numerador é um polinômio de grau ##n##, o denominador é um polinômio de grau ##n##, então a fração é constante e vale 1, já que os coeficientes iniciais ambos valem 1. O terceiro é a definição de ##e^{-\lambda}##. No quarto a segunda fração vai para zero, então o total tende a 1. Combinando tudo: $$\lim_{n\to\infty}p_s(n)= e^{-\lambda}\frac{\lambda^s}{s!}$$ Ahá, é uma distribuição de Poisson!

Não foi à toa que aquele histograma das figurinhas repetidas parecia familiar, elas seguem uma distribuição de Poisson! Em retrospecto até faz sentido. Se a distribuição das figurinhas é uniforme, então em um longo período de tempo, cada figurinha deve ter uma taxa de aparição constante (igual a ##\lambda=n/m##). Se a taxa é constante ao longo do tempo, então a quantidade de figurinhas em um pequeno intervalo tem que seguir a distribuição de Poisson.

Agora falta pouco para conseguir usar o teste do ##\chi^2## nas repetidas. Precisamos primeiro separar as amostras em categorias. O natural seria usar uma categoria por quantidade de repetidas, mas em teoria, um cara muito azarado poderia ter 424 figurinhas repetidas, e nesse caso o produto ##np## seria muito, muito menor que 5. Vamos então aglutinar todas as categorias de 6 repetidas em diante em um grupo só. A tabela fica assim:

Quantidade de repetições0123456 ou mais
Figurinhas, medido94150482453
Figurinhas, teórico17.140.247.437.221.910.35.9

Agora podemos calcular o ##\chi^2## medido e o teórico, lembrando que temos 6 graus de liberdade e ##\alpha## valendo 5%: $$\chi^2=\sum\frac{(Y_i-E_i)^2}{E_i}\sim 11.5275$$
InverseCDF[ChiSquareDistribution[6], 1 - 0.05]
>> 12.59
Conclusão: como 11.53 é menor que 12.59, então novamente não podemos rejeitar a hipótese de que as figurinhas são uniformes e independentes, e concluímos que sim, as figurinhas provavelmente são honestas. É verdade que tem muitas repetidas, mas o número de repetidas está dentro do esperado.

Quanto gastar para completar?


No pŕoximo post, vamos finalmente calcular qual o minimo de dinheiro que você precisa gastar para completar o álbum. Stay tuned!

segunda-feira, 5 de setembro de 2016

A Matemática do Pokémon GO

Você já capturou todos os pokémons do Pokémon GO? Eu ainda não consegui, e também não conheço ninguém que tenha conseguido. O motivo é simples: para cada Pikachu que você encontrar, vai passar por um centena de Zubats antes! Fica difícil pegar todos quando alguns são tão raros!

pokemon_zubat

A pergunta natural agora é: se eu quiser insistir e continuar jogando até pegar todos, quanto tempo vai levar? Ou, de outro modo, qual é o número esperado de capturas até completar a coleção? Esse problema é difícil de resolver no caso geral, mas fica mais simples usando Combinatória Analítica!

Combinatória Analítica Numerada


No último post do blog eu mostrei como funciona a Combinatória Analítica, mostrando uma aplicação de cálculo com sequências de cara ou coroa. Mas essa foi só metade da história! Aquele tipo de cálculo só funciona quando os objetos contados são indistinguíveis (um resultado cara não é diferente de outro resultado cara).

Se você quiser lidar com objetos que são todos diferentes entre si, então você precisa de outra abordagem. Por exemplo, quantas permutações de ##n## elementos existem? Nesse caso, os elementos precisam ser todos diferentes. Quando ##n=3##, existem ##6## permutações:
$$ \enclose{circle}{1} \enclose{circle}{2} \enclose{circle}{3} \longrightarrow \begin{cases} \enclose{circle}{1} \enclose{circle}{2} \enclose{circle}{3} \\ \enclose{circle}{1} \enclose{circle}{3} \enclose{circle}{2} \\ \enclose{circle}{2} \enclose{circle}{1} \enclose{circle}{3} \\ \enclose{circle}{2} \enclose{circle}{3} \enclose{circle}{1} \\ \enclose{circle}{3} \enclose{circle}{1} \enclose{circle}{2} \\ \enclose{circle}{3} \enclose{circle}{2} \enclose{circle}{1} \\ \end{cases}$$
Para esse tipo de cálculo nós precisamos da Combinatória Analítica Numerada. Existem duas diferenças principais entre as versões numeradas e não-numeradas da teoria. A primeira que é, ao invés de usar funções geradoras normais, na versão numerada nós usamos funções geradoras exponenciais (EGF), que são definidas como abaixo: 
$$EGF(F[n]) = \sum_{n\ge 0}{F[n]\frac{z^n}{n!}}$$
Quando a ordem dos objetos é importante, os resultados tendem a crescer muito rapidamente, então dividir por ##n!## termo a termo ajuda as funções a ficarem bem comportadas.

A outra diferença é em um dos operadores básicos. A adição continua representando a união de conjuntos, mas a multiplicação, que antes era concatenação, fica um pouco diferente. No caso numerado, você não pode simplesmente concatenar, porque os números precisam ser sempre de ##1## a ##n##, e a concatenação simples te deixaria com duplicatas:
$$ \enclose{circle}{1} \enclose{circle}{2} \star \enclose{circle}{2} \enclose{circle}{1} =\; ?$$
A solução é renomear as bolinhas após concatenar, mas você precisa renomear de todas as maneiras possíveis que deixem a ordem dentro dos subconjuntos originais intacta:
$$ \enclose{circle}{1} \enclose{circle}{2} \star \enclose{circle}{2} \enclose{circle}{1} = \begin{cases} \enclose{circle}{1} \enclose{circle}{2} \; \enclose{circle}{4} \enclose{circle}{3} \\ \enclose{circle}{1} \enclose{circle}{3} \; \enclose{circle}{4} \enclose{circle}{2} \\ \enclose{circle}{1} \enclose{circle}{4} \; \enclose{circle}{3} \enclose{circle}{2} \\ \enclose{circle}{2} \enclose{circle}{3} \; \enclose{circle}{4} \enclose{circle}{1} \\ \enclose{circle}{2} \enclose{circle}{4} \; \enclose{circle}{3} \enclose{circle}{1} \\ \enclose{circle}{3} \enclose{circle}{4} \; \enclose{circle}{2} \enclose{circle}{1} \\ \end{cases}$$
Com isso já conseguimos definir as permutações. Elas são formadas do conjunto vazio, ou de um elemento concatenado e renomeado com uma permutação:
$$P=\varepsilon + z\star P$$
Traduzimos direto para a EGF:
$$P=1 +z P\implies P=\frac{1}{1-z}$$
E extraímos os coeficientes por comparação com a definição de EGF:
$$\begin{align*} \sum_{n\ge 0}P[n]\frac{z^n}{n!} &= \frac{1}{1-z} \\ \sum_{n\ge 0}\frac{P[n]}{n!}z^n &= \sum_{n\ge 0}1\times z^n \\ \frac{P[n]}{n!} &= 1 \\ P[n] &= n! \end{align*}$$
A quantidade de permutações de tamanho ##n## é ##n!##, então o método funciona!

O operador SEQ


Também podemos definir o operador ##SEQ##, que enumera todas as sequências possíveis de um conjunto:
$$SEQ(B[z]) = \frac{1}{1-B[z]}$$
Note que as permutações também podem ser vistas como as sequências renomeadas de zero ou mais átomos, então a definição acima dá o mesmo resultado.

Outra definição útil é ##SEQ_k##, que enumera apenas as sequências com ##k## elementos:
$$SEQ_k(B[z]) = (B[z])^k$$

O operador SET


Todas as funções acima são para construções onde a ordem é importante. Quando a ordem não importa, você sempre pode voltar para a Combinatória Analítica Não-Numerada. Mas e se o seu problema envolver os dois tipos ao mesmo tempo?

Nesse caso, a solução é o operador ##SET##. Ele é um operador numerado que define um conjunto onde a ordem não importa. A definição é a abaixo:
$$SET(B[z]) = e^{B[z]}$$
Por exemplo, quantas sequências possíveis de números de ##1## a ##n## existem, se a ordem não for importante? Isso é dado por:
$$SET(z) = e^z$$
Agora basta abrir a exponencial em série de Taylor e comparar com a definição de EGF:
$$\begin{align*} \sum_{n\ge 0}F[n]\frac{z^n}{n!} &= e^z \\ \sum_{n\ge 0}\frac{F[n]}{n!}z^n &= \sum_{n\ge 0}\frac{1}{n!}z^n \\ \frac{F[n]}{n!} &= \frac{1}{n!} \\ F[n] &= 1 \end{align*}$$
Ou seja, quando a ordem não é importante, só existe uma única sequência de ##1## a ##n## que usa todos os números de ##1## a ##n##. Faz sentido né.

(É por causa desse caso que a EGF chama função geradora exponencial. A EGF de uma sequência de ##1##s é a função exponencial.)

Pokémons uniformes


Com as ferramentas em mãos, já conseguimos modelar o problema do Pokémon GO! Mas, antes de começar o problema real, vale a pena estudar um caso mais simples. Vamos supor que os pokémons são uniformes e resolver esse caso mais fácil primeiro, depois a gente generaliza para o caso onde alguns pokémons são mais raros que outros.

O truque para usar Combinatória Analítica Numerada é achar uma maneira de enxergar seu problema como uma sequência de ##1## de ##n##. No nosso caso podemos fazer da seguinte maneira: para cada pokémon nós temos uma caixinha, e cada vez que achamos um pokémon na rua, nós colocamos dentro da caixinha uma bola com um número indicando a ordem em que ele foi encontrado.

Por exemplo, suponha que andávamos pela rua e fizemos 6 capturas, nessa ordem: Bulbassauro, Squirtle, Bulbassauro de novo, Charmander, Squirtle de novo, Bulbassauro de novo. As caixinhas ficariam assim:

bulbassauro##\enclose{box}{\enclose{circle}{1} \,\enclose{circle}{3}\, \enclose{circle}{6}}##
squirtle##\enclose{box}{\enclose{circle}{2}\, \enclose{circle}{5} }##
charmander##\enclose{box}{\enclose{circle}{4}} ##

Se nenhuma caixinha está vazia, então em algum ponto nós pegamos todos os pokémons (nesse caso, foi na quarta captura). Então, para ##n## capturas, se calcularmos a quantidade de configurações em que nenhuma caixinha está vazia, dividida pela quantidade de configurações totais, nós vamos ter a probabilidade de que a coleção foi completada com ##n## capturas ou menos.

Vamos lá. Se nós temos ##k## pokémons e ##n## capturas, então a quantidade de configurações totais é ##k^n## (é só atribuir uma pokémon a cada captura).

Para a quantidade de configurações em que nenhuma caixinha está vazia, precisamos modelar o problema com os operadores que descrevemos antes. Os pokémons são diferentes entre si, e a ordem importa, então vamos usar o operador ##SEQ##. A ordem das bolinhas dentro da caixinha não importa (afinal, ##\enclose{circle}{1} \enclose{circle}{3} \enclose{circle}{6}## é a mesma coisa que ##\enclose{circle}{3} \enclose{circle}{6} \enclose{circle}{1}## nesse contexto: a ordem das capturas é determinada pelos números nas bolinhas, não pela posição delas). Então vamos usar o operador ##SET##. E não podemos esquecer de que nenhuma caixinha pode ser vazia. Então a construção combinatória que descreve o problema é:
$$C=SEQ_k(SET(z)-\varepsilon)$$
..ou seja, uma sequência de ##k## conjuntos não-vazios. A função geradora é, portanto: 
$$C[z]=(e^z-1)^k$$
Note que, se nós já temos a função geradora, então a parte combinatória do problema já acabou. Daqui em diante é só abrir a caixa azul e fazer as contas!

Nós queremos calcular o valor esperado do número de capturas até completar a coleção. Pela definição, isso é: $$E[X]=\sum_{n\ge 0}n P(X=n)$$ Onde ##P(X=n)## é a probabilidade de que nós completamos a coleção na captura ##n##. Mas nós não temos isso! O que temos é ##P(x\le n)##, a probabilidade de completamos a coleção com pelo menos ##n## capturas. Felizmente não tem problema, porque dá para converter um no outro: $$\begin{align*} E[X] &= \sum_{n\ge 0} n P(X=n) \\ &= \sum_{n\ge 0} \sum_{i\lt n} P(X=n) \\ &= \sum_{n\ge 0} \sum_{i\ge 0} P(X=n) [i\lt n] \\ &= \sum_{i\ge 0} \sum_{n\ge 0} P(X=n) [n\gt i] \\ &= \sum_{i\ge 0} \sum_{n\gt i} P(X=n) \\ &= \sum_{i\ge 0} P(X\gt i) \\ &= \sum_{i\ge 0} 1-P(X\le i) \\ \end{align*}$$ Esse truque funciona sempre que a distribuição é discreta, porque podemos abrir a probabilidade: $$P(X>n)=P(X=n+1)+P(X=n+2)+\cdots$$ Vamos substituir agora. A nossa probabilidade para ##X\le n## é a divisão da quantidade de configurações com caixinhas não-vazias (que é o coeficiente ##z^n## da função geradora, vezes ##n!## porque é uma função geradora exponencial) pela quantidade total (que é ##k^n##): $$P(X\le n)=\frac{n![z^n](e^z-1)^k}{k^n}$$ Logo a expectativa é: $$\begin{align*} E[X] &=\sum_{n\ge 0} 1-P(X\gt n) \\ &=\sum_{n\ge 0} 1-\frac{n![z^n](e^z-1)^k}{k^n} \\ &=\sum_{n\ge 0} 1-n![z^n](e^{z/k}-1)^k \end{align*}$$ Esse último passo usou uma propriedade de escala das funções geradoras. Se uma sequência ##f(n)## tem função geradora ##F(z)##, então: $$\begin{align*} \sum_{n\ge 0}f(n)z^n &= F(z) \\ \sum_{n\ge 0}k^n f(n)z^n &= \sum_{n\ge 0}f(n)(kz)^n = F(kz) \end{align*}$$ Voltando: nós vimos que a função geradora exponencial da sequência constante de ##1##s é ##e^z##, então dá pra jogar para dentro: $$\sum_{n\ge 0} 1-n![z^n](e^{z/k}-1)^k = \sum_{n\ge 0} n![z^n]\left(e^z -\left(e^{z/k}-1\right)^k\right)$$ Para prosseguir agora tem um truque muito bom, que funciona assim: se você tem uma função geradora ##F(z)## bem comportada, então vale que: $$\sum_{n\ge 0}n![z^n]F(z)=\int_0^\infty e^{-t}F(t)\,dt$$ Parece meio mágico, mas é só consequência da forma integral do fatorial: $$n!=\int_0^\infty t^n e^{-t}\,dt$$ Se nós chamarmos os coeficientes da função geradora de ##f_n##, então: $$\begin{align*} \sum_{n\ge 0}n![z^n]F(z) &=\sum_{n\ge 0}n! f_n = \sum_{n\ge 0} f_n n! \\ &=\sum_{n\ge 0}f_n\left(\int_0^\infty t^ne^{-t}\,dt\right) \\ &=\int_0^\infty \left(\sum_{n\ge 0}f_n t^n e^{-t}\right)\,dt\\ &=\int_0^\infty e^{-t}\left(\sum_{n\ge 0}f_n t^n \right)\,dt\\ &=\int_0^\infty e^{-t}F(t)\,dt\\ \end{align*}$$ Aplicando ao nosso caso: $$\begin{align*} \sum_{n\ge 0} n![z^n]\left(e^z -\left(e^{z/k}-1\right)^k\right) &= \int_0^\infty e^{-t} \left(e^t -\left(e^{t/k}-1\right)^k\right) \,dt \\ &= \int_0^\infty \left(1 -e^{-t}\left(e^{t/k}-1\right)^k\right) \,dt \\ &= \int_0^\infty \left(1 -\left(e^{-t/k}\right)^k\left(e^{t/k}-1\right)^k\right) \,dt \\ &= \int_0^\infty \left(1 -\left(1-e^{-t/k}\right)^k\right) \,dt \\ \end{align*}$$ Nesse ponto já daria para jogar num solver numérico e calcular a solução. Mas essa integral tem uma forma fechada! Para calcular, basta fazer a mudança de variável abaixo: $$\begin{align*} v &= 1-e^{-t/k} \\ 1-v &= e^{-t/k} \\ e^{t/k} &= \frac{1}{1-v} \end{align*}$$ Agora deriva para achar quanto vale ##dt##: $$\begin{align*} \frac{dv}{dt} &= \frac{e^{-t/k}}{k} \\ dt &= k \, e^{t/k}\, dv \\ dt &= \frac{k}{1-v} dv \end{align*}$$ Arruma os intervalos: $$\begin{align*} t=0 &\implies v=1-e^0=1-1=0 \\ t=\infty &\implies v=1-e^{-\infty}=1-0=1 \end{align*}$$ E por fim soca tudo na integral: $$\begin{align*} \int_0^\infty \left(1 -\left(1-e^{-t/k}\right)^k\right) \,dt &= \int_0^1 \left(1 -v^k\right) \frac{k}{1-v} \,dv \\ &= k\int_0^1 \frac{1-v^k}{1-v} \,dv \\ &= k\int_0^1 (1+v+v^2+\cdots+v^{k-1}) \,dv \\ &= k \left[v+\frac{v^2}{2}+\frac{v^3}{3}+\cdots+\frac{v^k}{k} \right]_0^1 \\ &= k \left[1+\frac{1}{2}+\frac{1}{3}+\cdots+\frac{1}{k} \right] \\ &= k H[k] \\ \end{align*}$$ Esse é o resultado final. O número esperado de capturas até você completar sua coleção de pokémons é ##k## vezes o harmônico de ##k##!

Qual seria o número esperado de capturas então? O jogo tem 151 pokémons, mas desses tem seis que ainda não podem ser pegos: Ditto, Articuno, Zapdos, Moltres, Mewtwo e Mew. Para pegar os 145 restantes, você precisaria de ##145\times H(145)\sim 806## capturas.

Mas isso só vale se os pokémons fossem uniformes! Como alguns são mais raros que outros, precisamos melhorar nossa estimativa.

Pokémons raros


A grande vantagem da Combinatória Analítica é que ela generaliza muito fácil. O problema de completar uma coleção com elementos não-uniformes é bem difícil com técnicas tradicionais, mas aqui fica simples: basta trocar os átomos ##Z_i## por ##p_iZ_i##, onde ##p_i## é a probabilidade de ocorrência dele.

Por que isso funciona? Suponha que você tem três pokémons, de probabilidades ##p_1=1/6##, ##p_2=1/3## e ##p_3=1/2##. Ao invés de fazer ##SEQ_3(SET(z)-\varepsilon)##, como no caso uniforme, você pode fazer:
$$(SET(Z_1)-\varepsilon)\star (SET(Z_2+Z_2)-\varepsilon)\star (SET(Z_3+Z_3+Z_3)-\varepsilon)$$
Isso vai dar um resultado seis vezes maior que o esperado, mas na proporção correta. Se você jogar esse valor ##1/6## para dentro, a conta fica:
$$\begin{align*} &(SET(\frac{Z_1}{6})-\varepsilon)\star (SET(\frac{2Z_2}{6})-\varepsilon)\star (SET(\frac{3Z_3}{6})-\varepsilon)=\\ &(SET(p_1Z_1)-\varepsilon)\star (SET(p_2Z_2)-\varepsilon)\star (SET(p_3Z_3)-\varepsilon) \end{align*}$$
Exatamente como mencionamos! Isso mostra que o resultado é válido para qualquer probabilidade racional, e você pode fazer um sanduíche se quiser reais também.

Daqui em diante as contas são as mesmas do caso uniforme, mudando somente o interior da integral: 
$$\begin{align*} E[X] &=\int_0^\infty\left( 1-\left(1-e^{-p_1t}\right)\left(1-e^{-p_2t}\right)\cdots\left(1-e^{-p_kt}\right) \right)\,dt\\ &=\int_0^\infty\left( 1- \prod_{1\le i\le k}\left(1-e^{-p_i t}\right) \right)\,dt\\ \end{align*}$$
Infelizmente essa integral não tem forma fechada, precisa ser calculada numericamente mesmo. Usando uma tabela com as probabilidades estimadas de cada pokémon, eu calculei a integral, e o resultado é que você precisaria de 51568 capturas para completar sua coleção!

Como de costume, eu fiz uma simulação para conferir os valores. Os scripts estão no github, e os resultados são:

Distribuição uniforme, valor teórico 805.8, computacional 805.1
Distribuição não-uniforme, valor teórico 51568, computacional 50931

Os resultados bateram como previsto! Vale a pena também estimar o tempo: se você pegar um pokémon a cada 3 minutos, então para fazer 51568 capturas você precisaria jogar por 107 dias seguidos, sem parar para dormir!

Na prática você pode completar em menos tempo que isso usando outros recursos do jogo, como evoluções e ovos, mas não tem muito segredo, para ser um mestre pokémon você precisa se dedicar bastante :)

PS: Como quase não tem literatura em português sobre o tema, as traduções ainda não são padrão. Eu traduzi Labelled Analytic Combinatorics como Combinatória Analítica Numerada, porque não achei nenhuma outra tradução que não fosse feia demais.

sábado, 30 de julho de 2016

Totorial de Combinatória Analítica

Se você jogar Cara ou Coroa várias vezes em sequência, qual a chance de nunca aparecer duas caras seguidas?

Esse é um problema que você resolve com ferramentas da Combinatória. Se você estudou para o vestibular, certamente deve lembrar que esse tipo de problema usualmente se resolvia com permutações, arranjos ou combinações.

Infelizmente, isso é verdade no vestibular, mas não na vida real. Os problemas avançados de Combinatória resistem a essas ferramentas básicas. Em geral, você vai precisar de recorrências complexas e somatórias complicadas; e até recentemente não tinha nenhum método que fosse simples e geral, para usar em qualquer situação.

A boa notícia é que isso mudou! Nos últimos anos, surgiu um novo ramo da matemática, chamado de Combinatória Analítica, que permite calcular esses problemas com facilidade. Como ainda não existe nenhum texto introdutório em português, resolvi escrever um totorial para mostrar como as coisas ficam bem simples com esse método!
totoro

Um exemplo simples


Como funciona esse método novo? A analogia mais simples é com desenho geométrico. Com as técnicas tradicionais, você precisa de régua, compasso e uma boa dose de insight para resolver os problemas. Ao invés disso, você pode usar geometria analítica: o problema é transformado em uma série de equações, aí você não precisa pensar, só resolver (e em muitos casos nem precisa resolver, é só jogar num sistema de computação automática como o Wolfram Alpha que ele resolve para você).

Com combinatória analítica é parecido: você vai transformar a descrição do seu problema de contagem em uma equação, e aí existem técnicas padrão para resolvê-la. Para mostrar como funciona, vamos pegar um problema simples: Quantas strings binárias de tamanho ##n## existem?

Para ##n=3##, por exemplo, existem oito strings: ##000##, ##001##, ##010##, ##011##, ##100##, ##101##, ##110## e ##111##.

A primeira coisa a ser feita é caracterizar o conjunto de todas as strings. Podemos recursivamente construir as strings binárias válidas com três casos:

Caso 1: A string vazia (##\varepsilon##), ou...

Caso 2: Uma string que começa com ##0##, seguida de uma string válida, ou...

Caso 3: Uma string que começa com ##1##, seguida de uma string válida.

Pronto, agora podemos agora escrever a equação que descreve o conjunto ##B## das strings binárias. Em combinatória analítica, a regra é que a descrição "ou" vira adição, e a operação "seguida de" vira multiplicação. Vou também substituir o caractere ##0## por ##Z_0## e ##1## por ##Z_1## só para deixar claro que nesse contexto são átomos, e não números:
$$B=\varepsilon+Z_0\times B+Z_1\times B$$
Essa é a construção combinatória correspondente às strings binárias. Agora o pulo do gato: eu vou trocar ##\varepsilon## pelo número ##1##, cada átomo individual pela variável ##z##, e resolver a equação:
$$\begin{align*} B&=\varepsilon+Z_0\times B+Z_1\times B\\ B&=1+zB+zB\\ B&=\frac{1}{1-2z} \end{align*}$$
Chegamos então no elemento principal da combinatória analítica: a função geradora do conjunto das strings binárias. E por que ela é tão importante? É que essa função simples tem escondida dentro dela a solução do problema para qualquer ##n##! Basta expandir a função em série (nesse caso, basta reconhecer que essa é a fórmula da soma infinita de uma PG com razão ##2z##): 
$$B=\frac{1}{1-2z}=1+2z+4z^2+8z^3+16z^4+\cdots$$
Olha só, a solução aparece na série! O termo ##8z^3## significa que, para ##n=3##, a solução é ##8##, exatamente como tínhamos visto enumerando as possbilidades.

O operador SEQ


O método acima foi simples né? Mas dá para ficar mais fácil ainda! O pessoal que estuda combinatória analítica faz bastante tempo inventou uma série de operadores que deixam a solução mais fácil de escrever.

O primeiro deles é o operador ##SEQ##, que serve para definir sequências. Ao invés da definição recursiva com três casos, podemos simplesmente dizer que uma string binária é formada por uma sequência de zero ou mais átomos, onde os átomos possíveis são ##0## e ##1##. Daí a construção combinatória sai direto:
$$B=SEQ(Z_0+Z_1)$$
E como transformar isso em função geradora? Se você sabe que um conjunto ##A## tem função geradora ##A(z)##, então a função geradora da sequência de ##A## é: 
$$SEQ(A)=\frac{1}{1-A(z)}$$ 
No nosso caso, a função geradora das strings binárias é:
$$SEQ(Z_0+Z_1)=SEQ(z+z)=\frac{1}{1-2z}$$
Pronto, a função geradora saiu direto, sem precisar resolver nenhuma equação!

O teorema de transferência


Ainda resta o problema de abrir a função geradora em série. Nesse caso foi fácil porque conseguimos ver de cara que a função era a soma de uma PG, mas e em casos mais complicados?

Bem, resolver de maneira exata continua sendo um problema difícil até hoje. Mas na maioria das vezes você não precisa da solução exata, um assintótico para ##n## grande já é suficiente. E para isso a combinatória analítica possui uma série de teoremas de transferência, que dão a resposta assintótica da função geradora sem precisar abrir a série!

O teorema de transferência usado nesse caso é simples de usar, embora sua definição seja meio assustadora:

Se a sua função geradora for uma função racional ##f(z)/g(z)##, onde ##f(z)## e ##g(z)## são relativamente primas, ##g(0)\neq 0##, ##g(z)## tem um único pólo ##1/\beta## de módulo mínimo, e ##\nu## é a multiplicidade desse pólo, então um assintótico para a função é dado por: $$[z^n]\frac{f(z)}{g(z)}=\nu\frac{(-\beta)^\nu f(1/\beta)}{g^{(\nu)}(1/\beta)}\beta^n n^{\nu-1}$$

Eu juro que é simples, apesar de assustador. Para o nosso caso, ##f(z)=1##, ##g(z)=1-2z##, ##g'(z)=-2##, ##\beta=2##, ##f(1/2)=1##, ##g'(1/2)=-2## e ##\nu=1##. Substituindo os valores, temos:
$$[z^n]B(z)=1\times\frac{(-2)^1\times 1}{-2}\times(2)^n n^{1-1} =2^n $$
Portanto, existem ##2^n## strings binárias de tamanho ##n##, o que bate com nosso resultado de ##8## para ##n=3##.

Outro exemplo


E se eu quiser calcular quantas strings binárias existem, mas só as que não possuem dois zeros seguidos? Esse problema é complicado com combinatória tradicional, mas com combinatória analítica fica simples!

Primeiro, vamos caracterizar essas strings. Podemos descrevê-las como uma sequência de ##1## ou ##01##, seguida de ##0## ou vazio:
$$C=SEQ(Z_1+Z_0Z_1)\times(Z_0 + \varepsilon)$$
Traduzindo para função geradora:
$$C=SEQ(z+z^2)\times(z+1)=\frac{z+1}{1-z-z^2}$$
Agora aplicamos o teorema de transferência. ##f(z)=z+1## e ##g(z)=1-z-z^2##. As raízes de ##g(z)## são ##0.618## e ##-1.618##, a de menor módulo é ##0.618## com multiplicidade ##1##. Então:
$$\begin{align*} \beta &= 1/0.618 = 1.618 \\ f(1/\beta)&=f(0.618)=1.618 \\ g'(z) &=-2z-1 \\ g'(1/\beta) &=g'(0.618)=-2.236 \\ C[n]&\sim \frac{-1.618\times 1.618}{-2.236}\times 1.618^n \\ C[n] &\sim 1.171\times 1.618^n \end{align*}$$
Prontinho, transformamos raciocínio, que é díficil, em álgebra, que qualquer chimpanzé bêbado consegue fazer!

A pergunta original


Podemos agora voltar à pergunta original: se você jogar Cara ou Coroa várias vezes em sequência, qual a chance de nunca aparecer duas caras seguidas?

Ora, essa probabilidade é exatamente igual ao número de strings binárias que não possuem dois zeros seguidos, dividida pelo número total de strings binárias. Esses são os dois casos que analisamos, logo:
$$p\sim\frac{1.171\times1.618^n}{2^n} = 1.171\times 0.809^n$$
Easy! E funciona mesmo? Eu fiz uma simulação computacional para comparar com o nosso assintótico, eis o resultado:

Para ##n=5##, simulação ##40.78\%##, assintótico ##40.57\%##.

Para ##n=10##, simulação ##14.17\%##, assintótico ##14.06\%##.

Para ##n=15##, simulação ##4.84\%##, assintótico ##4.87\%##.

Script com a simulação no github

Ou seja, funcionou super bem, e a precisão vai ficando melhor conforme o ##n## aumenta.

Para saber mais


Se você gostou do que viu, a referência definitiva sobre o tema é o livro Analytic Combinatorics, do Flajolet, que foi o inventor da técnica. Não é o livro mais didático do mundo, mas está inteiramente disponível na web. No Coursera, de tempos em tempos, o Sedgewick faz um curso online, que eu super recomendo, foi onde eu aprendi inclusive. E na wikipedia... o tópico é tão novo que não tem na wikipedia em português ainda! Fica como exercício para o leitor criar uma página na wikipedia sobre o tema :)

segunda-feira, 31 de março de 2014

Matemática x Tanques de Guerra

Quando eu era criança pequena lá na década de 80, eu tinha um brinquedo de tanque de guerra. Ele chamava Panzer, e era fabricado pela brinquedos Mimo. Para a época era bem divertido: tinha controle remoto, e, além de andar, ainda disparava mísseis!

panzer

Mas eu não percebia a ironia do brinquedo. O Panzer não era só um brinquedo, ele existiu de verdade. Na vida real, o Panzer era uma máquina nazista de destruição. A série Panzer de tanques alemães foi responsável pela morte direta de incontáveis aliados. Ver um Panzer surgindo no horizonte era como ver a morte vindo em sua direção. Era um mimo esse Panzer.

Por outro lado, o Panzer era forte, mas não era invencível. Os aliados tinham armas capazes de deter o Panzer, o problema era descobrir quantas dessas armas eles deveriam levar para o front. Poucas semanas antes do Dia D, os alemães começaram a usar o novo modelo Panzer V, e os aliados não tinham como estimar quantas armas levar sem saber quantos Panzer V os inimigos fabricaram.

Foi então que os aliados usaram a arma definitiva contra o Panzer: a Matemática! Usando um truque muito esperto, eles estimaram que os alemães estavam fabricando 270 Panzers por mês. Depois que a guerra acabou, eles foram nos arquivos alemães conferir os números, e na verdade estavam fazendo 276 Panzers por mês. A estimativa foi muito boa!


A Quantidade de Tanques


O truque dos aliados se baseava em um erro dos alemães e em uma hipótese estatística. O erro dos alemães era que os tanques tinham número de série sequencial. Quando saíam da fabrica, eram numerados como 1, 2, 3, e assim por diante. Ao fazer isso, eles deixaram uma brecha para os hackers aliados.

Suponha que você capturou um Panzer alemão, e o serial dele era 14. Você consegue estimar quantos tanques eles produziram baseado só nesse número? No caso geral não, mas você pode fazer uma hipótese simplificadora, que é supor que o tanque que você capturou é aleatório (ou seja, os alemães estavam distribuindo igualmente os tanques por toda a Europa).

Com essa hipótese a intuição é clara. Vamos supor, hipoteticamente, que os alemães fizeram 200 tanques. Qual chance de você capturar um tanque de serial menor ou igual a 14 nessa situação? É 14/200 né? E se os alemães fizeram só 20 tanques? Aí a chance é 14/20, dez vezes maior que 14/200. Com a hipótese de tanque aleatório, é mais provável que os alemães tenham feito 20 tanques do que 200 tanques.

Você pode formalizar esse argumento usando o teorema de Bayes, e calcular o valor esperado do número de tanques. Vamos fazer as contas na caixa azul:

Suponha que conseguimos capturar k tanques. Dentre esses k tanques, o maior serial encontrado tem o valor m. Com essas informações, vamos calcular qual é o valor esperado do número de tanques produzidos, que vamos chamar de n. Podemos abrir esse valor esperado pela definição:
$$E[n\;|\;m,k] = \sum_{n} n \; p(n\; |\; m,k)$$
O problema então é achar a probabilidade condicional de n, dados m e k. Nós ainda não sabemos quanto vale essa probabilidade, mas podemos abri-la usando o teorema de Bayes:
$$p(n\;|\;m,k)=\frac{p(m,k\;|\;n)\;p(n)}{p(m,k)}$$
Essas três probabilidades nós conseguimos calcular! A primeira delas, p(m,k|n), sai por combinatória. Dado um vetor de tamanho n, de quantas maneiras conseguimos escolher k elementos, de modo que o maior deles seja m?

Um dos elementos precisa necessariamente ser o m, então só precisamos calcular a posição dos outros k-1 restantes. E eles precisam ser menores que m, então só podem estar nas m-1 posições menores que m. Por isso, essa quantidade vale binomial(m-1k-1). Note, entretanto, que isso só vale quando n>=m, caso contrário a probabilidade é zero (não tem como o máximo ser m se o vetor tem um tamanho menor que m).

Para calcular a probabilidade desejada, temos que achar de quantas maneiras podemos escolher k elementos dentre n, sem considerar a restrição do maior deles ser m. Mas isso é o próprio binomial(n, k). Então, a probabilidade é:
$$p(m,k\;|\;n) = \frac{\displaystyle{m-1\choose k-1}} {\displaystyle{n\choose k}}[n\ge m]$$
Vamos agora para a segunda, p(n). Eu não tenho nenhuma informação sobre quantos tanques foram produzidos, então não sei qual é a probabilidade de ter n tanques. Mas eu posso chutar que a distribuição é uniforme. A chance de ter 15 tanques, ou de ter 500 tanques, a priori, deve ser a mesma. Por isso, vamos fazer p(n) ser uma constante independente de n:
$$p(n)=c$$
A última, p(m,k), é a probabilidade do máximo ser m em k escolhas, independente do valor de n. Para isso, podemos usar a lei da probabilidade total
$$p(m,k) = \sum_n p(m,k\;|\;n) \;p(n)$$ 
Essas duas já sabemos calcular! É só substituir:
 $$\begin{align*} p(m,k) &= \sum_n p(m,k\;|\;n) \;p(n) \\ &= \sum_n \frac{{m-1\choose k-1}} {{n\choose k}}[n\ge m]\;c\\ &= \sum_{n\ge m} c {m-1\choose k-1} {n\choose k}^{-1}\\ &= c {m-1 \choose k-1} \sum_{n\ge m} {n\choose k}^{-1}\\ \end{align*} $$
A somatória tem uma forma fechada, que você pode achar com o algoritmo de Gosper:
$$\sum_{n=m}^{\infty}{n\choose k}^{-1}=\frac{m}{k-1}{m\choose k}^{-1}$$ 
Substituindo:
$$\begin{align*} p(m,k) &= c {m-1 \choose k-1} \sum_{n\ge m} {n\choose k}^{-1}\\ &= c {m-1 \choose k-1} \times \frac{m}{k-1}{m\choose k}^{-1} \\ &= \frac{cm}{k-1} {m-1 \choose k-1} {m\choose k}^{-1} \\ \end{align*} $$
Agora podemos voltar e substituir as três probalidades na fórmula de Bayes:
$$ \begin{align*} p(n\;|\;m,k)&=\frac{p(m,k\;|\;n)\;p(n)}{p(m,k)} \\ &= \frac{\displaystyle{m-1\choose k-1}{n\choose k}^{-1}[n\ge m]\times c} {\displaystyle\frac{cm}{k-1}{m-1\choose k-1}{m\choose k}^{-1}}\\ &= \frac{k-1}{m}{m\choose k}{n\choose k}^{-1}[n\ge m] \end{align*}$$
Por fim, podemos substituir na fórmula do valor esperado:
 $$\begin{align*} E[n\;|\;m,k] &= \sum_{n} n \; p(n\; |\; m,k)\\ &= \sum_n n \frac{k-1}{m}{m\choose k}{n\choose k}^{-1}[n\ge m]\\ &= \frac{k-1}{m}{m\choose k}\sum_{n\ge m} n {n\choose k}^{-1}\\ \end{align*}$$
Essa somatória também tem uma forma fechada por Gosper:
$$ \sum_{n=m}^{\infty} n{n\choose k}^{-1}=\frac{m(m-1)}{k-2}{m\choose k}^{-1}$$
Chegamos então na substituição final:
$$\begin{align*} E[n\;|\;m,k] &= \frac{k-1}{m}{m\choose k}\sum_{n\ge m} n {n\choose k}^{-1}\\ &= \frac{k-1}{m}{m\choose k}\times\frac{m(m-1)}{k-2}{m\choose k}^{-1} \\ &= \frac{(k-1)(m-1)}{k-2} \\ \end{align*}$$

A Estimativa na Prática


A fórmula final é super simples, mas ela funciona mesmo na prática? Podemos testá-la fazendo uma simulação numérica. Para isso, eu fiz dez mil simulações. Em cada uma eu sorteio o número de tanques fabricados, e escolho aleatoriamente 3 deles. Aí eu calculo a estimativa pela fórmula, e normalizo o erro. O código está no github, e o histograma resultante é o abaixo. A normalização é tal que o valor 0 é uma estimativa totalmente correta:

germantank3

A estimativa não ficou muito boa não! Ao invés de ter uma gaussiana em torno do zero, o estimador tem um bias muito forte para os negativos. Felizmente, esse problema só acontece porque capturamos apenas 3 tanques. Se, ao invés disso, tivéssemos capturado 30 tanques, então o histograma seria bem melhor. Compare o histograma de 30 tanques sobreposto ao de 3 tanques:

germantank30

Bem melhor né? Agora quase todas as estimativas estão bem próximas do zero. Esse método é bastante sensível ao número de tanques capturados. Os aliados sabiam disso, por isso que a estimativa deles foi tão boa (estimaram 270 tanques/mês, quando o valor real era 276). Sabe quantos tanques eles capturaram para fazer essa estimativa? Só dois tanques!

Claro que tem um truque. A estimativa com dois tanques é muito ruim, mas eles notaram que os alemães colocavam serial em todas as peças do tanque. Em especial, cada tanque tinha 48 rodas, cada uma com um serial único. Por isso, eles conseguiram usar a fórmula com k=96, o que deu a precisão alta que queriam. E no fim a matemática ganhou a guerra :)

domingo, 23 de fevereiro de 2014

O Paradoxo do Pokémon Twitch

A última moda da internet é o Pokémon Twitch. Como essas modas evaporam tão rapidamente quanto surgem, vale a pena explicar do que se trata:

poketwitch

O twitch é um site para livestreaming de videogames. Você conecta seu videogame no site, e as pessoas do mundo todo podem te assistir enquanto você joga. Ele tem também uma janela de chat, então as pessoas podem comentar enquanto assistem.

Mas algum gênio aprimorou a idéia original. Ele ligou o twitch em um emulador de Game Boy, e configurou o emulador para pegar o texto do chat e interpretar como se fosse o joystick. Com isso, qualquer um que tenha uma conta lá pode controlar o personagem do jogo (no caso, a versão original de Pokémon). E melhor ainda, como todos estão controlando o mesmo jogo, então o que ele fez foi criar uma espécie de crowdgaming, onde a sabedoria das massas escolhe o melhor caminho do personagem.

Quer dizer, na teoria. Na prática, a internet será a internet. Para cada um que tenta jogar sério, tem um trollzinho que fica mandando o comando oposto. No pico de popularidade, ele chegou a ter mais de cem mil pessoas jogando ao mesmo tempo. Imagine metade disso sacaneando ao invés de jogando, e dá pra ter uma idéia da loucura que é o Pokémon Twitch (ou então leia o FAQ para entender tudo que já aconteceu até agora).

pikachu_muitcholoco

Quando eu fiquei sabendo do Pokémon Twitch, eu fiz o que qualquer pessoa normal faria: criei um modelo matemático do jogo. No meu modelo simplificado, o personagem só anda na vertical, um passo por iteração. O número de jogadores e o número de trolls é o mesmo, então a cada passo ele tem 50% de chance de andar uma unidade para cima ou para baixo. Assuma que ele começa da origem. Nesse modelo, eu faço duas perguntas:
  1. Depois de n iterações, onde você espera que o personagem esteja?
  2. Depois de n iterações, qual você espera que seja a distância do jogador à origem?
Eu sei, eu sei: "Ricbit, deixa de ser bobo! As duas perguntas são iguais! Suponha que a resposta da primeira pergunta é y. Então a resposta da segunda é y menos 0, ou seja, o mesmo valor y".

Paradoxalmente, isso está errado! A resposta das duas perguntas é diferente! Probabilidade é um tópico que escapa facilmente da nossa intuição, então vale a pena fazer as contas com cuidado para entender a solução.

O valor esperado da posição


Para a resposta da primeira pergunta, vamos calcular qual é o valor esperado da posição. Na iteração 0 ele ainda não executou nenhum movimento, então sabemos a posição deterministicamente: ele está na origem.


E na iteração n+1? Depende de onde ele estava na iteração n. Tem 50% de chance de ter subido, e 50% de ter descido:


Com isso podemos calcular o valor esperado de y(n+1) direto da lei da expectativa total:


Ou seja, o valor esperado em qualquer iteração é sempre zero, na média nós esperamos que ele fique andando em círculos e nunca fique muito longe da origem.

Entretanto, note que só porque o valor esperado é zero não quer dizer que ele vai sempre terminar na origem: a média é zero, mas a variância não é. Um exercício bastante curioso é calcular qual a probabilidade exata do personagem estar na origem após n iterações. Isso é fácil e dá pra resolver com combinatória de colégio. Imagine que você tem n caixinhas:


Cada uma dessas caixinhas você pode preencher com +1 ou -1; se a soma de todas elas for zero, então o personagem termina na origem. De quantas maneiras podemos fazer isso? Note que nós só precisamos escolher as posições dos +1:

+1 +1 +1 +1

Sabendo onde estão os +1, a posição dos -1 restantes fica unicamente determinada. Precisamos então descobrir de quantas maneiras podemos encaixar n/2 números +1 em n caixinhas. Mas isso é a definição do binomial:


Agora é só dividir pelo total. Como cada caixinha tem duas opções possíveis, +1 e -1, então o número total de caixinhas é 2^n. Portanto, a probabilidade dele terminar na origem é:


Essa fórmula tem um problema: embora ela seja exata, é muito ruim de calcular quando o n é grande (tente calcular para n=1000, o número de dígitos da sua calculadora vai acabar rapidinho). Podemos achar um assintótico usando a aproximação do coeficiente binomial central:


Mais fácil né? Agora podemos calcular facilmente que, para n=1000, a chance de terminar na origem é aproximadamente 2.5%.

Se você tem o olho bom, deve ter percebido uma pegadinha na derivação da fórmula acima. Ela só funciona quando n é par! De fato, quando n é ímpar, não tem como o número de +1 e -1 serem iguais, e portanto a chance dele terminar na origem é zero. Olha que curioso: quando o n é ímpar, ele nunca termina em zero, apesar do valor esperado ser zero!

O valor esperado da distância


Vamos agora à segunda pergunta, que é o valor esperado da distância. Por que afinal esse valor é diferente do anterior?

Para exemplificar, vamos supor que n=2. Em uma das amostras, digamos que os jogadores ganharam e ele andou 2 unidades para cima, então a distância para a origem é 2. Em outra amostra, suponha que os trolls ganharam, então ele andou 2 unidades para baixo. Qual a distância da origem? É 2 também! Distâncias são sempre positivas!

Em outras palavras, o valor esperado da distância é o valor esperado do módulo da posição! E quanto é esse valor? As contas aqui são um pouco mais complicadas, então vou precisar da caixa azul:

Dessa vez nós vamos ter que calcular o valor esperado pela definição:


Quando k=0 o somando vale zero também. Vamos então supor k positivo e usar a mesma idéia das caixinhas. A soma era zero quando a quantidade de +1 e de -1 era igual. Dessa vez, nós queremos que a soma seja k, então precisamos ter k números +1 a mais que o número de -1. Portanto, a chance dele terminar em um valor k positivo é:


Quando k é negativo, a fórmula é exatamente a mesma! Afinal, é só trocar todos os +1 por -1 e vice versa. Portanto, o valor esperado da distância é:

Essa fórmula é difícil de manipular, para deixar mais fácil vamos supor que n é par. Se n é par, podemos escrevê-lo como n=2m. Mas olha só, se n for par, então você nunca vai terminar numa distância k que seja ímpar. Se você tirar k ímpar de n par, o que sobra é ímpar, então não tem como ter soma zero. Portanto podemos supor k par, e fazer k=2q. Substituindo:


A somatória é uma soma hipergeométrica, então você pode resolvê-la com o algoritmo de Gosper, chegando assim na forma fechada. Mas se você não souber usar o algoritmo de Gosper, não tem problema, é só provar a identidade abaixo por indução finita:


Substituindo a fórmula temos:

Agora é só usar a aproximação do binomial central:

Chegamos então no resultado final, o valor esperado da distância, que é raiz de 2n sobre pi.


A matemática experimental


Probabilidade tem uma vantagem grande sobre outros ramos da matemática: é muito fácil conferir as contas. Basta escrever uma simulação computacional com uma quantidade razoável de amostras!
Eu escrevi um script em python (o código está no github), e os resultados para n=1000 e cem mil amostras foram os seguintes:

  • Posição: calculado 0.00, medido 0.05
  • Distância: calculado 25.23, medido 25.27
  • Término na origem: calculado 2.52%, medido 2.53%

Ou seja, a matemática experimental concordou com a matemática teórica.

A lição que fica dessa brincadeira é tomar muito cuidado com suas intuições sobre probabilidade, porque ela tem uma tendência a nos enganar bem facilmente :)