BayArea

BayArea (chamado por alguns autores de Bayes Area para evitar confusão na busca de palavras-chave, já que Bay Area é o nome da área costeira da cidade de São Francisco, CA) é um método proposto por Landis et al. (2013) que visa inferir cenários acerca da biogeografia histórica de um grupo de interesse. Segundo os autores, biogeografia histórica é “estudo da distribuição passada das espécies e dos processos que a influenciam”. Devido a problemas de tradução e para simplificação, deixemos claros alguns conceitos e termos usados nesse texto:

• Área: área discreta que descreve apenas os locais em que os terminais da árvore estão presentes (em inglês: area).
• Área de distribuição: conjunto de uma ou mais áreas discretas que visa modelar uma área contínua, i.e., o padrão de presença e ausência através do conjunto de áreas (em inglês: range).
• BayArea: foco desse artigo. Quando falamos em BayArea, estamos nos referindo ao modelo proposto por Landis et al. (2013).
• Pré-DEC: o BayArea é baseado no método proposto em 2005 por Ree et al. – que mais tarde daria origem ao método conhecido como DEC (dispersal-extinction-cladogenesis) (Ree & Smith 2008). Por isso nos referiremos a ele como “pré-DEC”.

Figura 1: Num Processo de Markov, há a mudança de estado de uma variável através do tempo de forma que a cada tempo t a variável tem um valor xt (a). No entanto, para a seleção de cada mudança de estado, o passado é ignorado (b) e o futuro é decidido com base unicamente no valor presente da variável (c).

O pré-DEC é um método baseado em verossimilhança que se utiliza de processos de Markov de tempo contínuo. O BayArea, além de modificar algumas premissas e acrescentar parâmetros, trabalha dentro da inferência bayesiana. Para entender tanto o pré-DEC quanto o BayArea, precisamos antes visitar alguns dos conceitos matemáticos e estatísticos envolvidos com essa inferência.

Conceitos matemáticos envolvidos com o método

O primeiro conceito é o de Processo de Markov. Para quem já trabalho com estatística Bayesiana, o nome Markov deve soar familiar. Um processo de Markov é um processo do tipo “sem-memória”, isto é, ele desconsidera o passado ao fazer decisões acerca do futuro. Digamos que estamos tentando encontrar o valor mais provável para uma determinada variável. Numa busca que se encaixe no conceito do processo de Markov, esse valor muda de estado de tal forma que o próximo estado dela depende apenas do estado presente dela. A figura 1 mostra uma representação disso: existe uma trilha formada por uma série de estados que levou a variável a chegar ao estado x (a); num processo de Markov, essa trilha é ignorada para ir do estado x ao estado x + 1 (b) e apenas o estado presente x é considerado (c). Um jeito bonitinho de colocar essa premissa é que o passado e o futuro, enquanto condicionais ao presente, são independentes.

Uma série de estados obtidas de maneira obediente ao processo de Markov (isto é, a “trilha” da fig. 1) é denominada cadeia de Markov. Na inferência bayesiana de programas como o MrBayes, a cadeia de Markov é obtida através de um método chamado “Monte Carlo random walk” (o velho conhecido MCMC). O método de Monte Carlo é uma “caminhada aleatória” e isso quer dizer que ela é medida em passos discretos. Por isso, esse método obedece a um processo de Markov de tempo discreto. No caso do processo de Markov de tempo contínuo, as coisas ficam um pouco mais complicadas e o BayArea usa as duas coisas em momentos diferentes. (Uma série de vídeos bastante fáceis para entender o processo de Markov está disponível em: https://www.youtube.com/playlist?list=PLX2wbP1VmGem7E2T_fC6co9oxQ1SF3FAC)

Imagine que você é um gerente de banco num universo extremamente fictício em que os bancos querem realmente melhorar o atendimento. Nesse universo fictício, filas não existem porque quando as pessoas chegam ao banco e não podem ser atendidas imediatamente, elas vão embora. Se você quer melhorar o atendimento, o início do expediente não importa para você, pois assim que o banco abre todos os caixas estão livres; estamos interessados nos horários de pico, em que os caixas não conseguem ficar disponíveis por muito tempo.

Figura 2: No exemplo do banco (ver texto), temos um número x de caixas cada um capaz de atender somente um cliente por vez. Quando um cliente chega e encontra um caixa vazio, ele se dirige imediatamente a ele. No entanto, na falta de caixas disponíveis, o cliente dá meia-volta e sai do banco sem ser atendido.

Estabelecido esse cenário (ver Fig. 2), vamos estipular alguns parâmetros e variáveis. Cada caixa tem atendentes diferentes, resolvendo problemas diferentes para pessoas diferentes. Cada um deles vai demorar um tempo τ que podem ser diferentes entre os caixas, mas também podem ser iguais. Dessa forma, se fôssemos oniscientes, veríamos um pequeno reloginho despertador que dispara ao fim do atendimento. Entretanto, existe um outro relógio com que devemos nos preocupar: os clientes também chegam a uma certa taxa e, por isso, é possível que imaginemos um reloginho à porta do banco também. Nossa preocupação é o estado x de caixas disponíveis em que o banco se encontra num tempo t, pois, toda vez que o relógio de chegada de clientes dispara, é preciso ter caixas disponíveis, senão o cliente simplesmente vai embora e, toda vez que o relógio de um caixa dispara, é preciso que clientes cheguem para que ele não fique ocioso (Fig. 3).

Figura 3: Baseado no exemplo do banco, colocamos agora relógios que medem o tempo probabilístico de cada atendimento além de um outro relógio que mede o tempo probabilístico de chegada de clientes. Quando o relógio do cliente dispara, um cliente chega no banco e segue a regra normal do cliente: ir embora na ausência de caixas disponíveis (cenário 1 na animação). Quando o relógio do caixa dispara, o atendimento termina e o cliente sai do banco. Assim o caixa fica disponível para caso o relógio do cliente dispare novamente (cenário 2 na animação). (Lembre-se de que os relógios são reiniciados a cada disparo, então, um outro cenário possível é o disparo do relógio de outro caixa e o consequente aumento do número de caixas disponíveis.)

Existem duas grandes sacadas desse modelo matemático. A primeira é que a taxa da probabilidade de salto (“reloginhos”) não é fixa. Ela é uma curva de probabilidade de acordo com os seus dados – mais especificamente, ela é dada pela exponenciação de uma matriz que você insere no método. A segunda grande sacada jaz no fato de que esse é um processo de Markov, e isso diz muito sobre como esses reloginhos funcionam: em um banco real, o tempo é o mesmo para todos os caixas, de forma que se o caixa 1 leva 10 minutos para ficar livre e o caixa 2 leva 20 minutos, quando o relógio do caixa 1 dispara, faltarão apenas 10 minutos para o relógio do caixa 2 disparar. O estado do caixa 2 se torna correlacionado ao estado do caixa 1 e isso não condiz com um processo de Markov – lembre-se: o estado futuro de uma variável só depende do estado presente dela. Para manter a premissa de Markov, a cada disparo do relógio, todos são reiniciados.

Existem vários processos de tempo de contínuo com os quais você pode já estar acostumado: birth-death e processos de Poisson fazem parte desse tipo de algoritmo. Agora que já estamos familiarizados com a matemática da coisa, vamos ver como isso funciona para o BayArea.

Sobre o método

Sobre o pré-DEC

Um dos principais objetivos do pré-DEC é corrigir as falhas do método mais avançado da época, o DIVA (Ronquist 1996, 1997). A primeira falha é o chamado paradigma área-cladograma: a biogeografia histórica visava obter árvores para a história das áreas a partir das árvores representantes da história dos organismos. Esse conceito levou à aplicação de métodos que se utilizavam do conceito de “árvore-dentro-de-árvore” (Ree et al. 2005) que tentavam encaixar as árvores de organismos dentro das árvores de área baseados em critérios de otimalidade. O DIVA surgiu da tentativa de se desvencilhar desse problema, mas a história de áreas e de organismos permaneceram conceitualmente ligadas. O segundo problema é o fato de o DIVA explicitamente favorecer eventos de vicariância.

Ree et al. (2005) argumentam, então, que a abordagem do DIVA só pode ser útil quando a única informação disponível é a relação filogenética entre os organismos de interesse e propõem um novo modelo para inferência da história biogeográfica. O modelo do pré-DEC é composto por dois componentes teórico-práticos principais: um paleogeográfico (que consiste de um conjunto de áreas com conexões (rotas de dispersão) entre pares de áreas) e um biológico (que consiste em uma filogenia com os comprimentos de ramos calibrados e parâmetros que ditam as taxas de dispersão e extinção local). Portanto, este é um método paramétrico, i.e., depende dos parâmetros dos dados que lhe servem de input, e separa consistentemente a evolução das áreas da evolução das linhagens.

A partir desse modelo, ele propõe dois níveis para a evolução da distribuição dentro das linhagens: as dispersões e extinções locais ocorrem em nível anagenético, ou seja, ao longo dos ramos da árvore; em nível cladogenético, as áreas podem apenas ser herdadas pelas espécies filhas. Para modelar a evolução anagenética, o pré-DEC se utiliza de uma função de Poisson (que é uma função de tempo contínuo) com a exponenciação da matriz de uma série de valores que o usuário insere para informar sobre a disponibilidade de rotas de dispersão para uma determinada área (1 para rotas bem definidas e/ou áreas conhecidamente próximas; 0 para áreas muito distantes ou com barreiras intransponíveis entre elas; valores intermediários são aceitos). Simplificando para o modelo do banco: existe um relógio para a probabilidade de dispersão, um relógio para a extinção local e um relógio para nenhum evento; o primeiro que disparar é o mais provável de acontecer. Para modelar a evolução cladogenética, o método enumera cenários de herança da área de distribuição da espécie ancestral para as espécies-filhas e usa do modelo de poda de Felsenstein (1981) o mais provável, isto é, é um modelo baseado em máxima verossimilhança. Discutiremos as premissas desse método junto com a apresentação do BayArea em si.

Sobre o modelo de Landis et al. (2013)

Landis et al. (2013) introduz o problema do método de Ree et al. (2005) falando das limitações computacionais. As dimensões matriz utilizada pelo pré-DEC é exponencialmente proporcional ao número de áreas utilizadas: o número de estados (áreas de distribuição) possíveis é dado por 2ⁿ – 1 porque, além dos locais de origem de cada terminal, são contabilizadas combinações entre esses locais para que seja possível considerar extinção local ou expansão da distribuição por dispersão. Ou seja, para 10 áreas, existem 1023 estados possíveis. Esse número de estados para uma única característica (distribuição) acaba requerendo um poder computacional muito grande*.

Os autores do BayArea ressaltam ainda que a aproximação matemática feita pelo pré-DEC permite reconhecer todas as possibilidades de colonização de novas áreas e extinção local e que a primeira motivação de reescrever o método foi a necessidade de aumentar o número de áreas suportado pelo algoritmo. A primeira mudança que o BayArea apresenta é a implementação de um modelo bayesiano no lugar do modelo de máxima verossimilhança de Ree et al. (2005). Esse modelo bayesiano foi inspirado pelo método de análise de sequências de aminoácidos introduzido por Robinson et al. (2003). Em vez de integrar todos os estados possíveis (integração analítica), os autores implementam uma Cadeia de Monte Carlo Markov (integração numérica) para descrever as probabilidades finais da história biogeográfica e de cada evento.

A segunda mudança que Landis et al. (2013) propõe é a reestruturação dos modelos de evolução anagenética e cladogenética. Primeiramente, o BayArea requer coordenadas geográficas para as áreas inseridas. Isso porque a matriz utilizada para calcular a probabilidade de dispersão vs. extinção não é mais apenas uma matriz de probabilidades definidas a priori, mas uma matriz de probabilidades dependentes da distância entre esses pontos: pontos mais próximos têm maior probabilidade de ser conectados por rotas de dispersão do que pontos mais distantes (os relógios da dispersão que se referem a áreas mais próximas disparam mais rápido do que os que se referem a áreas distantes). Dessa forma, a evolução das áreas é descrita por dois modelos: um modelo nulo em que todas as áreas têm igual taxa de colonização (dispersão) e extinção e um segundo modelo em que a colonização é dependente da distância.

Dessa forma, a probabilidade de uma determinada história biogeográfica é definida pelo produto:

*- Existem modos de reduzir o número de estados possível para o número de áreas e assim aumentar o número de áreas, mas eles envolvem diminuir os arranjos e configurações entre as áreas e devem ser usados com cautela, baseando-se em noções empíricas

Premissas

Premissas do pré-DEC mantidas no BayArea

• Cada par de áreas possui uma rota de dispersão que pode ser descrita por uma função de probabilidade.
• A árvore inserida deve ser calibrada e tanto a sua topologia quanto sua datação serão consideradas corretas.
• Exclusão de hipóteses de vicariância ad hoc (cenários de vicariância são permitidos, mas não são nem favorecidos a priori nem explicitamente associados com os eventos das áreas).
• Quebra do paradigma área-cladograma: ao contrário de uma simples análise de caráter, em que a especiação gera táxons inicialmente com os mesmos caracteres, há a enumeração dos possíveis cenários de subdivisão e herança da área de distribuição ancestral, incluindo o cenário em que as espécies-filhas não herdam a mesma distribuição.
• Cladogênese é um evento instantâneo e independente de dispersão e/ou extinção local.
• Anagênese descreve o tempo de duração de uma linhagem e, portanto, o tempo disponível para o aparecimento de eventos de dispersão e extinção local.
• Áreas são unidades geográficas discretas que mantêm sua identidade dentro do tempo de interesse.

Premissas apresentadas pelo BayArea:

• As rotas de dispersão são descritas com base em sua distância geográfica partindo-se de um modelo nulo em que todas as áreas têm igual taxa de extinção e colonização.
• Substituição da integração analítica (i.e., contabilizando todos os estados possíveis) por uma integração numérica (baseada em MCMC): isso permite uma maior resolução pela possibilidade do aumento do número de áreas, mas leva a (1) uma maior sensibilidade aos dados geográficos inseridos e (2) à necessidade dos cuidados que se têm normalmente com uma inferência bayesiana qualquer (por exemplo, observar convergência e fazer múltiplas corridas para avaliar a concordância entre elas), pois esse tipo de análise está sujeito às ilhas de probabilidade.**

**- Landis et al. (2013), propõem uma extensão de seu método que permitiria a construção de novas árvores através da integração da verossimilhança marginal. Uma aproximação dessa abordagem é a inserção de várias árvores ao invés de uma só. Trataremos dessa possibilidade abaixo.

Tutorial

Em 2014 (Matzke 2014) e em 2015 (Yu et al. 2015) houve uma tentativa de reunir as alternativas entre os métodos de inferência biogeográfica em um único programa ou script, de forma a facilitar a comparação entre elas. A abordagem de Nicholas Matzke (2014) é interessante por permitir um teste de modelo, dar valores comparáveis aos resultados obtidos e facilitar a escolha do usuário. No entanto, essa abordagem, feita no R através do pacote BioGeoBEARS, é bastante primária em termos de implementação: por exemplo, o autor retirou a possibilidade de inserir coordenadas geográficas no BayArea e, portanto, obriga o usuário a retirar o fator distância da análise (ver discussão sobre isso no artigo de Matzke).

A abordagem de Yu et al. (2015) deu origem ao software RASP. O mais interessante do RASP é a possibilidade de inserir conjuntos de árvores e, logo, reduzir a incerteza contida na filogenia do input. O software possui algumas limitações computacionais, mas, para a maioria dos casos, essas limitações não devem ser tão críticas. Abaixo encontram-se tutoriais para o BioGeoBEARS (1) e para o RASP (2) (em inglês). (O tutorial para o RASP é incluído no próprio download do programa.)

1. http://phylo.wikidot.com/biogeobears
2. http://www.ib.usp.br/~silvionihei/pdf/A%20Rough%20Guide%20to%20RASP.pdf