Nos últimos dias parte da Twittosfera ficou quente devido a um polêmico estudo divulgado pelo economista Samy Dana, do qual ele é um dos autores. A principal conclusão do estudo é que seguindo o padrão de confinamento até então, o número de mortes no Brasil e São Paulo seria muito menor do que estudos anteriores, especialmente comparado ao estudo do Imperial College. Tanto é que a chamada diz “Novo estudo contesta previsões”. Este artigo faz uma análise pra ver se é isso mesmo.
No vídeo que ele divulgou num programa semanal, Samy afirma que o estado de São Paulo teria no total entre 5900 e 15600 mortes com mediana de 8557. O pico de óbitos diários seria entre os dias 29 de abril e 4 de maio com a mediana de óbitos diários no pico de 219 mortes.
Muitos demonstram ceticismo com relação aos resultados do artigo. Alguns, inclusive, deram pouca credibilidade aos resultados uma vez que não há nenhum especialista na área específica (epidemiologia). A maioria das críticas certamente não olhava os detalhes do método. Ainda que o Preprint esteja disponível, vamos lembrar que está escrito em inglês com um linguajar bem técnico. Há uma apresentação em PowerPoint disponível em português, porém ela não contém todos os detalhes.
Confesso de antemão que também estava cético dos resultados apresentados. Eu não tinha visto previsões de que o número de mortes cairia a zero já em junho deste ano. Também não sou especialista na área específica, mas tenho alguma experiência com modelagem matemática e em estimação Bayesiana ainda que a minha experiência neste último seja bem menor.
Entendendo o Modelo
Justiça seja feita, o Pre-print do artigo está acessível a todos. No entanto, em nenhuma das suas apresentações (vídeos) é claro o mecanismo que explicaria a queda abrupta no número de casos. De tanto ler sobre Covid19 nos últimos meses, minha pergunta sempre foi “o que faz no modelo deles ter o R menor do que 1?”. Esse R ao qual me refiro é o falado “número de reprodução” que é definido como o número médio de pessoas que uma pessoa infectada transmite a doença (há um vídeo da Angela Merkel que ela explica isto muito bem).
Um R de 1,2 quer dizer que cada pessoa infecta outras 1,2 pessoas. Esses números quebrados não ajudam uma vez que não infectamos 0,2 pessoas, mas digamos que numa cidade por aí existam 10 pessoas infectadas. Estas 10 pessoas vão levar a infecção de outras 12 pessoas. Estas 12 pessoas, infectarão entre 14 e 15 pessoas e assim sucessivamente. No décimo ciclo isto já seriam umas 70 pessoas. No vigésimo umas 5000. Deste jeito chegamos há mais de 4 milhões de infectados no mundo a partir de um pequeno grupo de pessoas em poucos meses. Essa Thread no Twitter discute um pouco o significado e os desafios de medir esse R.
Agora imagina que este R é de 0,8 e iniciamos com essas mesmas 10 pessoas infectadas. No ciclo seguinte 8 novas pessoas seriam infectadas. No segundo ciclo seriam 6,4 e por aí vai. Veja que se este R quando é menor que 1, o número de casos diminui a cada ciclo em vez de aumentar. Portanto, se o modelo prevê que o número de mortes cairá a zero isto significa necessariamente que o valor de R em algum momento tornou-se menor do que 1 (Para uma discussão sensacional, leia o texto do Marcel aqui).
Não sou epidemiologista e também não sei se estou usando esta notação corretamente (ficaria grato com sugestões), mas há três grandes fatores que influenciam o número de reprodução R. Primeiro são os nossos cuidados diários para evitar o contágio. Lembra quando no início o número de casos dobrava a cada 3 ou 4 dias e depois passou a dobrar muito mais lentamente? Isto é reflexo das medidas de higiene e distanciamento social que todos nós realizamos. Neste caso o R pode ter saído de algo como 3 para algum valor próximo de 1.
O segundo fator é a imunidade obtida seja se infectando ou por vacina. Digamos que em um local o valor do R é 1,25 com nenhuma pessoa infectada. Se 20% da população desenvolve imunidade, o R passa a ser efetivamente 1,25×0,8=1. A ideia é simples, basicamente temos 20% a menos de pessoas para infectar. Este mecanismo que explica a tão falada imunidade de rebanho. No caso do Covid, até onde eu sei, acredita-se que a imunidade de rebanho deve ser alcançada com algo em torno de 60% de pessoas infectadas.
Por último, o R também é influenciado pelo próprio vírus em si e dos sintomas das pessoas infectadas. Por exemplo, esse R seria menor caso infectados assintomáticos não transmitem o vírus. E o R seria maior se o período de incubação fosse maior.
Só para ter uma ideia inicial de qual tem sido o valor do R, eu estimei o valor do R (variando com o tempo…) com base no número de casos e mortes do estado de São Paulo. Deixo claro aqui que utilizei um método simplificado adaptado de uma sugestão de um epidemiologista no Twitter (link). Basicamente, como o ciclo do vírus é em torno de 5 dias, basta dividir o número de casos no dia t pelo número de casos no dia t-5. Como há bastante variação diária, principalmente no final de semana, eu fiz esta estimação com base em média móvel de uma semana (isto é, o valor que eu considerei do tempo t é sendo a média dos 7 dias anteriores).
Número de reprodução R estimado para o estado de São Paulo
Veja que o R começa substancialmente alto (o primeiro dia do gráfico corresponde ao dia 10 de março logo após a primeira morte). Ele posteriormente cai, há um pico novamente e depois volta a cair. Veja que, apesar de boa parte da população não estar saindo de casa, o valor de R raramente ficou abaixo de 1. Isto quer dizer que o número de casos tem crescido.
Esse valor de R é tanto um dado de entrada do modelo quanto um dado de saída. No estudo, eles utilizam o valor de R no passado para estimar o valor futuro de R. Eles realizam a estimação do R com um modelo na literatura. Aparentemente é um método de estimação Bayesiano e relativamente sofisticado até onde eu entendi. Para projetar os dados no futuro, eles usam esta estimação de R com base na seguinte função:
Aproximação do R com o tempo proposto pelos autores.
onde beta é obtido pelo método de calibração Bayesiano que utilizam, porém o detalhe importante é em como eles estimam gamma, delta, e lambda. Eles ajustam estes valores separadamente com base no R que eles mediram pelos dados históricos.
As equações podem parecer intimidadoras, mas isto quer dizer o seguinte: “vamos tentar estimar uma linha que consegue replicar de forma satisfatória o R como função do tempo t”. A premissa básica é que podemos identificar uma certa tendência, por exemplo, uma tendência de baixa e podemos utilizar este padrão para estimar o R no futuro.
Para dar uma ideia do conceito, eu vou usar os dados do PIB do Brasil para exemplificar. A linha preta na figura abaixo mostra a taxa anualizada do PIB no Brasil entre o início de 2014 e o fim de 2018. Imagina que no início de 2018 (marcada com uma linha vertical preta no gráfico) eu quisesse estimar o PIB do Brasil para o ano de 2018. Como eu não sou economista, eu não teria ideia por onde começar. Por que não tentar dar uma olhada na “cara” da curva de PIB nos últimos 4 anos (curva preta) para fazer a estimativa do PIB? Como essa curva preta tem duas concavidades, eu achei que uma função do terceiro grau seria boa para descrever essa curva.
Dito e feito, coloquei meu computador para calcular os parâmetros de uma função do terceiro grau que melhor se aproxime da reta preta (lembre-se, estou em 2018). O resultado desta curva é a linha contínua azul no gráfico. Veja que ela aproxima muito bem a linha preta até 2018! Mas lembre-se eu quero calcular o que acontecerá depois de 2018 e a linha azul só está definida entre 2014 e 2018. A solução é simples, eu tenho que apenas “continuar” a curva na mesma tendência, ou seja, “extrapolar” essa curva. Essa continuação da curva está plotada com uma linha tracejada azul. Se eu tivesse no início de 2018 eu acharia que o PIB iria crescer bastante naquele ano. Faria sentido, afinal o PIB foi de -5% para 1% em um ano por que não continuaria subindo? Aí o Felipe do fim de 2018 descobriu que a previsão estava muito errada. O PIB simplesmente parou de subir, ou seja, o PIB não manteve a tendência dos trimestres anteriores.
Exemplo de interpolação polinomial
Esta é a ideia de calcular linhas de tendências. Olhamos para o passado para identificar a tendência e depois projetamos para o futuro. Esta abordagem é muito comum em economia e ciências em geral e por vezes esta aproximação é ótima, porém tem que ser utilizada com cuidado para não ser enganado como no caso da Figura acima.
Pois bem, com base no R mostrado na primeira figura, eu estimei os parâmetros gamma, delta e lambda como eles propõem no artigo, com base no perfil de óbitos no estado de São Paulo. O resultado que eu cheguei está ilustrado na figura abaixo. A linha vermelha mostra o valor do R até onde temos dados. A linha azul continua mostra essa função que descreve a tendência do R como a equação acima. O valor da extrapolação (isto é, os valores previsto para dias futuros) é mostrado em linhas tracejadas. Veja que nos primeiros dias quando o R era alto (não estávamos implementando isolamento social) o modelo subestima os valores. Logo após a tendência se reverte. A grosso modo o R estimado é razoável até o ponto em que temos dados disponíveis.
Resultado da estimação do R com a melhor combinação de constantes para a curva que eles propõem.
Veja, no entanto, o que o modelo prevê para além do dia atual: uma queda consistente no valor do R. Basicamente isto quer dizer que o número de casos vai começar reduzir de forma drástica.
Só para me certificar que o que eu implementei está consistente com os resultados deles, abaixo eu ploto o número de mortes observados e o número de mortes previsto por este modelo. Se eu somo o número de mortes, chego a um valor de 8545 mortes que é bem próxima da mediana da previsão no artigo deles (8557). Pode haver uma diferença porque simplifiquei partes. Não estou dizendo que os números que reporto aqui são exatamente iguais aos que eles reportam no artigo deles. No entanto, qualquer diferença não parece ser significativa. Quem tiver interesse, pode me contactar no Twitter (@fas_souza) que eu envio todo o código para gerar todas as figuras. Utilizei dados do BrasilIO para isso baixando os dados em formato CSV para o Estado de São Paulo.
Estimação do número de mortes com base na estimação do R como na figura anterior.
O que eu descrevo acima é a base do modelo deles. Eles modelam uma série de outras coisas como utilização de leitos de UTI e número de hospitalizações. Não vi esses detalhes, mas pela apresentação em vídeo que eu vi me pareceu razoável.
Discussão e Crítica
Deixo claro a partir de agora que tudo que eu comentar é baseado somente na minha opinião. Podem discordar à vontade.
A minha crítica central é que eles não apresentam nenhuma justificativa lógica para esta interpolação do R. Como comentado, o R depende da capacidade do vírus de se transmitir, da política de isolamento social e do número de pessoas imunes. Já que não podemos controlar a transmissibilidade inerente ao vírus e considerando que eles assumem que a política de quarentena se mantém, só sobra confiar que há uma subestimação significativa (mais do que já conhecemos) no número de casos e que estaríamos próximos da imunidade de rebanho.
Faria sentido? Bom, se eles estimam 10 mil mortes e a população do Estado de São Paulo sendo de 44 milhões daria em torno de 0,02% de mortes no total. Isto me parece muito otimista se considerarmos que no estado de Nova York já morreu 0.14% da população (não é 0,14% da mortalidade do vírus, é 0,14% da população em geral!).
Outro argumento possível seria confiar que haveria queda nas transmissões devido a novas descobertas ou implementação de políticas eficientes de rastreamento e isolamento de casos suspeitos. Até faz algum sentido. Começamos com o R em torno de 2 e logo após a quarentena foi implementada, o R caiu para um valor próximo de 1,5. Depois passou-se a utilizar a máscara e isso provavelmente tenha contribuído para uma nova queda do valor de R. O padrão em geral faz sentido, porém os efeitos dessas melhorias tendem a ser menores ao longo do tempo. Não faz sentido considerar que o R pode ir a zero em tão pouco tempo. Portanto, acho improvável que isso aconteça e no todo me parece que a ideia de aproximar o R do jeito que eles propõem é o perfeito análogo de aproximar a curva de crescimento de PIB brasileiro com um polinômio de terceira ordem.
Mas vamos relevar isso por enquanto e digamos que a proposta de aproximar o R daquela forma seja plausível. Quem sou eu para dizer que não? Vamos olhar a predição de óbitos pelo modelo deles. A previsão é que basicamente não haverá óbitos já em junho, ou seja, algo em torno de 3 semanas a partir de agora. Com exceção de talvez Coréia do Sul e China, vocês viram algum país que conseguiu levar o número de casos a zero em 3 semanas? Veja os gráficos de países como Alemanha, Itália, Espanha e Estados Unidos e observe que a queda é sempre mais lenta do que o crescimento de casos. A Itália, por exemplo, estava no dia 20 de março com 6500 casos diários aproximadamente. Transcorridos dois meses de um confinamento agressivo, em meados de maio eles ainda reportam algo em torno de 1000 casos diários.
E aqui fica a crítica principal, se os resultados são incríveis, precisamos fornecer evidências incríveis (alô, Carl Sagan). Se estamos prevendo que a queda de casos comparada com a Itália será ao menos duas vezes mais rápida, teríamos que prover evidências adicionais da plausibilidade deste resultado. Por mais sofisticada que seja a análise e o método que utilizamos, temos que estar cético dos nosso resultados.
A minha crítica quanto ao estudo em si é essa. Acho os resultados pouco plausíveis ainda que que eu torça que eles estejam corretos, uma vez que significaria um número de mortes muito menor do que esperaríamos. Agora eu discuto dois argumentos que foram utilizados pelo Samy nos vídeos que ele apresentou o estudo.
Primeiro, eles argumentam que propuseram um modelo adaptado à realidade nacional uma vez que os dados no Brasil não são confiáveis e etc. Segundo ele, “ninguém tem modelo no Brasil, o pessoal pega os modelos importados e adapta”. Essa linha de argumento é motivo de discussões na academia. Por vezes, este argumento faz bastante sentido e em outros casos nem tanto. Via de regra, quanto mais “dura” é a ciência, menos esse argumento faz sentido. Por exemplo, a propagação da luz ocorre exatamente da mesma forma numa cidadezinha do interior de Santa Catarina e na Times Square em Nova York. Por outro lado, a dinâmica de desigualdade social é, provavelmente (não é minha área…), muito diferente no Brasil e na Noruega.
No entanto, há que ter um cuidado quando se afirma que um modelo não se aplica em determinado contexto. É o modelo que não se aplica ao Brasil ou os parâmetros do modelo que são diferentes? Por exemplo, a aceleração gravitacional é diferente na superfície da Terra e na superfície de Marte, porém a lei da gravitação dos corpos se aplica em ambos. Como sabemos, a aceleração da gravidade depende da massa e diâmetro do planeta. Os modelos importados que o Samy se refere são os SIR e SEIR que, da mesma forma que o modelo gravitacional, precisa de parâmetros de entrada. Quando o Samy diz que o modelo SIR e SEIR superestimam o número de mortos, é uma característica do modelo em si ou uma falha dos dados de entrada? Se é uma falha do modelo, o que exatamente faz o modelo superestimar o número de mortes? Em nenhum momento ele explica.
Eu até concordo que as estimativas para o Brasil talvez estivessem superestimadas, mas temos que considerar que: (i) parte destas estimativas foram realizadas antes do confinamento e portanto a velocidade de transmissão do vírus claramente diminuiu; e (ii) tínhamos pouquíssimas informações sobre a doença em si e a letalidade do vírus. Na época, acreditava-se que a mortalidade era em torno de 3% (baseado na mortalidade de China e Coréia do Sul), mas estudos recentes apontam que a taxa de mortalidade está provavelmente entre 0,6% e 0,9% dependendo das características demográficas de cada país ou região. Se ele não ajustou as estimativas do modelo com base nesses dois aspectos, ele não pode afirmar que o modelo importado não se adapta a realidade nacional.
A segunda crítica que eu teria seria em focar demais no método (Hamiltonian Monte-Carlo) como se fosse a peça mais importante do quebra-cabeça e respondesse todas as questões. Métodos de Monte-Carlo têm inúmeros usos, mas via de regra esta metodologia é utilizada para refinar os resultados de um modelo existente que provavelmente já seria útil sem a utilização deste método. Neste caso específico, utilizar Monte-Carlo permite prover uma faixa de incertezas nas previsões e não apenas prever um único valor. Por exemplo, eles afirmam que a mediana das previsões é de 8557 mortes, porém podem estar entre 6400 e 14400. Isto com certeza é interessante. Porém o fato é que sem utilizar esta metodologia eu consegui replicar a grosso a mediana da previsão deles. Portanto, não vejo que este seja um caso que o uso de métodos de Monte Carlo é fundamental. Por outro lado, isto não quer dizer que esta metodologia não possa ser usada ou que os resultados não são plausíveis por causa disso.
Utilizar o “Hamiltonian Monte Carlo” (ou qualquer algoritmo com nome difícil e que parece ser sofisticado) com base em um modelo ruim e esperar ótimos resultados não faz muito sentido. Seria a mesma coisa que colocar rodas mais largas e trocar o carburador de um carro da década de 60 e esperar que ele seja mais estável e eficiente que um carro produzido nesta década. Por outro lado, este método associado a um modelo adequado pode aprofundar o entendimento do problema e as incertezas do modelo. A utilização desta técnica seria como ajustar cada peça de um carro para uma corrida específica para ganhar 0,5s por volta, como na Fórmula 1, é uma ótima ideia. Meio segundo por volta é pouco em uma volta de 90 segundos, mas pode ser o detalhe que define o vencedor da corrida.
Agora vou a um ponto polêmico: O Samy pode publicar seus resultados mesmo não sendo especialista na área?
Claro que pode e deve ser incentivado a fazer isso. Porém como qualquer coisa que fazemos pela primeira vez, temos que ter cuidado quadruplicado. A analogia básica seria dirigir em uma rodovia movimentada após tirar a carteira de motorista. Se temos carteira de motorista, claro que podemos dirigir. Porém nas primeiras vezes sempre deixamos mais “folga” nas manobras e dirigimos um pouco mais devagar. Uma vez que temos maior noção de velocidade, distância e do comportamento dos outros motoristas, realizamos as manobras com uma “folga” menor – ainda que sempre tomando cuidado!
A minha maior crítica nesse caso é que ele não tomou nem o cuidado do motorista experiente e muito menos do motorista novato. Ele chegou a um resultado que prevê número de mortes ZERO por junho (ao menos é isso que eu vejo no gráfico no preprint e na primeira apresentação dele). Sendo um resultado que é bem diferente do que tem ocorrido em outros países e sem ter uma explicação sólida do mecanismo que explicaria a queda de casos, divulgar esses resultados seria como realizar uma ultrapassagem perigosa.
Por último, baseado neste último parágrafo nem preciso falar o quão temerário é utilizar resultados que não parecem plausíveis para orientar políticas públicas. Nos vídeos o Samy fala que contactou governadores e etc. Ele pode fazer isso? Sim. Deveria? Na minha modesta opinião não. Ou só faria depois de, no mínimo, do artigo passar por revisão por pares e depois de mostrar os resultados para muitos epidemiologistas (Para ler mais sobre revisão por pares, tem esse texto aqui). Caso minha análise esteja certa (não estou afirmando com 100% de certeza do que ocorrerá) e a queda não seja como o modelo proposto por eles, vários estados poderiam decidir não aumentar o número de leitos de UTI, por exemplo. Acho os resultados deles muito preliminar para basear decisões de políticas públicas. Isto obviamente não quer dizer que os resultados de outros modelos não podem ter erros. Portanto, o fato de eu fazer uma crítica do modelo não quer dizer que eu tudo que o Imperial College publica devemos aceitar sem criticar também.
Muito da minha discussão acima foi influenciada diretamente ou indiretamente por essas 4 arrobas no Twitter que são sensacionais para acompanhar uma decisão técnica sobre o assunto:
@ Nataliexdean
Por último, fiz algumas críticas pontuais do que foi falado em relação a divulgação deste artigo relacionado ao Samy Dana – por ser a pessoa de maior exposição. As críticas específicas a forma como foi divulgado portanto não deve se aplicar aos outros autores. Além disso, nenhum dos autores, incluindo o Samy, afirma que a quarentena deve acabar e portanto não estou acusando de forma alguma que o artigo seja apenas um pretexto para tal.
PS: Agradeço ao Alexandre Simas, um dos coautores do artigo. Ele foi solicito e confirmou que a tendência do R apresentada aqui bate, a grosso modo, com os resultados dele. Diga-se, todos no Twitter tem confirmado que eles tem sido solícitos quanto a isso.
Referências
- Preprint do artigo deles
- Chamada da apresentação do modelo deles que contestam estimativam anteriores. Este link contém uma apresentação em slides em português e para um dos vídeos ao qual me refiro.
- Chamada para a segunda apresentação quando afirmam que estamos no pico e que o modelo é com características brasileiras.
- Caso queiram ter acesso ao código para gerar estas figuras, podem utilizar os comentários ou me contactar no Twitter (@fas_souza)
Agradecimentos
Alguns amigos deram uma lida neste texto antes e deram boas sugestões, especialmente o Eduardo Hulse. Minha prima Patrícia e irmã Mariana leram e não entenderam boa parte da discussão sobre interpolação o que me motivou a simplificar esta parte. Espero que esteja mais fácil de entender!
