A Marcha dos Lemmings

Lemmings são pequenos roedores que vivem na tundra ártica. Eles possuem pelos longos e macios, rabo curto, são herbívoros e explodem.

Sim, eles explodem! Os lemmings possuem a capacidade de ligar um contador luminoso sobre suas cabeças. A cada passo que dão, o contador decrementa, até que finalmente explodem quando o contador chega a zero.

O que você faria se soubesse que vai explodir em seis passos? Parece que não tem muitos passeios possíveis com exatamente seis passos, mas a intuição é meio falha aqui. Na verdade, o número de passeios possíveis é exponencial no número de passos! Como fazer um algoritmo que calcula o número de caminhos possíveis que o lemming pode fazer antes de explodir?

As Regras do Jogo

Para modelar o problema, vamos supor que os passos que lemming pode fazer são descritos por um grafo orientado. Por exemplo, no grafo abaixo o lemming começa no nó A. Do A ele pode ir para B ou C, de B ele pode ir para C, e de C ele pode ir para A.

Se o nosso lemming pode dar seis passos, então existem sete caminhos diferentes que ele pode fazer:

ABCABCA

ABCACAB

ABCACAC

ACABCAB

ACABCAC

ACACABC

ACACACA

Logo ##p(6)=7##, e podemos até tabelar ##p(n)## para todos os ##n## de 0 a 6:

$$p(0..6)=1,2,2,3,4,5,7$$
Dado ##n##, como fazer um algoritmo que calcule ##p(n)##?

Eu vou mostrar três algoritmos possíveis usando a linguagem do Wolfram Mathematica. Para isso, o grafo precisa ser descrito de alguma maneira que os algoritmos entendam. Eu vou usar uma matriz de adjacências: para um grafo com ##k## nós, a matriz é ##k\times k## e o elemento na linha ##i## da coluna ##j## é 1 se existe uma ligação saindo do nó ##j## e indo em direção ao nó ##i##, e 0 em caso contrário. Para o nosso exemplo, a matriz é:

$$M=\left[\begin{matrix}0 && 0 && 1 \\ 1 && 0 && 0 \\ 1 && 1 && 0\end{matrix}\right]$$

Força Bruta

Sempre que eu tenho que fazer um algoritmo do zero, prefiro começar com a força bruta. É fácil e rápido de implementar, e embora não seja eficiente, serve para conferir o resultado dos algoritmos mais avançados.

Nesse caso, para implementar a força bruta você passa o grafo, o nó em que o lemming está, quantos passos faltam para terminar, e o resto é uma recursão simples. Como toda recursão, não esqueça de fazer o caso base! Para esse problema, o caso base é ##p(0)=1##, ou seja, tem um único caminho possível sem dar passo nenhum, que é não sair do lugar.

No fim, a solução em força bruta fica assim:

count[graph_, pos_, 0] := 1
count[graph_, pos_, size_] := Sum[
    If[graph[[i, pos]] == 1, count[graph, i, size - 1], 0],
    {i, 1, Length[graph]}]

Vamos conferir o resultado:

graph := {{0, 0, 1}, {1, 0, 0}, {1, 1, 0}}
Print[Table[count[graph, 1, n], {n, 0, 6}]]

{1,2,2,3,4,5,7}

Perfeito! Mas qual é a complexidade desse algoritmo? Bem, ele precisa visitar cada um dos caminhos possíveis, e o número total de caminhos é exponencial, logo a complexidade vai ser do tipo ##O(c^n)##, onde a constante ##c## depende do grafo. Ou seja, você não quer usar esse algoritmo quando o ##n## é grande.

Multiplicação de Matrizes

Como melhorar então? Tem uma maneira calcular ##p(n)## mais rápido usando multiplicação de matrizes. É mais fácil mostrar como funciona por exemplos: no estado inicial, nós temos um lemming no nó A, então vamos montar um vetor correspondente: ##V=[1 \;0\; 0]^T##. Se multiplicarmos ##M## por ##V##, qual o resultado?

$$M V=\left[\begin{matrix}0 && 0 && 1 \\ 1 && 0 && 0 \\ 1 && 1 && 0\end{matrix}\right] \left[\begin{matrix}1\\0\\0\end{matrix}\right]=\left[\begin{matrix}0\\1\\1\end{matrix}\right]$$

Ou seja, um lemming que está em A vai parar em B ou C. Continuando:

$$M V=\left[\begin{matrix}0 && 0 && 1 \\ 1 && 0 && 0 \\ 1 && 1 && 0\end{matrix}\right] \left[\begin{matrix}0\\1\\1\end{matrix}\right]=\left[\begin{matrix}1\\0\\1\end{matrix}\right]$$
Dado um lemming em B ou C, o resultado é um lemming em A ou C. Se o vetor de entrada mostra quantos lemmings tem em cada nó, multiplicar pela matriz M vai mostrar onde eles podem estar depois de dar um passo.

Agora é só estender o raciocínio. Se multiplicar por M correponde a um passo, então multiplicar ##n## vezes por M é o equivalente a ##n## passos. E se você começou o processo com um vetor unitário, então a soma do elementos de V é o número de caminhos possíveis partindo do nó inicial!

Uma implementação possível é a abaixo:

count[graph_, pos_, size_] :=
  Total[MatrixPower[graph, size] . UnitVector[Length[graph], pos]]

Conferindo:

graph := {{0, 0, 1}, {1, 0, 0}, {1, 1, 0}}
Print[Table[count[graph, 1, n], {n, 0, 6}]]

{1,2,2,3,4,5,7}

Bateu certinho novamente! E qual é a complexidade? Se você implementou a multiplicação e a exponenciação de matrizes do jeito mais simples, vai ser da ordem de ##O(k^3 n)##, que certamente é melhor que exponencial mas ainda não é boa o suficiente para usar na prática.

No entanto, você pode melhorar isso com implementações melhores! A multiplicação de matrizes pode usar o algoritmo de Strassen, e a exponenciação pode usar o método binário, aí o total fica ##O(k^{2.8} \log_2 n)##, bem melhor; e se a matriz for esparsa dá para ficar ainda mais rápido.

Função Geradora

Mas tem um método que supera o anterior, que é usar funções geradoras! Nós começamos tudo com um grafo orientado, e todo grafo orientado é equivalente a um autômato finito. Por sua vez, todo autômato finito é equivalente a uma gramática regular, então deve ter um jeito de descrever os caminhos do lemming usando uma gramática. Para o exemplo dado, a gramática fica assim:

$$\begin{align*}\alpha &\to A \;|\; A \beta \;|\; A \gamma \\ \beta &\to B \;|\;  B \gamma \\ \gamma &\to C \;|\; C \alpha\end{align*}$$
Todos os caminhos possíveis produzem strings que são aceitas por essa gramática, e a gramática não aceita nenhuma string que não seja um caminho possível.

O truque agora é usar combinatória analítica, que para o problema de enumerar strings de uma gramática é bem simples: basta trocar cada token por um ##z##, e cada OR por uma adição. A gramática vira um sistema de equações:

$$\begin{align*}\alpha &= z + z\beta +z \gamma \\ \beta &= z+z \gamma \\ \gamma &= z+z \alpha\end{align*}$$
Resolvendo para ##\alpha##, chegamos na função geradora:

$$\alpha=\frac{z + 2 z^2 + z^3}{1 - z^2 - z^3}$$
E para que serve a função geradora? Olhe só o que acontece quando abrimos a função em série de potências:

$$\frac{z + 2 z^2 + z^3}{1 - z^2 - z^3}=z+2z^2+2z^3+3z^4+4z^5+5z^6+7z^7+\dots$$

Ahá! Os coeficientes de ##z## são exatamente o número de caminhos no grafo! Para calcular ##p(6)##, basta olhar o coeficiente de ##z^7## (você precisa somar um porque a função geradora está contando o número de tokens na strings, e nós queremos o número de passos).

Podemos ir além agora. Quando a função geradora é uma função racional, você sempre consegue inverter e conseguir uma fórmula explícita para ##p(n)##:

$$p(n)=0.957\times 1.325^n+0.426\times 0.869^n\cos(1.469+ 2.438 n)$$
Essa fórmula é exata (ou seria, se eu não tivesse truncado os números para caber na tela). Note que o resultado sempre vai ser um inteiro, e que a parte oscilante da fórmula sempre vai ser menor que 0.4, então dá para simplificar mais ainda:

$$p(n)=\left\lfloor 0.5+0.957\times 1.325^n\right\rfloor$$
Basta pegar a parte exponencial e arredondar para o inteiro mais próximo. Podemos observar duas coisas aqui. Primeiro, agora nós temos um algoritmo onde o número de operações independe do valor de ##n##, então efetivamente esse algoritmo é ##O(1)##. Na prática você vai implementar usando bignum, aí não fica ##O(1)## de verdade; mas os outros algoritmos mostrados também vão ficar mais lentos na mesma proporção.

Segundo, nós provamos que aquele primeiro algoritmo era de fato ##O(c^n)##! Conseguimos até calcular quem é o tal do ##c##: é aproximadamente 1.325; ou, se você quiser o valor exato com radicais:

$$c=\frac{\sqrt[3]{9+\sqrt{69}}+\sqrt[3]{9-\sqrt{69}}}{\sqrt[3]{18}}$$
Para implementar esse método computacionalmente, é só reescrever aquele sistema de equações em forma matricial:

$$\left[\begin{matrix}\alpha\\ \beta \\ \gamma\end{matrix}\right] = z \left[\begin{matrix}0 && 1 && 1\\0 && 0 && 1\\ 1 &&0 && 0\end{matrix}\right] \left[\begin{matrix}\alpha\\ \beta \\ \gamma\end{matrix}\right]  + z\left[\begin{matrix}1\\ 1\\1\end{matrix}\right] $$

A matriz ali no meio é a transposta da matriz de adjacências! Agora já podemos implementar:

count[graph_, pos_, size_] := With[{n = Length[graph]},
  SeriesCoefficient[LinearSolve[
      IdentityMatrix[n] - z Transpose[graph],
      ConstantArray[z, n]][[pos]],
    {z, 0, size + 1}]]

Como esperado, a resposta é correta:

graph := {{0, 0, 1}, {1, 0, 0}, {1, 1, 0}}
Print[Table[count[graph, 1, n], {n, 0, 6}]]

{1,2,2,3,4,5,7}

Um Problema em Aberto

Nós chegamos em um algoritmo que é ##O(1)##, dá para otimizar mais? Não! Esse é um limite teórico.

Quer dizer, existem algoritmos que são ##O(0)##. Por exemplo, "crie um algoritmo que, dado uma circunferência e um diâmetro, calcule a razão entre os dois". Esse é algoritmo é ##O(0)## porque a resposta é sempre ##\pi##, independe do input. Mas, se a resposta não é constante, ##O(1)## é o melhor algoritmo possível. (A não ser que você crie um algoritmo que dê a resposta antes do input chegar. Nesse caso me ensine como faz, porque você inventou a viagem no tempo).

Mas existe uma pequena variação no problema dos lemmings que o torna muito, muito mais difícil. Suponha que agora o lemming não pode passar duas vezes pelo mesmo nó. A força bruta para essa variação ainda é ##O(c^n)##. Mas ninguém bolou um algoritmo que seja melhor que a força bruta nesse caso!

Quando o grafo é finito, você pode usar a força bruta para enumerar todas as soluções, já que eventualmente o grafo acaba. Mas se o grafo for infinito, aí danou-se. Para esse caso, tem gente procurando a solução desde a década de 50. Já são setenta anos procurando o algoritmo e ninguém achou nada para um grafo genérico, nem mesmo um assintótico! O melhor resultado até hoje é um de 2011 onde acharam o ##c## para um grafo específico, o lattice hexagonal. Para todos os outros grafos, o melhor que temos são heurísticas e aproximações numéricas.

Se você estiver sem nada o que fazer no próximo fim de semana, essa pode ser uma boa maneira de ocupar o tempo!

O dia que o Knuth ganhou meu cheque

No último dia 10, o Donald Knuth completou 80 anos de vida. Para comemorar, ele fez uma festinha na remota e gelada cidadezinha de Piteå, no norte da Suécia, e eu tive a sorte de poder participar. A festa foi em formato de simpósio de matemática: ao invés de presentes, cada amigo apresentou uma palestra sobre um tema relacionado ao trabalho do Knuth. (E os amigos dele são os caras mais famosos da computação, estavam lá o Sedgewick, o Karp, o Tarjan, entre outros).

Eu também fiz minha contribuição, presenteei o Knuth com vários artigos contendo idéias para futuras revisões dos livros dele. Mas em um dos artigos aconteceu uma coisa engraçada: o Knuth notou que uma das minhas contas podia ser simplificada, e como agradecimento pela sugestão eu tive a oportunidade de inverter expectativas, e dar um cheque para o Knuth :D

A Soma dos Quadrados

O artigo em questão era uma sugestão para o Concrete Mathematics, que é o meu livro de matemática predileto. Sou fã desse livro desde os 17 anos, quando o usava como referência para estudar para a Olimpíada de Matemática.

Quem já leu sabe que o livro contém mais de dez demonstrações diferentes da fórmula para a soma dos quadrados. Para ##n\ge 0##, a fórmula é:
$$\sum_{1\le x\le n} x^2=\frac{n(n+1)(2n+1)}{6}$$
As demonstrações vão desde as mais simples (usando indução finita), até as mais complexas (usando a fórmula de Euler-MacLaurin). Porém, uns anos atrás, eu bolei uma demonstração nova que é provavelmente a mais curta demonstração de todas! Ela está na caixa azul abaixo:

Testando valores pequenos, notamos que a fórmula é verdadeira para ##n=0##, ##n=1##, ##n=2##, ##n=3##. Logo, é verdadeira para qualquer ##n##.

Sim, eu sei o que você está pensando: "Ricbit, você tá doido? Só porque funciona para alguns números não quer dizer que funciona para todos os números!".

A preocupação é válida: o que mais tem na matemática são fórmulas que funcionam para números pequenos mas falham para números grandes. Um exemplo famoso é o polinômio sortudo de Euler, ##k^2-k+41##, que produz apenas números primos para ##k=1##, ##k=2##, ##k=3##, etc., mas falha lá na frente quando ##k=41##.

Porém, nesse caso, testar números pequenos é suficiente para demonstrar a fórmula! Para entender o motivo precisamos de dois fatos sobre cálculo e polinômios.

O Cálculo Discreto

Quem já estudou cálculo com certeza sabe algumas integrais de cabeça. Por exemplo, a integral do monômio:
$$\int x^n\;dx=\frac{x^{n+1}}{n+1}$$
O que você talvez não saiba é que o cálculo que a gente aprende nas faculdades de exatas é só um tipo de cálculo: o Cálculo Contínuo. Existem outros tipos, como o Cálculo Discreto, que lida com somatórias ao invés de integrais. Esse cálculo possui várias fórmulas análogas ao cálculo mais comum. Em especial, ele tem uma fórmula análoga à integral do monômio:
$$\sum x^{\underline{n}}\;\delta x=\frac{x^{\underline{n+1}}}{n+1}$$
Nesse fórmula o monômio não usa a potência simples, ele usa a potência fatorial, que é definida da seguinte maneira:
$$x^{\underline n}=x(x-1)(x-2)\dots(x-n+1)$$
É como se fosse um fatorial, mas ao invés de ir até o fim, você pega só os ##n## primeiros termos. (Para pronunciar ##x^{\underline n}##,  você fala "x elevado a n caindo").

Para converter uma potência fatorial em uma potência normal, você pode abrir a definição e multiplicar todos os termos, mas isso dá trabalho. É mais fácil ter à mão uma tabela com os números de Stirling (que vêm em dois tipos, o primeiro tipo ##{n\brack k}##, e o segundo tipo ##{n\brace k}##). Esses números são fáceis de calcular porque eles obedecem uma regra similar ao triângulo de Pascal.  Tendo a tabela, as fórmulas abaixo fazem a conversão:
$$\begin{align*}
x^{\underline{n}}&=\sum_k{n\brack k}(-1)^{n-k}x^k\\
x^{n}&=\sum_k{n\brace k}x^{\underline{k}}\end{align*}$$
Usando as fórmulas acima e alguma paciência, você consegue demonstrar a fórmula da soma dos quadrados (converta ##x^2## em potências fatoriais, use a fórmula de somatória do cálculo discreto, depois converta de volta as potências fatoriais em potências tradicionais).

Mas isso dá muito trabalho. Ao invés disso, note que as fórmulas acima tem uma consequência interessante: com elas é possível ver que a somatória de um polinômio de grau ##n## sempre vai ser um polinômio de grau ##n+1##, e em especial a somatória de ##x^2## vai ser um polinômio de grau 3, a gente só não sabe qual polinômio ainda. Mas aí temos outro truque nas mangas!

Interpolação de polinômios

Certo, eu não sei qual é o polinômio de grau 3 que resolve a somatória. Vamos então escrever esse polinômio na forma ##P(n)=an^3+bn^2+cn+d##. Apesar de não conhecermos o polinômio, é fácil descobrir os pontos por onde ele passa, basta calcular a somatória!

Para valores pequenos de ##n##, podemos tabelar ##P(0)=0##, ##P(1)=1##, ##P(2)=1+4=5##, ##P(3)=1+4+9=14##, agora podemos achar os coeficientes usando álgebra linear:
$$\left[\begin{array}{c}0\\1\\5\\14\end{array}\right]=
\left[\begin{array}{cccc}0^3&0^2&0^1&1\\1^3&1^2&1^1&1\\2^3&2^2&2^1&1\\3^3&3^2&3^1&1\end{array}\right]
\left[\begin{array}{c}a\\b\\c\\d\end{array}\right]$$
Novamente, dá muito trabalho. Ao invés disso, notamos que um polinômio de grau ##n## fica completamente determinado por ##n+1## pontos, e agora podemos finalmente entender a demonstração inicial!

Olhe novamente para o que queremos demonstrar:
$$\sum_{1\le x\le n} x^2=\frac{n(n+1)(2n+1)}{6}$$
De um lado, temos uma somatória de polinômios de grau 2, que sabemos que é um polinômio de grau 3. Do outro lado, temos um polinômio cujo grau também é 3. Para provar que esses dois polinômios são iguais, é suficiente testar os dois polinômios em quatro pontos. Escolhemos os mais fáceis que são ##n=0##, ##n=1##, ##n=2##, ##n=3##, como a fórmula funciona para esses quatro valores, então funciona para todos os valores, QED.

O cheque do Knuth

Quando eu presenteei o artigo com essa demonstração para o Knuth, em segundos após a leitura ele comentou: "Você podia ter usado -1 aqui!".

É verdade! Essa demonstração é tão curta que dá para fazer de cabeça, mas as contas para ##n=3## ficam meio chatinhas. Ao invés disso, você pode usar -1 no lugar de 3, e aí fica bem mais fácil: do lado esquerdo a soma é vazia e dá 0, do lado direito o termo ##(n+1)## zera, logo o resultado dá zero também. Você precisa expandir o domínio da fórmula, de ##n\ge 0## para ##n\ge -1##, o que significa que a fórmula provada com -1 é até melhor que a fórmula provada com 0!

Tradicionalmente, o Knuth oferece $0x1 hexadólar para quem acha um erro em um artigo dele, e $0x0.20 para sugestões. Pela simplificação que ele sugeriu, eu achei pertinente oferecer um cheque também!

Aqui o detalhe do cheque. O Knuth emite cheques pelo banco fictício de San Serriffe (eu tenho $0x5.80 lá atualmente). Os meus cheques eu resolvi emitir pelo Banco de Dados:

O Knuth adorou o cheque! Ele ficou impressionado em como eu consegui fazer tão rápido um cheque que parece mesmo um cheque, mas o segredo é simples: quem fez fez a arte foi a Ila Fox, eu só imprimi no hotel onde estava, usando papel reciclado.

Sem dúvida foi a melhor festa de aniversário que já participei! O Knuth já tinha feito uma festa dessas quando completou 64 anos, e para manter a simetria ele disse que vai fazer outra quando completar 96 anos. Eu certamente estarei lá mais uma vez :)

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!



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:

	##\enclose{box}{\enclose{circle}{1} \,\enclose{circle}{3}\, \enclose{circle}{6}}##
	##\enclose{box}{\enclose{circle}{2}\, \enclose{circle}{5} }##
	##\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^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^k}{k^n} \\ &=\sum_{n\ge 0} 1-n^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^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.

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!

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 :)

A Busca Por Bozo Binário

Se você trabalhar por tempo o suficiente com Computação, eventualmente uma das suas atribuições vai ser entrevistar candidatos para a sua equipe. Entrevistar é uma arte: existem inúmeros truques e macetes para descobrir se é esse sujeito na sua frente é mesmo o cara que você quer como seu parceiro. O Joel tem uma boa introdução sobre o tema que vale a pena ler.

Uma das habilidades que você precisa desenvolver rápido é aprender que o candidato não pensa como você. Em especial, quando ele responde uma pergunta de um jeito diferente daquele que você esperava, não quer dizer necessariamente que ele esteja errado. E esse caso eu já vi de perto, com o algoritmo que eu apelidei de Busca Por Bozo Binário.

Essa história aconteceu faz algum tempo. Eu estava na cantina conversando sobre entrevistas, quando alguém comentou que rejeitou um candidato porque ele não conseguiu implementar uma busca binária. Poxa, isso é grave mesmo! E qual o algoritmo errado que ele fez? A resposta foi essa:

"Sorteie um número de 1 até o tamanho do vetor, e verifique o elemento que tinha esse índice. Se for igual você achou, então retorna. Se for menor, repete para o lado esquerdo; se for maior, repete para o lado direito."

Isso é errado, me disseram, porque o pior caso é ##O(n)##. De fato, se o elemento que você quer achar é o primeiro, e o gerador de aleatórios retornar sempre o último, ele vai passar pelo vetor inteiro.

Mas nesse ponto eu tive que interromper: "Peraí! Você explicitou que queria ##O(\log n)## no pior caso? Se você não falou nada no enunciado, ele pode ter entendido que era ##O(\log n)## no caso médio, e aí o algoritmo dele está certo!"

Na hora ninguém botou fé que esse algoritmo de Busca por Bozo Binário era de fato ##O(\log n)##, então eu tive que voltar para casa e escrever a demonstração. Se você tem medo de Matemática Discreta, pule a caixa azul:

Para a Busca por Bozo Binário, nós queremos calcular qual é o valor médio do número de comparações que ele realiza. E como você começa uma análise de caso médio?

Primeiro você define as características da entrada: eu vou assumir que todos os casos estão distribuídos uniformemente. Depois, é só usar a fórmula do valor esperado: você soma o número de comparações em cada caso, multiplicado pela probabilidade daquele caso ocorrer. $$ F[x] = \sum_i p[i]C[i] $$ Antes de começar, vejamos os casos extremos. Quando o vetor é vazio, não precisa de nenhuma comparação, então ##F[0]=0##. Se o vetor tiver tamanho um, então só precisa de uma comparação, ##F[1]=1##. É sempre bom ter uns casos pequenos para conferir no final.

Vamos ver agora o caso de um vetor de tamanho ##n##. Eu não sei qual número vai ser escolhido para cortar o vetor em dois, é o gerador de números aleatórios que escolhe. Então eu vou tirar a média sobre todas as escolhas possíveis. Como os ##n## valores são equiprováveis, então a probabilidade individual é ##1/n##: $$F[n] = \sum_{0\le k\lt n}\frac{1}{n}F[k,n] = \frac{1}{n}\sum_{0\le k\lt n}F[k,n] $$ Na fórmula acima, ##F[n]## é o número médio de comparações para um vetor de tamanho ##n##, e ##F[k,n]## é o número médio de comparações para o vetor de tamanho ##n##, no caso em que o corte foi feito no elemento ##k##.

Vamos calcular ##F[k,n]## agora. Esse valor depende de onde está o número que estamos procurando. Eu não sei onde está o número; como estamos calculando o caso médio, precisamos tirar a média sobre todos as posições possíveis. $$F[k,n]=\sum_{0\le i\lt n}\frac{1}{n}F[i,k,n]=\frac{1}{n}\sum_{0\le i\lt n}F[i,k,n]$$ Agora ##F[i,k,n]## significa o valor médio para um vetor de tamanho ##n##, cortado na posição ##k##, e com número final na posição ##i##.

Nós temos três casos agora. Se ##i=k##, então você só precisa comparar uma vez, a comparação vai ser verdadeira e o algoritmo termina. Se ##i\lt k##, então você compara uma vez e faz recursão para o lado esquerdo. Se ##i\gt k##, então você compara uma vez e faz recursão para o lado direito. Resumindo: $$F[i,k,n]=\begin{cases} 1+F[k] &(i \lt k) \\ 1 &(i = k) \\ 1+F[n-k-1] &(i \gt k) \\ \end{cases}$$ Podemos sintetizar tudo que vimos em uma única recorrência. $$\begin{align*} F[0] &= 0 \\ F[n] &= \frac{1}{n} \sum_{0\le k \lt n}\left( \frac{1}{n}\sum_{0\le i\lt k}\left(1+ F[k]\right) + \frac{1}{n} + \frac{1}{n}\sum_{k\lt i\lt n}\left(1+ F[n-k-1]\right) \right) \end{align*}$$ Agora é só resolver!

Antes de mais nada, vamos coletar quem é constante em relação ao índice da somatória. Todos os ##1/n## coletam, e as somatórias mais internas viram constantes. $$\begin{align*} F[n] &= \frac{1}{n^2} \sum_{0\le k \lt n}\left( \sum_{0\le i\lt k}\left(1+ F[k]\right) + 1 + \sum_{k\lt i\lt n}\left(1+ F[n-k-1]\right) \right) \\ &= \frac{1}{n^2} \sum_{0\le k \lt n} k \left(1+ F[k]\right) + 1 + \left(n-k-1\right)\left(1+ F[n-k-1]\right) \\ &= \frac{1}{n^2} \sum_{0\le k \lt n} k + k F[k] + 1 + n-k-1 +\left(n-k-1\right)F[n-k-1] \\ &= \frac{1}{n^2} \sum_{0\le k \lt n} n + k F[k] +\left(n-k-1\right)F[n-k-1] \\ &= 1+\frac{1}{n^2} \left(\sum_{0\le k \lt n} k F[k] +\sum_{0\le k \lt n} \left(n-k-1\right)F[n-k-1] \right) \\ \end{align*}$$ Vamos focar naquela segunda somatória. Primeiro fazemos ##q=n-k-1##: $$ \sum_{0\le k \lt n} \left(n-k-1\right)F[n-k-1] =\sum_{0\le n-q-1 \lt n} qF[q] $$ Agora é só manipular os índices: $$\begin{align*} 0\le n-&q-1 \lt n \\ 1-n\le -&q \lt 1 \\ -1\lt &q \le n-1 \\ 0\le &q \lt n \\ \end{align*}$$ Olha só, a segunda somatória é igual à primeira! $$ \sum_{0\le n-q-1 \lt n} qF[q] =\sum_{0\le q \lt n} qF[q] =\sum_{0\le k \lt n} kF[k] $$ Substituindo, chegamos em uma fórmula curta para ##F[n]##: $$\begin{align*} F[n] &= 1+\frac{1}{n^2} \left(\sum_{0\le k \lt n} k F[k] +\sum_{0\le k \lt n} \left(n-k-1\right)F[n-k-1] \right) \\ &= 1+\frac{1}{n^2} \left(\sum_{0\le k \lt n} k F[k] +\sum_{0\le k \lt n} kF[k] \right) \\ &= 1+\frac{2}{n^2} \sum_{0\le k \lt n} k F[k] \\ \end{align*}$$ A somatória ainda está atrapalhando, o ideal seria sumir com ela. Uma maneira de fazer isso é isolando: $$\begin{align*} F[n] &= 1+\frac{2}{n^2} \sum_{0\le k \lt n} k F[k] \\ \sum_{0\le k \lt n} k F[k] &= \frac{n^2}{2}\left(F[n]-1\right) \end{align*}$$ Essa última fórmula vale para todo ##n##, em especial vale também para ##n-1##: $$\begin{align*} \sum_{0\le k \lt n} k F[k] &= \frac{n^2}{2}\left(F[n]-1\right) \\ \sum_{0\le k \lt n-1} k F[k] &= \frac{(n-1)^2}{2}\left(F[n-1]-1\right) \\ \end{align*}$$ Ahá! Agora é só subtrair uma da outra que a somatória desaparece! $$\begin{align*} \sum_{0\le k \lt n} k F[k] - \sum_{0\le k \lt n-1} k F[k] &= \frac{n^2}{2}\left(F[n]-1\right) - \frac{(n-1)^2}{2}\left(F[n-1]-1\right)\\ (n-1) F[n-1] &= \frac{n^2}{2}\left(F[n]-1\right) - \frac{(n-1)^2}{2}\left(F[n-1]-1\right)\\ (n-1) F[n-1] &= \frac{n^2}{2}F[n]-\frac{n^2}{2} - \frac{(n-1)^2}{2}F[n-1]+\frac{(n-1)^2}{2}\\ \frac{n^2}{2}F[n] &= \left(\frac{(n-1)^2}{2}+(n-1)\right)F[n-1]+\left(\frac{n^2}{2} -\frac{(n-1)^2}{2}\right)\\ n^2F[n] &= \left(n^2-1\right)F[n-1]+(2n-1)\\ \end{align*}$$ Chegamos finalmente em uma recorrência sem somatórias. Melhor ainda, essa recorrência está no ponto certo para usar a técnica do summation factor! Como esse é um truque muito útil de se conhecer, eu vou fazer em câmera lenta. Você pode usar um summation factor sempre que sua recorrência for da forma abaixo: $$a_n F[n] = b_n F[n-1]+c_n$$ Esse é o nosso caso, se usarmos a correspondência abaixo: $$\begin{align*} a_n &= n^2 \\ b_n &= n^2-1 \\ c_n &= 2n-1 \end{align*}$$ O segredo do summation faction é tentar achar uma sequência mágica ##s_n## que tenha a seguinte propriedade: $$s_n b_n=s_{n-1}a_{n-1}$$ E por que isso é bom? Olha só o que acontece quando você multiplica a recorrência original por ##s_n##: $$\begin{align*} a_n F[n] &= b_n F[n-1]+c_n \\ s_n a_n F[n] &= s_n b_n F[n-1]+s_n c_n \\ s_n a_n F[n] &= s_{n-1}a_{n-1}F[n-1]+s_n c_n \\ \end{align*}$$ Agora, se você substituir ##T[n]=s_n a_n F[n]##, temos: $$\begin{align*} s_n a_n F[n] &= s_{n-1}a_{n-1}F[n-1]+s_n c_n \\ T[n] &= T[n-1]+s_n c_n \end{align*}$$ Opa, essa é fácil! Dá pra resolver de cabeça, é praticamente a definição de somatória! $$\begin{align*} T[n] &= T[n-1]+s_n c_n\\ T[n] &= T[0] + \sum_{1\le k \le n} s_k c_k \\ s_n a_n F[n] &= s_0 a_0 F[0] + \sum_{1\le k \le n} s_k c_k \\ F[n] &= \frac{1}{s_n a_n} \left(s_0 a_0 F[0] + \sum_{1\le k \le n} s_k c_k \right) \\ \end{align*}$$ Pronto, agora não é mais uma recorrência, é uma fórmula simples. A dificuldade do método é encontrar o tal do ##s_n##. Para achá-lo você pode usar a intuição, chutar, fazer um sacrifício aos deuses; ou então você pode usar o método do porradão. Nesse método você começa com ##s_n=s_{n-1}\left(a_{n-1}/b_n\right)## e abre recursivamente, chegando no monstrinho abaixo: $$s_n = \frac{a_{n-1}a_{n-2}\ldots a_2 a_1}{b_n b_{n-1}\ldots b_3 b_2}$$ Essa fórmula parece grande, e é grande mesmo. Mas na maioria das vezes tem um monte de cancelamentos, o que deixa o resultado final pequeno. Vamos ver no nosso caso como fica: $$\begin{align*} s_n &= \frac{a_{n-1}a_{n-2}\ldots a_2 a_1}{b_n b_{n-1}\ldots b_3 b_2} \\ &= \frac{(n-1)^2 (n-2)^2 \ldots 2^2 1^2} {\left(n^2-1\right) \left(\left(n-1\right)^2-1\right) \ldots (3^2-1) (2^2-1)} \\ &= \frac{(n-1)^2 (n-2)^2 \ldots 2^2 1^2} {\left((n+1)(n-1)\right) \left((n-1+1)(n-1-1)\right) \ldots (3+1)(3-1) (2+1)(2-1)} \\ &= \frac{(n-1) (n-2)(n-3) \ldots3\cdot 2\cdot 1} {(n+1) (n) (n-1)\ldots (4+1) (3+1) (2+1)} \\ &= \frac{2} {(n+1) (n) } \\ \end{align*} $$ Cancelou mesmo! Bastou lembrar que ##(a^2-b^2)=(a+b)(a-b)## que o resto saiu. Agora é só lembrar que ##F[0]=0## e substituir na fórmula: $$\begin{align*} F[n] &= \frac{1}{s_n a_n} \left(s_0 a_0 F[0] + \sum_{1\le k \le n} s_k c_k \right) \\ F[n] &= \frac{n(n+1)}{2n^2} \sum_{1\le k \le n} \frac{2(2k-1)}{k(k+1)}\\ \end{align*}$$ Essa é a solução da recorrência. Dá para achar uma forma fechada? Nesse caso dá, desde que você tope uma forma fechada em função dos números harmônicos ##H[n]##, que são definidos como: $$H[n]=\sum_{1\le k \le n}\frac{1}{k}$$ Você começa abrindo a fração no somatório: $$\begin{align*} \frac{2(2k-1)}{k(k+1)} &= \frac{A}{k} + \frac{B}{k+1}\\ \frac{4k-2}{k(k+1)} &= \frac{A(k+1)+B k}{k(k+1)}\\ 4k-2 &= Ak+A +B k \\ 4k-2 &= (A+B)k+A \\ A &= -2 \\ B &= 6 \end{align*}$$ Agora você divide o somatório em dois: $$\begin{align*} F[n] &= \frac{n(n+1)}{2n^2} \sum_{1\le k \le n} \frac{2(2k-1)}{k(k+1)}\\ F[n] &= \frac{n(n+1)}{2n^2} \sum_{1\le k \le n} \frac{6}{k+1} - \frac{2}{k}\\ F[n] &= \frac{n(n+1)}{2n^2} \left( \sum_{1\le k \le n} \frac{6}{k+1} - \sum_{1\le k \le n} \frac{2}{k} \right) \\ F[n] &= \frac{n(n+1)}{2n^2} \left( \sum_{1\le k \le n} \frac{6}{k+1} - 2\sum_{1\le k \le n} \frac{1}{k} \right) \\ F[n] &= \frac{n(n+1)}{2n^2} \left( \sum_{1\le k \le n} \frac{6}{k+1} - 2H[n] \right)\\ \end{align*}$$ O primeiro somatório precisa de um pouco de massagem: $$\begin{align*} \sum_{1\le k \le n} \frac{6}{k+1} &= \sum_{1\le q-1 \le n} \frac{6}{q} \\ &= \sum_{2\le q \le n+1} \frac{6}{q} \\ &= -6+\frac{6}{n+1}+6\sum_{1\le q \le n} \frac{1}{q} \\ &= -6+\frac{6}{n+1}+6H[n] \\ \end{align*}$$ Por fim: $$\begin{align*} F[n] &= \frac{n(n+1)}{2n^2} \left( \sum_{1\le k \le n} \frac{6}{k+1} - 2H[n] \right)\\ F[n] &= \frac{n(n+1)}{2n^2} \left( -6+\frac{6}{n+1}+6H[n] - 2H[n] \right)\\ F[n] &= -3 +\frac{2(n+1)}{n}H[n]\\ \end{align*}$$ Ufa, deu trabalho! E isso é ##O(\log n)## mesmo? Bem, podemos verificar usando um assintótico para o harmônico, por exemplo: $$ H[n] = \ln n + \gamma + O(1/n) $$ Agora você substitui: $$\begin{align*} F[n] &= -3 +\frac{2(n+1)}{n}H[n]\\ &= -3 +2H[n] +\frac{H[n]}{n} \\ &= -3 +2\ln n +2\gamma +O(1/n)+ 2\frac{\ln n}{n}+2\frac{\gamma}{n}+O(1/n^2) \\ &= 2\ln n + O\left(\frac{\ln n}{n}\right) \\ &= O(\log n) \end{align*}$$ Como queríamos demonstrar!

Ou seja, realmente a Busca por Bozo Binário é ##O(\log n)##. Em demonstrações grandes assim eu sempre gosto de fazer uma simulação por Monte Carlo para conferir as contas. A minha simulação está no github:

Busca por Bozo Binário, simulação por Monte Carlo

Simulando 50000 vezes em um vetor de tamanho 3000, o número médio de comparações foi de ##14.1769##. O valor teórico previsto pela fórmula é de ##14.1732##, nada mau!

Eu também fiz uma simulação da busca binária tradicional como comparação. Com os mesmos parâmetros, a busca tradicional usou ##10.6356## comparações no caso médio, ou seja, foi mais eficiente que o Bozo Binário.

Por que isso acontece? Analisando o assintótico fica claro. Para ##n## grande, o valor médio da Busca por Bozo é de ##2 \ln n## (logaritmo natural), enquanto que na busca binária é ##\log_2 n## (logaritmo de base 2). Então, em média, esperamos que Bozo Binário seja mais lento por um fator de:

$$\begin{align*} \frac{2\ln n}{\log_2 n} &= {2\frac{\log n}{\log e}} \div {\frac{\log n}{\log 2}} \\ &= 2\frac{\log n}{\log e} \times \frac{\log 2}{\log n} \\ &= 2\frac{\log 2}{\log e} \\ &= 2\ln 2 \\ &\simeq 1.386 \end{align*}$$ Ou seja, mais ou menos uns 38% mais lento, o que bate com a simulação.

Depois da análise do algoritmo, a dúvida que resta é: e o candidato? Bem, se fosse eu a entrevistar, precisaria coletar mais dados para saber se ele realmente sabia do tempo médio, ou se estava perdido e implementou a primeira coisa que veio na cabeça.

Um possível follow-up é "Certo, esse método é mais complicado que a busca binária tradicional, mas funciona. Você pode me falar alguma situação onde ele é preferível sobre o tradicional?". Curiosamente, existe sim um caso onde a Busca por Bozo Binário é melhor que a busca binária tradicional! Mas esse caso fica para a próxima :)