Data Mining, Python

Quatro Módulos Essenciais para Desenvolver Machine Learning em Python

set 14, 2015 por Mario Filho

Ter as ferramentas corretas é essencial para o sucesso de qualquer projeto de machine learning. A maior parte do tempo gasto em um projeto que envolva extrair valor de dados é preparar os mesmos para análise. Neste artigo recomendo módulos em Python que são essenciais para cientistas de dados atingirem uma produtividade maior.

Antes de começarmos: PIP

https://pip.pypa.io/en/latest/installing.html

Esta ferramenta permite que você instale e atualize módulos Python através da linha de comando. Isto facilita bastante, pois evita dores de cabeça e perda de tempo na instalação destes módulos.

Este é um gerenciador de pacotes como Homebrew e Apt-get, com uma extensa biblioteca de módulos Python. A instalação dos módulos, inclusive de quaisquer dependências necessárias, acaba se tornando muito simples, o que nos permite concentrar em utilizá-los.

Todos os módulos do artigo estão disponíveis no repositório do PIP.

Numpy/Scipy

http://www.numpy.org/
http://www.scipy.org/

Estes módulos são indispensáveis para qualquer usuário de Python que trabalhe com projetos que exijam computação científica. Além de facilitar o armazenamento dos dados em estruturas (como arrays) que possibilitam fazer várias operações matemáticas, eles implementam diversas funções úteis de estatística, álgebra linear e computação numérica.

Melhor ainda, por possuir seus fundamentos em operações de álgebra linear, eles utilizam bibliotecas otimizadas em C por baixo do código Python, o que faz uma diferença substancial na performance, tornando as computações mais rápidas do que em Python puro.

Pandas

http://pandas.pydata.org/

Este módulo oferece estruturas e funções para facilitar a análise de dados em Python. Quem já usou R vai perceber algumas similaridades como o “data.frame”. Originalmente suas funções eram concentradas para análise de séries temporais, mas com o tempo foi evoluindo e hoje possui suporte para diversos tipos de dados.

Outra grande vantagem é que boa parte do código é implementada em C para melhorar a performance. Ou seja, temos a facilidade da linguagem Python com a performance similar à da linguagem C.

Scikit-learn

http://scikit-learn.org/stable/

Essencial para qualquer cientista de dados que use Python e Machine Learning. É o módulo de Machine Learning mais completo disponível na linguagem. Possui algoritmos para aprendizado supervisionado, semi-supervisionado, e não supervisionado.

A documentação é bastante clara, com diversos exemplos. Possui funções para facilitar a transformação dos dados, além de várias ferramentas para ajudar na validação e avaliação do modelo. Sua interface é simples, e em poucas linhas de código é possível treinar modelos bastante complexos.

Mais uma vez temos um módulo que utiliza código em C para tornar as operações mais rápidas.

Matplotlib

http://matplotlib.org/

Este módulo contém ferramentas para criar gráficos em 2d. É possivelmente o módulo de visualização mais utilizado em Python, e funciona muito bem com Numpy e Scipy. Existem muitos tipos de gráficos, adequados à maioria dos dados.

A única desvantagem é que, apesar de ser possível fazer boas visualizações com ele, falta um cuidado maior com a estética. É bastante útil para visualizações claras, mas deixa a desejar em ferramentas para torná-las mais bonitas. Ele é mais indicado para criar visualizações para trabalhos acadêmicos, ou durante o processo de exploração de um banco de dados, antes da análise.

Bônus: Keras

https://github.com/fchollet/keras

O scikit-learn possui uma diversidade grande de modelos, mas não disponibiliza redes neurais. Apesar de ser um modelo bastante conhecido, suas características acabam exigindo uma implementação própria para agregar toda a complexidade e variedade.

O Keras é um módulo que descobri recentemente que facilita demais o processo de criação e treino de uma rede neural. Ele disponibiliza vários algoritmos de treino, tipos de layers, regularização, e modos de operação, como redes recorrentes ou sequenciais.

Redes neurais são complexas por natureza, e ter um módulo com uma interface simples ajuda bastante na hora de utilizá-las.

Ele é construído sobre o módulo Theano, o que possibilita personalizar as funções otimizadas pela rede neural, além de permitir o treino da mesma em uma GPU.

Classificação, Kaggle, Natural Language Processing, Regressão

Determinando a Relevância de Produtos num Sistema de Buscas Usando Machine Learning

ago 24, 2015 por Mario Filho

Uma das tarefas mais importantes para um site de comércio eletrônico é garantir que o usuário encontre o produto que está buscando. Isto inclui um sistema de buscas que retorne produtos relevantes. Nesta competição a CrowdFlower, uma plataforma especializada na coleta de dados, disponibilizou dados sobre buscas por produtos, e a mediana da pontuação de relevância dada por três pessoas para cada um.

O objetivo era criar um modelo que pudesse calcular a relevância de um produto para uma determinada frase de busca.

Dados

Foram disponibilizados cerca de 10 mil registros para treino, e 20 mil para teste. Nestes registros havia a frase de busca, o título e a descrição do produto, a mediana da relevância e a variância da relevância atribuída pelos usuários.

Nos exemplos de teste havia uma parcela de registros “falsos”, que não eram utilizados para calcular a performance, mas foram colocados para evitar classificação humana.

Métrica

A métrica de avaliação foi o Quadratic Weighted Kappa. Ela mede a discordância entre dois “avaliadores”. Neste caso, media a discordância entre a pontuação dada pelo modelo e a pontuação dada por humanos. Ela é sensível não apenas à precisão, mas também à distribuição da pontuação. Em tese, uma atribuição de pontos aleatória resultaria num Kappa de 0, e uma concordância total, em Kappa igual a 1. Existem casos em que o Kappa se torna negativo.

Processamento dos Dados

Texto

A transformação básica do texto do título e da descrição é fazer uma bag of words. Em sua forma mais simples é uma matriz cujas linhas são os exemplos (documentos) e as colunas são palavras (termos). O valor de cada célula é o número de vezes que a palavra aparece naquele documento.

Uma transformação mais avançada é pegar esta matriz e fazer algumas operações utilizando a frequência dos termos em cada documento, e também no corpo de documentos. Esta é chamada de TF-IDF. Neste caso, os termos com maior presença dentro de um documento específico recebem um peso maior, mas que é descontado por um fator inversamente proporcional à presença dele em outros documentos.

Simplificando: palavras raras que aparecem bastante em um documento recebem maior peso, e palavras frequentes, que aparecem em todos os documentos, se tornam menos importantes.

Fora isso, é importante determinar um limite mínimo de vezes que a palavra deve aparecer no corpo de documentos, para diminuir o ruído.

Neste caso, a melhor representação que encontrei foi treinar uma matriz baseada apenas nos títulos, e outra apenas na descrição. Alguns participantes sugeriram treinar nos dois ao mesmo tempo, mas como se trata de um problema de busca, acredito que seja importante diferenciá-los.

SVD – Latent Semantic Analysis

Como você já deve imaginar, a matriz de documentos fica enorme, e é bastante esparsa. Isso dificulta o uso com alguns modelos não lineares. Neste caso existe a opção de reduzir a dimensionalidade usando o SVD. Este algoritmo de álgebra linear basicamente encontra os componentes principais da matriz de documentos, onde há maior variação. Isso reduz a dimensionalidade do problema, mas pode causar perda de informação.

Nestes casos de NLP eles acabam funcionando bem, porque filtram um pouco o ruído dos dados. Aplicar o SVD à matriz TF-IDF te dá uma representação conhecida como Latent Semantic Analysis. A interpretação é que cada componente obtido se refere a um “conceito” representado pelos documentos. Ele também é bastante usado para visualizar os documentos num espaço com menores dimensões.

Na competição utilizei SVD para reduzir a dimensionalidade da matriz TF-IDF antes de alimentá-la aos modelos.

Outros atributos

Além dos atributos acima, criei outros que pareciam fazer sentido, e alguns que encontrei em trabalhos acadêmicos relacionados. O mais importante deles foi a porcentagem de palavras da busca que estavam no título. Além disso, criei atributos baseados na similaridade média entre um produto e os outros da mesma frase de busca, e também a contagem de vezes que a frase de busca aparece nos dados.

Uma parte importante foi codificar cada frase de busca com one hot. Ou seja, para cada uma delas criei uma coluna que possuía o número 1 caso ela estivesse presente, e o número 0 caso contrário. Fiz isso para ajudar o modelo a capturar variações próprias de cada frase de busca.

Regressão ou Classificação?

Existiam duas maneiras de encarar essa tarefa: considerar cada nível da mediana como uma classe, sendo assim um problema de classificação, ou considerar como um problema de regressão e arredondar as previsões.

Num primeiro momento achei mais adequado tratar como regressão, mas como o Kaggle queria as previsões arredondadas, acabei criando vários modelos de classificação e não tive tempo de testar a outra alternativa.

Com a classificação temos duas opções: usar a classe prevista pelo modelo, ou utilizar as probabilidades de cada classe para fazer uma média ponderada para obter a pontuação. A primeira alternativa se mostrou melhor neste caso.

Modelos

SVM

Os melhores modelos foram obtidos usando SVM. Ele é um dos modelos mais utilizados para classificação de texto, então não é surpresa que tenha boa performance.

O melhor modelo individual, utilizando os componentes do SVD e os outros atributos, atingia um Kappa de 0,672 na validação cruzada, e 0,654 na LB.

Random Forest

Um dos meus modelos preferidos, mas normalmente não são a melhor opção para dados esparsos em altas dimensões. De qualquer maneira, após vários testes, consegui 0,6695 na validação cruzada, e 0,6508 na LB usando uma matriz TF-IDF junto com os outros atributos.

Uma razão que pode ajudar a explicar o fato dela se comportar melhor com a matriz TF-IDF do que com a matriz do SVD é que este é um método que reduz a variância, então acaba sendo robusto diante de muitos atributos, desde que uma parcela deles seja relevante para a previsão.

Gradient Boosted Trees

Apesar deste modelo também não ser muito recomendado para uso com dados esparsos em altas dimensões, decidi treinar vários modelos “fracos” e fazer o stacking usando as previsões deles e os outros atributos, que não envolviam a representação das palavras.

Os modelos do stacking tinham uma média de 0,60 Kappa na validação.

A pontuação deste modelo ficava por volta de 0,624 na LB. Mas a intenção real com ele era ter um modelo diferente do SVM para fazer o ensemble.

Pós Processamento

Em todas as soluções de Top 10 que eu li, havia um detalhe em comum: o pós processamento das previsões. As equipes usaram métodos especiais para fazer arredondamento, multiplicando por coeficientes, ou determinando pontos de corte diferentes para a previsão de cada “classe”.

Em um dos casos, o time ordenou as previsões e, de acordo com a proporção de cada classe nos dados de treino, arredondou as previsões. Kappa é uma métrica preocupada com a distribuição das previsões, então isso acabou ajudando.

O mais próximo que cheguei de um ajuste nas previsões foi utilizar a média geométrica em vez da média simples para o ensemble. Os modelos costumavam superestimar a relevância dos produtos, e ela puxava as previsões para a classes mais baixas.

Ensemble e Overfitting

No fim, fiz uma média geométrica com meus melhores modelos. O ensemble era formado por vários SVMs, uma Random Forest e um GBT.

Apesar de não ser uma prática recomendável, principalmente em casos com uma quantidade de dados pequena como esta competição, utilizei a LB para validar o ensemble. Isso acabou fazendo com que eu escolhesse o meu segundo melhor ensemble, mesmo assim acabei no Top 5% da competição.

Classificação, Kaggle, Python

Como Consegui o 6º lugar Entre 985 Participantes na Competição de Machine Learning do Facebook

ago 3, 2015 por Mario Filho

Nesta competição de recrutamento, o Facebook disponibilizou dados anônimos sobre usuários de um site de leilões. A tarefa era criar um modelo que determinasse qual usuário é humano e qual é robô. Os usuários humanos de um site de leilões não querem ter que competir com máquinas. Um site com muitos robôs pode sofrer um êxodo de usuários humanos, então é importante ter um sistema robusto de detecção de usuários suspeitos.

Segundo os administradores, estes dados não vinham do Facebook, mas acredito que o interesse era criar um ambiente parecido com o sistema de leilões de publicidade deles.

Esta competição tornou-se bastante especial para mim, já que foi meu primeiro Top 10 no Kaggle, e acabei me tornando Kaggle Master. Este era meu primeiro objetivo, mas também tem o bônus de ser considerado para uma entrevista para a vaga de Data Scientist no Facebook.

Dados

Uma das partes mais interessantes desta competição, e que a tornava bastante próxima de um caso real de ciência de dados, era o fato de não termos um arquivo bruto com variáveis relacionadas ao objetivo, que pudesse ser utilizado diretamente com algum modelo de machine learning.

Foram disponibilizados dois arquivos: o primeiro continha cerca de 2000 linhas com informações anônimas dos usuários, como ID e endereço. O segundo continha mais de 7 milhões de lances dados no site.

As variáveis do segundo arquivo incluíam: ID do lance, ID do usuário, ID do leilão, categoria do leilão, dispositivo do usuário, tempo, país, IP e URL. O tempo estava codificado, de maneira que não era possível descobrir a data exata dos lances.

Validação

Ter um ambiente de validação robusto é essencial para qualquer tarefa de machine learning. Neste caso havia o desafio de termos apenas 2000 exemplos para treinar e validar. Após testar várias combinações, decidi utilizar uma validação cruzada estratificada, com 5 divisões, e repetir essa validação por 5 vezes.

A repetição acaba fazendo com que os exemplos sejam colocados, aleatoriamente, em divisões diferentes, o que ajuda a estabilizar a estimativa de erro, já que um ou outro exemplo particular terá menor influência na média.

Criação de variáveis

A parte mais importante de qualquer trabalho de ciência de dados é garantir que você tenha os dados corretos para fazer a previsão. Um erro muito comum é achar que é só colocar os dados no algoritmo, e ele vai aprender o que você quiser. Esta foi, de longe, a parte mais importante da minha solução.

Contagens

O primeiro passo foi criar variáveis com a contagem de quantos lances cada usuário tinha feito. Estendendo esta ideia, criei variáveis para a contagem de quantos dispositivos, IPs, categorias, países e URLs diferentes um mesmo usuário tinha utilizado. Com estas variáveis simples, uma Random Forest, sem tuning, atingia 0,86 AUC na validação cruzada com 5 divisões.

No caso dos países, além de integrar a variável própria para contar o número de países diferentes, fiz variáveis individuais para cada um deles, contando quantos lances o usuário tinha dado vindos de cada país. Isto ajuda o modelo a identificar “países de risco”, onde pode ser mais fácil hospedar robôs para dar os lances.

Estatísticas simples

Depois de ver que as contagens possuíam um grande valor informativo, decidi testar estatísticas simples, como médias e desvios. Exemplos de variáveis importantes:

– Média e desvio padrão de lances por leilão;
– Número máximo de lances dados num mesmo leilão, com a mesma timestamp;
– Média de tempo entre lances.

Temporais

A variável que indicava o tempo nos dados originais estava codificada de maneira que não dava para saber a data exata em que os lances ocorreram. Mas, após analisar os dados, percebi que mesmo não sabendo a data exata, era possível identificar o que seriam dias, horas e minutos.

Isto me permitiu criar variáveis relativas a médias, desvios, proporções e contagens baseadas em unidades de tempo. Algumas variáveis importantes:

– Desvio padrão de lances por dia;
– Média de lances por dia;
– Média de lances na mesma timestamp;
– Contagem de lances por hora do dia;
– Média de leilões que o usuário participou no mesmo dia.

Parece que todos os participantes que se posicionaram no topo da tabela conseguiram identificar estes padrões temporais.

Proporções

Algumas proporções baseadas nas contagens e no tempo também foram importantes. Alguns exemplos são:

– Proporção máxima do total de lances em um dia atribuídos ao usuário;
– Densidade de lances por hora do dia

Modelos

Quando você tem variáveis que contêm bastante informação sobre o evento, a parte de seleção e ajuste de parâmetros do modelo se torna mais simples. Neste caso, apesar de testar brevemente modelos lineares, acabei optando por explorar melhor modelos baseados em ensembles de decision trees.

Random Forest

O primeiro modelo que testei, e o que se provou mais útil neste caso.

Acredito que sua performance superior se deu por ser um modelo mais voltado para a diminuição da variância. Por ter poucos dados para treinar, era importante que o modelo conseguisse lidar bem com valores extremos, e estabilizasse as previsões.

Uma Random Forest com 500 árvores e parâmetros ajustados através da validação cruzada apresentava 0,9112 AUC na medição local, e 0,9203 na LB.

Gradient Boosted Trees

Este é um dos modelos mais utilizados em competições de ciência de dados porque reconhecidamente oferece uma boa performance em diversos tipos de tarefas. Neste caso não foi assim. Os meus melhores modelos ficaram por volta dos 0,90-0,91 AUC na LB apesar de alcançarem a 0,94 na validação.

Isso levanta a seguinte questão: temos um modelo baseado em decision trees que nos dá uma boa performance (Random Forest) e, em tese, boosting melhoraria a pontuação, por que este não é o caso aqui?

A minha resposta é: Random Forest aumenta o bias, ou seja, torna o modelo mais rígido, com menor variância, enquanto o GBT aumenta a variância, sendo mais sensível ao ruído presente nos dados. Neste caso temos menos de 2000 exemplos, e um forte desequilíbrio entre as classes. Estas duas razões combinadas geram dados com uma boa quantidade de variância, aumentá-la ainda mais vai causar overfitting.

A solução seria fazer o tuning e regularizar mais o GBT, mas com poucos dados para validar fica difícil confiar nos resultados. Por isso, decidi confiar nas razões teóricas e continuar trabalhando com a Random Forest para ter um modelo mais estável.

Ensembles

O toque final de boa parte das melhores soluções de competições é fazer um ensemble com os melhores modelos.

Ensemble Aleatório

Devido ao pequeno número de exemplos, era difícil confiar nos parâmetros encontrados através da validação. Por isso, decidi criar um ensemble com Random Forests “aleatórias” nos parâmetros. Cada uma delas tinha uma seed e alguns parâmetros diferentes. Este ensemble atingiu uma AUC de 0,92 na estimativa out-of-bag e na LB.

AdaBoost Random Forest

Apesar da instabilidade do GBT, decidi aplicar boosting à Random Forest. Talvez aplicando este método em modelos mais estáveis que as decision trees pudesse ajudar, mesmo com a pequena quantidade de exemplos. E foi isso mesmo que aconteceu.

Na validação obtive 0,925 de AUC, e 0,928 na LB.

Ainda assim havia uma diferença entre a pontuação da validação e da LB quando eu adicionava novas variáveis.

Estabilizando as seeds

Para estabilizar as previsões decidi rodar o modelo de Boosted Random Forest usando várias seeds diferentes, e fazer a média. Isso não deu um bom resultado na LB, mas apresentava estabilidade na validação cruzada, o que me fez confiar no modelo.

Conclusão

No fim, utilizei um modelo de Boosted Random Forest com um conjunto de atributos que apresentava AUC 0,94 na validação cruzada, e foi basicamente minha pontuação no test set final da competição, que me garantiu o 6º lugar.

Big Data, Kaggle, Python, Regressão, Séries Temporais

Machine Learning Para Previsão de Vendas Usando Dados Meteorológicos

jul 20, 2015 por Mario Filho

O WalMart é uma rede com milhares de lojas em 27 países. É possível encontrar vários artigos sobre os mecanismos tecnológicos utilizados para gerenciar a logística e distribuição dos produtos. É a segunda vez que eles oferecem uma competição no Kaggle com a intenção de encontrar candidatos para entrevistas para vagas de cientistas de dados.

Uma grande vantagem deste tipo de competição é termos acesso a dados de grandes companhias, e entender quais são os problemas que eles estão tentando resolver com modelos probabilísticos.

O objetivo da competição era criar um modelo que pudesse prever a quantidade de venda de alguns produtos, em lojas específicas, nos dias antes e depois de nevascas e tempestades. O exemplo dado por eles na descrição da tarefa foi a venda de guarda-chuvas, que intuitivamente deve ver um aumento antes de uma grande tempestade.

Dados

Para treinar o modelo foram oferecidos dois arquivos: um deles continha informações sobre a identificação das lojas, produtos, e as estações meteorológicas mais próximas. O outro continha dados meteorológicos de cada estação.

No total foram disponibilizados dados de 111 produtos cujas vendas podem ser afetadas pelas condições climáticas, 45 lojas, e 20 estações meteorológicas. O objetivo era prever a quantidade de cada produto, em cada loja, que seria vendida 3 dias antes, 3 dias depois, e no dia do evento climático.

Um evento climático era considerado caso tivesse sido registrada mais que uma polegada de chuva, ou mais que duas polegadas de neve naquele dia.

A combinação entre lojas e produtos nos dava cerca de 4,6 milhões de exemplos para treino, e cerca de 500 mil para teste. Cada exemplo se referia a um dia, loja e produto.

Métrica

A métrica utilizada foi o Root Mean Square Log Error. É basicamente o RMSE aplicado à transformação log(Y + 1) das previsões. Isso significa que erros em previsões que deveriam ser próximas de 0 seriam punidos mais severamente do que erros em previsões com números mais altos. Por exemplo, prever 5 itens quando deveria ser 0, tem uma punição maior do que prever 105 quando deveria ser 100.

Transformação dos dados

Como eu só tinha 10 dias para trabalhar nesta competição, decidi verificar até onde era possível chegar com um modelo baseado apenas nas variáveis originais e, quem sabe, fazer um ensemble.

Uma diferença desta para outras competições é que, apesar dos dados virem organizados, você era responsável por unir as variáveis climáticas com os dados identificados de produtos e lojas. Faz todo o sentido, pois o Walmart não iria querer um data scientist que não saiba manipular dados.

Para isto utilizei Pandas. Esta biblioteca Python para manipulação de dados é uma das mais utilizadas, e lembra bastante as estruturas de dados disponíveis em R.

Num primeiro momento usei todas as variáveis como numéricas, e treinei um XGBoost com leve tuning, excluindo as variáveis que codificavam alertas especiais, e usando 10% do dataset para determinar o número de árvores. Como já era de se esperar, o resultado foi ruim, cerca de 0,1643 na LB.

Binarizando as variáveis

Depois de testar o primeiro modelo, codifiquei as variáveis categóricas com one-hot encoding. Ou seja, para cada nível foi criada uma coluna com o indicador 0 ou 1, caso a variável estivesse presente naquele exemplo. Normalmente o número de colunas deve ser o número de níveis menos um, para não ter problemas com colinearidade. Como eu planejava usar modelos que não eram sensíveis a esse problema, eu não me preocupei em excluir uma das colunas.

Após fazer tuning usando um subset com 10% dos dados, e obter o RMSLE de 0,1216 na validação, enviei a solução, e obtive o valor de 0,1204 na LB, uma boa melhora sobre o anterior.

Imputação

Muitos dados meteorológicos não estavam presentes, então decidi testar um método simples de imputação: substituir os NaNs pela média dos valores da coluna. Após fazer o tuning novamente, agora adequado a esses novos valores, obtive o RMSLE de 0,1140 nos 10% de validação e 0,1095 na LB.

Variáveis temporais

Não explorei muito a dimensão temporal dos dados, em uma das tentativas adicionei informações meteorológicas do dia anterior, o que reduziu o erro para 0,1083.

Subdivisões dos dados

Um método que funcionou muito bem na minha primeira competição, e que sempre acabo tentando, é dividir o training set em pequenos subsets relacionados a alguma variável e treinar um modelo específico para cada um. Neste caso, decidi criar 19 modelos, um para cada uma das 19 estações meteorológicas presentes no test set. O RMSLE na LB foi de 0,1101. Pior do que com um modelo que trata as estações como variáveis no dataset inteiro.

Um problema grave com esta abordagem é tentar usar o mesmo modelo, com os mesmos parâmetros, para datasets diferentes. Sabendo disso, decidi fazer um pequeno tuning dos parâmetros para cada dataset, o que reduziu o RMSLE da LB para 0,1069.

Apesar da pequena diferença, parece que a divisão em modelos individuais para cada estação capturava algumas informações que não estavam presentes no modelo que considerava todas juntas.

Modelos

Dos modelos que testei, dois se destacaram: Gradient Boosted Trees (XGBoost) e Random Forests.

Random Forest

Eu tinha usado Random Forest para regressão apenas uma vez, em um trabalho, mas nunca com uma quantidade grande de dados. Após fazer o tuning dos parâmetros, aplicando o modelo nos dados imputados, ela resultou num RMSLE de 0,1116 na LB.

XGBoost

O XGBoost apresentou o melhor erro de um modelo individual. Além de ajustar parâmetros, como a profundidade das árvores, usando um subset dos dados, era necessário ajustar a quantidade de árvores e a learning rate. Normalmente uma learning rate pequena e uma grande quantidade de árvores é a receita segura para melhorar a performance, em troca de mais tempo para treinar o modelo.

Como o XGBoost é sensível à seed do RNG, decidi fazer um ensemble de XGBs mudando apenas este valor. Este método melhorou marginalmente a minha pontuação em outras competições, mas nesta o impacto dele foi maior pelo seguinte fato: o XGBoost possui uma função que permite deixar dados separados para que ele determine o número de árvores que minimiza o erro. Neste caso decidi usar a função, e deixei 5% dos dados separados. Além de variar a seed para o próprio XGB, eu variei a seed para fazer a divisão dos dados, isso fez com que os modelos ficassem mais diversos, o que é essencial para um bom ensemble.

O máximo que cheguei a treinar foram 10 XGBoosts neste esquema. Apesar de ser um modelo bastante estável, o RMSLE do ensemble foi de 0,1041, apresentando uma redução comparado a 0,1095 do modelo individual.

Solução Final

No fim, juntei todas as soluções que eu havia enviado, e acabei obtendo um RMSLE de 0,1028, garantindo uma posição entre as 20% melhores.

Possíveis Melhorias

Um dia após o término da competição eu revisei as variáveis do meu XGBoost, e descobri que as variáveis que identificavam os produtos (item_nbr) não estavam em formato binário, e eram consideradas por ele as mais importantes. Com a codificação correta acredito que seria possível diminuir mais o erro, e alcançar uma posição final melhor.

Classificação, Kaggle

Classificando 150 Mil Produtos em Categorias Diferentes Usando Machine Learning

jul 6, 2015 por Mario Filho

O grupo Otto é uma das maiores empresas de comércio online do mundo.
Segundo eles, devido à diversidade da infraestrutura global da empresa, muitos produtos similares acabam classificados em categorias incorretas. Por isso eles disponibilizaram dados sobre aproximadamente 200 mil produtos, pertencentes a nove categorias. O objetivo era criar um modelo probabilístico que classificasse corretamente os produtos dentro de cada categoria.

Dados

Para treino foram disponibilizados cerca de 60 mil produtos e, para teste, cerca de 150 mil.

Os 93 atributos não foram identificados, a única informação dada sobre eles é o fato de serem variáveis numéricas. Isso dificulta bastante para entendermos a melhor maneira de trabalhar com os dados, e também o trabalho de criação de novos atributos.

Além disso, a tarefa envolve múltiplas classes, ou seja, temos que gerar probabilidades sobre 9 classes.

Métrica

A métrica escolhida para esta competição foi a log loss multiclasse.

Esta métrica pune probabilidades muito confiantes em previsões erradas, e é sensível ao desequilíbrio entre as classes. Algumas delas representavam 20% dos dados de treino, enquanto outras, 9%. Neste caso, tentar equilibrar as classes, seja com subsampling ou penalizando as classes menores, não vai ajudar.

Criação de Atributos

Mesmo com as variáveis anônimas, decidi explorar a possibilidade de integrar interações entre elas, que talvez não fossem capturadas pelos modelos.

Básico

Procurei criar novos atributos baseados em somas, subtrações, razões e multiplicações entre os originais. O único atributo que contribuiu significativamente foi a soma de todos os atributos para um determinado produto. A maioria dos atributos era bastante esparsa (tinha mais zeros do que outros valores), então talvez isso tenha contribuído para a falta de sucesso.

Usando um GBM/RF para Selecionar Atributos

Como nós estamos falando de um espaço gigantesco de possibilidades de combinações de atributos, decidi usar uma técnica usada por um participante de outra competição e publicada neste link: http://trevorstephens.com/post/98233123324/armchair-particle-physicist

Basicamente consiste em criar alguns datasets com as interações desejadas e treinar um modelo de Gradient Boosted Trees em cada um. Após isso, verifica-se quais são os atributos mais importantes em cada dataset. A lógica é: se um atributo for importante em vários datasets diferentes, ele deve ser importante no geral.

Eu fiz isso nestes dados, usando 2 divisões estratificadas, mas apesar deles concordarem na importância de algumas interações, elas não melhoraram meu modelo base. Como eu já havia passado dias me dedicando a essa parte, eu decidi que ia me focar mais em ajustar modelos para um ensemble.

Caso eu tivesse continuado, talvez selecionar as X melhores interações e reajustar os hiperparâmetros de algum modelo pudesse extrair valor destas variáveis.

Modelos

Gradient Boosted Trees (XGBoost)

O meu modelo individual, com a menor log loss, foi criado com o XGBoost, que é uma implementação rápida e paralela das poderosas Gradient Boosted Trees. Esse modelo já esteve em soluções vencedoras de boa parte das competições e normalmente apresenta uma performance superior.

Para encontrar o melhor grupo de atributos usei uma busca aleatória. Fixei a learning rate em 0,05 e variei atributos como a profundidade das árvores, mínimo de exemplos que deveriam compor um nó, e a proporção de exemplos e variáveis que o algoritmo deveria selecionar aleatoriamente para criar cada árvore.

Normalmente quanto menor a learning rate, melhor a precisão, mas são necessárias mais árvores, o que aumenta o tempo de treino. No fim, deixei este valor em 0,01 e treinei 2000 árvores.

Estes atributos servem para controlar o overfitting. Após encontrar os melhores atributos, a log loss da validação cruzada em 3 divisões era 0,4656, e na leaderboard, 0,4368.

Random Forests

Apesar de não ter muito sucesso inicial com as Random Forests nesta competição. Uma sugestão dada no fórum foi de calibrar as previsões. Recentemente a equipe desenvolvedora do scikit-learn disponibilizou uma nova ferramenta que nos permite ajustar os valores de saída de um modelo de maneira que eles se tornem mais próximos das probabilidades reais.

Uma Random Forest normalmente faz uma votação entre seus componentes, e a proporção de cada classe é dada como a probabilidade do exemplo pertencer àquela classe. Essas proporções não correspondem às probabilidades reais dos eventos, então só teremos probabilidades relevantes se fizermos o ajuste.

Como esta competição pedia para prevermos a probabilidade de um produto pertencer a uma categoria, era muito importante ter as probabilidades ajustadas corretamente. Neste caso, utilizei a nova ferramenta do scikit-learn com 5 divisões dos dados. Em 4 delas ele treinava a Random Forest e, na restante, uma regressão isotônica. Isso era repetido 5 vezes, e a média da probabilidade das 5 regressões era o resultado.

Isso melhorou bastante a previsão da Random Forest e, apesar de também usarem decision trees, como o GBM, elas acabaram contribuindo significativamente para um ensemble.

Redes Neurais

Tive a oportunidade de conhecer dois módulos ótimos para treinar redes neurais de maneira simples em Python. Eles são Lasagne e NoLearn. Com eles é simples criar e treinar redes neurais de última geração, usando inclusive GPUs para processamento.

Mas não se engane, apesar deles facilitarem a implementação, uma rede neural exige uma grande quantidade de decisões para se tornar útil. Temos que determinar desde a arquitetura, até os métodos de inicialização dos pesos e regularização. Alguns pesquisadores sugerem fazer uma busca aleatória por parâmetros iniciais e continuar a ajustar manualmente.

Neste caso, a arquitetura que me serviu bem foi a seguinte: 3 layers de neurônios com o seguinte número de unidades em cada uma 768-512-256. Eu treinei duas versões, uma delas com dropout de 0.5 entre os hidden layers. E outra com dropout de 0.5 apenas no primeiro hidden layer. As unidades dos hidden layers eram ReLu, e a unidade de saída, softmax.

Uma interpretação interessante, de um outro competidor que chegou a uma arquitetura parecida, é de que o primeiro hidden layer permite que a rede faça projeções aleatórias dos dados, e os layers menores, que vêm após este, buscam uma sintetização das informações. Redes neurais são naturalmente difíceis de interpretar, mas eu achei este um ponto de vista bem interessante.

Por fim, como a mínima de uma rede neural é local, e dependente dos pesos iniciais (e aí temos muitos pesos iniciais), decidi fazer a média de 5 redes neurais com pesos iniciais diferentes. Esta média me deu um score de 0.4390 na LB, comparável ao XGBoost.

Outros

Ainda tentei treinar: SVM, Nearest Neighbors, Regressão Logística, mas nenhum deles apresentava uma boa performance sozinho, ou contribuía significativamente para o ensemble.

Ajustando os Hiperparâmetros

Já que não tínhamos acesso ao conteúdo das variáveis, era vital ajustarmos os parâmetros para termos os melhores modelos. No início fiz uma busca aleatória com valores razoáveis, mas que permitissem aos modelos explorar um pouco o espaço.

Normalmente os hiperparâmetros de um modelo não são independentes, ou seja, ajustar um de cada vez não vai te levar à melhor combinação. Por isso é importante fazer uma exploração com combinações de parâmetros. Infelizmente na maioria das vezes é impossível testar todas as combinações possíveis, mas quando fazemos a busca aleatória podemos ter uma ideia do espaço de soluções, sem ter que explorá-lo exaustivamente.

Depois de feita esta exploração, é bom variar manualmente alguns deles, e ver se é possível melhorar a performance em alguma combinação vizinha.

Resultado

No fim minha solução era formada por uma média ponderada dos seguintes modelos:

– Gradient Boosted Trees (XGBoost) sobre as variáveis originais e a soma das colunas.
– Random Forest sobre as variáveis originais, log(X+1) delas, e soma das colunas.
– Rede neural com 3 layers e dropout de 0.5 em cada um, sobre as variáveis originais, log(X+1) delas, e soma das colunas.
– Rede neural com 3 layers com dropout de 0.5 no primeiro hidden layer, sem dropout no segundo, sobre as variáveis originais, log(X+1) delas, e soma das colunas.

Além desses modelos, adicionei um bias, simplesmente prevendo a proporção de cada classe como probabilidade de um exemplo pertencer a elas.

Esta solução atingiu a log loss de 0,4080, e me garantiu a posição 35 entre 3514 times (Top 1%).

Soluções vencedoras

A solução do time vencedor consistia de um ensemble de três camadas:

Na primeira camada foram treinados 33 modelos diferentes, variando tanto o tipo do modelo, quanto as variáveis utilizadas para treinar cada um. Alguns deles eram múltiplos modelos treinados com bagging.

Na segunda camada, as previsões destes modelos eram utilizadas para alimentar um XGBoost, uma rede neural, e um ExtraTrees com Adaboost.

Na terceira, e última, camada, foi feito o bagging dos três modelos da segunda camada (totalizando 1100 modelos) e depois uma média ponderada entre os três.

Além disso, foram criadas novas variáveis, baseadas na distância do exemplo para os exemplos mais próximos de cada classe, e também outras baseadas em transformações TF-IDF e T-SNE do dataset original.

Big Data, Classificação, Kaggle, Python

Usando Machine Learning Para Identificar Motoristas Através de Dados GPS

jun 22, 2015 por Mario Filho

Nos últimos anos a indústria de seguros tem buscado maneiras de aprimorar seus modelos usando Machine Learning. Uma delas é utilizar dados que vão além de um formulário preenchido pelo segurado para determinar o risco de acidentes.

Um dos métodos utilizados é usar dados comportamentais do motorista, obtidos através de rastreamento via GPS. Desta maneira acredita-se ser possível capturar informações e padrões de perfil que vão além dos métodos tradicionais.

O que torna um motorista diferente do outro? Identificar o motorista que está conduzindo um veículo durante uma viagem é o primeiro passo para conseguirmos um bom modelo.

Descrição da Competição

A seguradora AXA decidiu disponibilizar viagens anônimas de aproximadamente 2700 motoristas. Para cada motorista são 200 viagens, sendo que algumas delas não são do motorista em questão, e a tarefa era criar um modelo que identificasse quais eram as viagens atribuídas a este motorista incorretamente.

Para evitar o vazamento de informações, os dados foram centralizados e rotacionados de maneira aleatória, além de ter alguns trechos das viagens removidos do início e do fim.

Dados

No total foram disponibilizadas 200 viagens de 2736 motoristas. Dentro da pasta de cada motorista havia uma quantidade aleatória de viagens falsamente atribuídas ao motorista em questão. Uma informação importante é que tínhamos a garantia que a maioria das viagens pertencia verdadeiramente ao motorista.

Cada viagem era descrita em um arquivo CSV com um valor numérico para a posição x e y, que seria a distância da origem, em metros, e cada linha correspondia a 1 segundo de viagem.

Métrica

A métrica utilizada para avaliação foi a área sob a curva ROC, abreviada por AUC.

Atributos

Os atributos que utilizei se baseiam principalmente nas três primeiras derivadas do deslocamento: velocidade, aceleração e arrancada. Outras derivadas acabavam por prejudicar a performance.

Frequências – FFT

Procurei extrair frequências, tanto dos dados originais brutos, quanto das derivadas. Fiz isto com a função FFT (Fast Fourier Transform) disponível no Numpy. Dividi o intervalo de frequências por 20, e o valor dos atributos era a soma das frequências do respectivo intervalo.

Neste caso eu extraí as frequências isoladas para cada componente (x e y). Talvez obter as frequências usando a velocidade e aceleração, considerando os dois componentes juntos, seja melhor, mas não testei esta hipótese.

O AUC máximo que consegui, usando frequências das três derivadas, foi 0,7058, com uma regressão logística sem otimização de parâmetros. Apesar de não ter testado, acredito que otimizando os parâmetros do modelo para cada motorista este número poderia ser aumentado até um máximo de 0,73.

Histogramas

Esta é a transformação que encontrei em todos os trabalhos acadêmicos que li sobre o assunto. Consiste em fazer um histograma dos atributos. Neste caso, ganhamos um novo parâmetro, o tamanho dos intervalos.

A melhor maneira que encontrei foi dividir a magnitude da velocidade em 50 intervalos, e a aceleração, em 10.

Dois detalhes que melhoraram a performance: utilizar o atributo “density” da função do Numpy, que calcula a função densidade da probabilidade, e normalizar os resultados para que sua soma fosse 1.

Valores em Porcentagens

Outra sugestão bastante frequente era utilizar o valor dos atributos em determinadas porcentagens. Isso ajuda a descobrir quão extremos são os valores dos mesmos.

Um exemplo pode facilitar o entendimento: no caso da velocidade, para obter o valor presente no percentil 5%, ordenei os valores de maneira crescente e localizei o ponto da lista que representava a fronteira dos primeiros 5% dos dados ordenados.

Se tivermos uma lista com 100 valores de velocidade ordenados, com um índice que vai de 1 a 100, o valor presente na posição 5 seria o desejado. O valor presente na posição 50, por exemplo, seria a mediana.

Isso é importante porque diferencia motoristas que possuem um perfil mais rápido ou mais lento.

Eu utilizei os percentis 5% a 95%, em intervalos de 5%.

Viagens Repetidas

Uma abordagem que foi importante para que os melhores colocados subissem na tabela foi criar um algoritmo para encontrar viagens similares. Era notável que alguns trechos, apesar de possuírem um ângulo diferente com relação à origem, se tratavam do mesmo caminho.

Eu fiz um algoritmo simples, baseado na velocidade do motorista durante a viagem mas, apesar de ter encontrado algumas viagens similares, acabava prejudicando a pontuação.

Os algoritmos que tiveram sucesso eram mais complexos, e em alguns casos envolviam a transformação das variáveis, redução de dimensões, e limites para determinar a similaridade. Dentre todos, dois métodos pareceram essenciais:

Ramer–Douglas–Peucker: como cada viagem tinha uma duração diferente, este método foi utilizado para reduzir, e igualar, o número de pontos.

Dynamic Time Warping: este método foi utilizado para comparar os trechos e determinar a similaridade entre eles.

Amostragem

Como não era possível saber realmente quais eram as viagens de cada motorista. Ou seja, existiam viagens marcadas como positivas, quando na realidade eram negativas, havia a oportunidade de se decidir a melhor maneira de fazer a amostragem dos dados.

Aleatório

Num primeiro momento, com métodos não supervisionados, usei apenas as 200 viagens de cada motorista. Uma sugestão no fórum da competição foi usar as 200 viagens de um motorista como exemplos positivos e usar viagens de outros motoristas como exemplos negativos.

Em geral esta foi uma abordagem bastante positiva. Dentre as maneiras de fazê-la estavam:

– Escolher, aleatoriamente, uma quantidade N de motoristas e tratar todas as viagens deles como negativas.

– Escolher, aleatoriamente, uma quantidade V de viagens de N motoristas diferentes, para usar como negativas.

A abordagem que me serviu melhor foi baseada na segunda opção: selecionar um pequeno número de viagens de vários motoristas diferentes.

Amostras similares/diferentes

Outra alternativa que testei foi selecionar amostras negativas de motoristas que eram mais similares ou diferentes das amostras positivas, de acordo com uma medida como distância euclidiana.

Criei atributos como a média de velocidade, distância percorrida, e duração das viagens para cada motorista e determinei quais eram os mais parecidos, e quais eram os mais diferentes do motorista em questão.

Apesar da abordagem utilizando os motoristas “mais diferentes” mostrar uma melhora durante a validação cruzada, ela não melhorou a pontuação no test set. A razão que atribuo a isto é o fato dos dados da validação conterem viagens que são falsas, mas estão marcadas como verdadeiras, enquanto que no test set temos a realidade absoluta.

Reforçando as previsões

Procurei reforçar as previsões do modelo treinando mais de uma vez, mas agora usando as previsões do modelo anterior no lugar da variável dependente original, da seguinte maneira:

No primeiro passe fiz o procedimento normal, usando as viagens do motorista como positivas, e outras aleatoriamente selecionadas de outros motoristas como negativas para obter previsões.

Nos próximos passes eu utilizei as previsões do modelo anterior como variável dependente. Isso melhorou bastante na validação cruzada, mas apresentou uma melhora muito pequena nos dados de validação do Kaggle.

Modelos e Métodos

Clustering

Num primeiro momento decidi agrupar as viagens de acordo com características similares. Essa ideia foi baseada na suposição que, como a maioria das viagens pertencia ao motorista, seria possível separar as viagens em dois clusters, e aquele que tivesse a minoria das viagens seria o “negativo”.

Outra ideia era aumentar o número de clusters, já que existiam diferenças das viagens entre si, mesmo pertencentes ao mesmo motorista.

De qualquer maneira, estas ideias se mostraram inúteis, com AUC próximo de 0,5.

Outra maneira não supervisionada que tentei foi calcular a similaridade de cosseno entre a média dos atributos das viagens e cada uma delas. Neste caso a “probabilidade” de uma viagem pertencer ao motorista seria um menos este valor. O que também se mostrou inútil.

Regressão Logística

O meu modelo com maior sucesso foi a regressão logística. Isso me surpreendeu, já que o sucesso descrito por vários participantes no fórum envolvia boosted trees.

Utilizei cerca de 1000 exemplos, 200 positivos e 800 negativos. Um detalhe importante foi ativar a opção “class_weight=auto”. Esta opção atribui uma penalidade diferente aos exemplos, de acordo com a distribuição das classes. Neste caso, um erro ou acerto na classe positiva valia mais do que o mesmo na classe negativa. Desta maneira, apesar de termos mais exemplos positivos, a penalização por erros era equilibrada.

Treinei uma regressão para cada motorista, personalizando a penalidade de regularização (L1 ou L2) e o coeficiente da mesma. Cada modelo fazia a validação cruzada de valores para este coeficiente utilizando os dois tipos de regularização, e no fim selecionava o tipo de regularização e o coeficiente que apresentasse a melhor AUC.

O melhor modelo individual atingiu a AUC de 0,86.

Random Forests/Boosted Trees

De acordo com os fóruns, o modelo mais utilizado foi Gradient Boosted Trees. Basicamente ele consiste em treinar uma sequência de decision trees pouco complexas, atribuindo um peso maior a exemplos que foram classificados incorretamente pelas árvores anteriores. A probabilidade de um exemplo específico é uma combinação linear das previsões de cada árvore.

No meu caso as tentativas de utilizar este modelo não foram bem sucedidas. Dois motivos podem ter contribuído:

– Por desconfiança da validação, eu acabei não dando muita atenção à otimização dos parâmetros dele para conseguir uma boa AUC.
– O processo de validação cruzada para encontrar bons parâmetros ia demorar muito, e testes preliminares não haviam apresentado bons resultados.

De qualquer maneira treinei um modelo de Random Forests, com os parâmetros padrões, para utilizar num ensemble. Ele sozinho possuía AUC de 0,847.

Ensemble

Próximo do fim da competição decidi partir para a forma mais “fácil” de se melhorar a pontuação: ensembles.

Para uma parte do ensemble treinei os seguintes modelos:

– Regressão Logística apenas com regularização L1
– Regressão Logística apenas com regularização L2
– SVM RBF
– SVM Linear
– Random Forests

Apesar do SVM não oferecer probabilidades, acreditei que ele seria uma boa contribuição para o ensemble, o que se mostrou verdadeiro.

A média simples da previsão destes modelos produziu uma AUC de 0,8825.

Por fim, juntei estes com todos os outros que já havia produzido desde o início da competição, para uma AUC final de 0,8915. Isto me garantiu uma colocação entre as 15% melhores soluções da competição.

Métodos das soluções vencedoras

As soluções vencedoras usaram 400 exemplos para cada motorista, selecionando viagens aleatórias de outros motoristas como exemplos negativos. O modelo mais comum foi Gradient Boosted Trees. Além disso, algoritmos complexos para encontrar viagens similares e atribui-las diretamente ao motorista conquistaram os pontos finais necessários para se posicionar entre as 10 melhores soluções.

Otimização, Python

Usando Otimização Para Aproximar a Menor Rota Entre Mais de 5.500 Municípios Brasileiros

jun 3, 2015 por Mario Filho

Viajar é a opção de lazer de muitas pessoas, e hoje os aplicativos que calculam rotas e nos ajudam a chegar nos lugares utilizando GPS são bastante úteis. Eles resolvem um dos problemas clássicos da otimização: o problema do caixeiro viajante (Traveling Salesman Problem). É bem provável que este seja o problema mais estudado na área de otimização discreta.

Por que ele é importante?

Imagine que você seja um representante comercial e precise visitar clientes em 5 cidades diferentes. Cada cidade deve ser visitada apenas uma vez e, ao final, você volta para a primeira. Você precisa escolher a menor rota possível entre elas. Quantas opções de rotas existem?

Você precisa escolher uma das cinco cidades para começar, depois uma das quatro para ser a próxima, aí ficam mais três… e assim vai. Matematicamente, você tem 5 x 4 x 3 x 2 = 120 opções. Neste caso pequeno seria possível pedir que o computador calculasse a distância percorrida em todas as combinações e retornasse a menor.

Mas e quando temos mais cidades? Se falarmos apenas das capitais brasileiras, temos 27 cidades (26 estados e um Distrito Federal). A quantidade de combinações é um número de 28 zeros.

Isso significa que precisamos de métodos que calculem estas rotas de maneira mais inteligente, sem precisarmos verificar todas as combinações.

Eu andei pesquisando na internet, mas não encontrei algum lugar em que houvesse o cálculo de uma rota entre todos os municípios brasileiros. Então a tarefa aqui será aproximar a menor rota entre mais de 5.000 municípios brasileiros.

O número de municípios brasileiros não é fixo, mas a melhor lista que consegui, com as posições geográficas de cada um, contém 5.509 municípios. Ela foi retirada originalmente deste link.

Existem alguns métodos para calcular esta rota, sem eles seria impossível calcular, por força bruta, qual é a melhor, mesmo que tivéssemos todos os computadores do mundo à nossa disposição.

Os arquivos utilizados neste artigo estão disponíveis neste repositório GitHub. O script usa Python puro, então recomendo rodar o script usando PyPy para ir mais rápido.

Dados

Para cada município temos o nome, estado, e a latitude e longitude. A distância entre eles será calculada baseada nessas informações. Apesar de haver diferença com a realidade, já que as estradas de uma cidade a outra não formam, necessariamente, uma linha reta, esta será uma boa aproximação.

Dificilmente chegaremos ao mínimo absoluto com o método demonstrado aqui, mas ele é uma maneira de se aproximar bastante deste número. Para efeitos práticos, aproximar-se do mínimo absoluto sem precisar de um grande poder computacional, e em pouco tempo, é o bastante.

Preparando os dados

Com o XLS transformado em CSV, vamos criar funções dentro do script Python para processar os dados, colocando-os em estruturas adequadas para rodarmos os algoritmos.

Este script é baseado no código disponibilizado no curso Discrete Optimization, e na implementação que utilizei para resolver uma das tarefas, que são alguns TSPs.

Primeiro criamos uma estrutura de dados para armazenar as informações de cada cidade, usaremos a named tuple. Esta estrutura cria uma subclasse da tuple original cujos dados podem ser acessados usando nomes, em vez dos índices da tuple padrão.

Agora criamos uma função para carregar o arquivo CSV e armazenar cada linha numa named tuple contendo o ID, latitude, longitude, nome e estado onde fica a cidade. Estas serão armazenadas numa lista. Retornamos a lista das cidades, e a quantidade.

town = namedtuple("town", ['id','name','lat', 'long'])

def parseInput():
    
    with open('MunicipiosBrasil.csv','r') as input_data_file:
        town_list = [l.split(';')[:4] for l in input_data_file.read().split('\r')[1:]]
        
        town_list = [town(int(t[0])-1,t[3] ,math.radians(float(t[1].replace(',','.'))), math.radians(float(t[2].replace(',','.')))) for t in town_list]
        
        townCount = int(len(town_list))

        return town_list, townCount

town = namedtuple("town", ['id','name','lat', 'long'])

def parseInput():

with open('MunicipiosBrasil.csv','r') as input_data_file:

town_list = [l.split(';')[:4] for l in input_data_file.read().split('\r')[1:]]

town_list = [town(int(t[0])-1,t[3] ,math.radians(float(t[1].replace(',','.'))), math.radians(float(t[2].replace(',','.')))) for t in town_list]

townCount = int(len(town_list))

return town_list, townCount

Note que, decidi subtrair 1 de cada ID para termos um índice começando em 0, e além de formatar a latitude e longitude de maneira que o Python pudesse transformá-las em floats, já aproveitei para transformá-las em radianos. Originalmente elas estão em graus decimais. Isso facilitará os cálculos de distância.

Resolvendo o problema

Mixed Integer Program

Uma maneira de obter uma solução bastante próxima da melhor é usar um software que resolva MIPs. Este softwares resolvem problemas de otimização discreta aproximando-os com Programação Linear, de maneira contínua, e depois explorando árvores binárias cujos nós são soluções com algumas variáveis discretas fixas.

Existem várias formulações para resolver este problema, e é possível saber quão próximos estamos do ponto mínimo (o intervalo entre o ponto mínimo e a nossa solução) durante o processo de solução.

O melhor software gratuito para resolver MIPs que encontrei foi o SCIP, que pode ser utilizado em Python com o módulo Numberjack. Mas fique atento, dependendo da formulação, um problema deste tamanho exigirá bastante memória e tempo para chegar à solução.

Cálculo da distância

Para calcularmos a distância entre dois pontos originalmente identificados pela latitude e longitude decimais, vamos transformá-los em radianos e usar a Aproximação Equiretangular. Como as distâncias entre os pontos, em geral, são pequenas, não teremos problemas com erros muito grandes, e esta acaba sendo uma fórmula bastante simples.

def length(town1, town2):
    R = 6371.
    x = (town2.long - town1.long) * math.cos((town2.lat + town1.lat)/2.)
    y = town2.lat - town1.lat
    return math.sqrt(x**2 + y**2) * R

def length(town1, town2):

R = 6371.

x = (town2.long - town1.long) * math.cos((town2.lat + town1.lat)/2.)

y = town2.lat - town1.lat

return math.sqrt(x**2 + y**2) * R

Você pode obter mais informações sobre esta fórmula, e outras maneiras de calcular esta distância neste link, que é a fonte que utilizei para decidir qual fórmula utilizar.

(Meta) Heurísticas

Greedy

def greedy(sol, townCount,town_list):
    
    for i in xrange(townCount-1):
        cost = length(town_list[sol[i]], town_list[sol[i+1]])
        for j in xrange(i+2,townCount):
            cost2 = length(town_list[sol[i]], town_list[sol[j]])
            
            if cost2 < cost:
                old = sol[i+1]
                sol[i+1] = sol[j]
                sol[j] = old
                cost = cost2

    return sol

def greedy(sol, townCount,town_list):

for i in xrange(townCount-1):

cost = length(town_list[sol[i]], town_list[sol[i+1]])

for j in xrange(i+2,townCount):

cost2 = length(town_list[sol[i]], town_list[sol[j]])

if cost2 < cost:

old = sol[i+1]

sol[i+1] = sol[j]

sol[j] = old

cost = cost2

return sol

A solução que normalmente nos vêm à cabeça primeiro é simplesmente selecionar a cidade com a menor distância como a próxima a ser visitada no percurso. Esta é uma solução rápida e em boa parte dos casos vai atingir um resultado satisfatório. Além disso, podemos utilizá-la como solução inicial para ajudar algoritmos mais complexos a chegar mais rápido a uma solução melhor.

No nosso caso, esta solução nos dá uma distância total de 156.760 km. Esta será a base de comparação das outras soluções.

Iterated Local Search 2-Opt

Outra maneira de aproximar a distância mínima é usar uma busca local. Este método, diferente do MIP, funciona melhor com problemas grandes, mas normalmente não nos dá a menor rota.

Basicamente consiste em encontrar soluções vizinhas à nossa atual. Por exemplo, trocando duas cidades de posição. Dentre as maneiras de encontrar soluções vizinhas, uma que funciona muito bem, e é bastante utilizada é a K-OPT.

Se visualizarmos o problema como uma rede (ou um grafo): cada cidade sendo um nó, conectado à cidade anterior e à próxima, e o valor desta conexão é a distância entre elas. Neste caso, K é o número de nós que vamos trocar de posição para ver se melhoramos a solução.

Neste exemplo será utilizado o 2-OPT, então buscaremos vizinhos trocando 2 nós. Em código isso é feito simplesmente selecionado duas cidades, e invertendo a ordem das cidades entre elas.

def selectNextTour2Opt(solution, townCount):
    first = randint(0,townCount)
    last = randint(first,townCount)
    
    slice = solution[first:last+1]
    slice.reverse()
    new_solution = solution[:first] + slice + solution[last+1:]

    return new_solution

def selectNextTour2Opt(solution, townCount):

first = randint(0,townCount)

last = randint(first,townCount)

slice = solution[first:last+1]

slice.reverse()

new_solution = solution[:first] + slice + solution[last+1:]

return new_solution

Num exemplo simples: imagine que as cidades A – B – C – D estão conectadas, e resolvemos aplicar o 2-OPT, selecionando as cidades B e C como início e fim do trajeto. Vamos inverter o trajeto B – C, que se tornará C – B, e teremos uma nova solução A – C – B – D. Se a distância do trajeto desta solução for menor do que a anterior, aceitamos, caso contrário, descartamos.

Na prática este procedimento desfaz trajetos que contenham rotas que cruzem sobre as outras, e por isso, são ineficientes. Em nosso exemplo vamos utilizar um 2-OPT aleatório devido ao tamanho do problema, mas uma outra maneira de decidir os elos que serão trocados é selecionar o primeiro elo que diminua a distância, verificando cada um deles, ou verificar exaustivamente e encontrar o que apresente a maior redução entre todos eles.

Simulated Annealing

Um dos problemas do método acima é que ele pode ficar parado em pontos mínimos locais. Existem vários métodos para “escapar” desses pontos, e buscar uma mínima melhor, mesmo sem garantias que ela seja global.

Uma destas maneiras é o Simulated Annealing. Baseado no processo de aquecimento e resfriamento da indústria de metalurgia, este procedimento possui dois parâmetros: temperatura e coeficiente de resfriamento. Basicamente, em vez de aceitar apenas trajetos que possuam uma distância menor, nós temos a probabilidade de aceitar um trajeto que piore, temporariamente, a distância, na busca por uma solução melhor no futuro.

    for i in xrange(1000000):

        new_solution = selectNextTour2Opt(solution,townCount)
        
        newCost = costSolution(new_solution,town_list,townCount)
        if newCost < obj:
            obj = newCost
            solution = new_solution[:]
        elif expCoinFlip(newCost,obj,T):
            solution = new_solution[:]
            obj = newCost


        if obj < min_obj:
            
            min_obj = obj
            print 'Epoch',i,'Objective', min_obj

        T = 0.999 * T



        if T <= final_T:
            T = init_T

for i in xrange(1000000):

new_solution = selectNextTour2Opt(solution,townCount)

newCost = costSolution(new_solution,town_list,townCount)

if newCost < obj:

obj = newCost

solution = new_solution[:]

elif expCoinFlip(newCost,obj,T):

solution = new_solution[:]

obj = newCost

if obj < min_obj:

min_obj = obj

print 'Epoch',i,'Objective', min_obj

T = 0.999 * T

if T <= final_T:

T = init_T

Esta probabilidade é proporcional à temperatura. Pense da seguinte maneira: no início, com a temperatura alta, nós estamos abertos a explorar várias soluções, mesmo que elas não melhorem nosso objetivo, conforme o tempo vai passando, decidimos caminhar em direção à mínima local de uma destas soluções.

No nosso caso, ainda teremos o reaquecimento, que é simplesmente mandar o algoritmo voltar à temperatura inicial após atingirmos uma temperatura muito baixa.

Resultados

Este algoritmo requer o ajuste de alguns parâmetros para que possamos atingir o melhor resultado possível. O valor ideal para estes parâmetros depende da estrutura do problema, e uma das maneiras mais simples de se encontrá-los é testando algumas combinações.

Não quero complicar muito, então os parâmetros ajustados, e os respectivos valores testados serão:

Temperatura inicial: 1, 10
Coeficiente de resfriamento: 0,99; 0,999; 0,9999

Vamos ver quais são os valores que nos ajudam a alcançar a distância mínima em 1 milhão de iterações. Assim como este problema, encontrar os parâmetros é, normalmente, um problema com pontos ótimos locais. Se tivéssemos muitos valores para visitar, seria recomendado fazer uma busca aleatória, mas como estamos restritos a 6 combinações, e a avaliação leva pouco tempo, vou testar todas elas.

No fim, a combinação com a menor distância em 1 milhão de iterações foi a temperatura inicial 1, e o coeficiente de resfriamento 0,999. Mas como a diferença entre as opções foi pequena, acredito que tenha sido uma questão de sorte. A distância final foi 151.125 km, uma redução aproximada de 5.635 km, ou 3,59%.

Decidi rodar o mesmo script, com os parâmetros anteriores, por 10 milhões de iterações. Após cerca de 40 minutos a distância atingida foi 142.794 km, uma redução aproximada de 13.966 km, ou 8,91%. Por se tratar de um problema com muitas cidades, e cada uma delas bem próximas de outras, acredito que rodar por mais tempo ainda possa diminuir bastante esse número.

Classificação, Regressão

As Métricas Mais Populares para Avaliar Modelos de Machine Learning

maio 25, 2015 por Mario Filho

Durante o processo de criação de um modelo de machine learning nós precisamos medir a qualidade dele de acordo com o objetivo da tarefa. Existem funções matemáticas que nos ajudam a avaliar a capacidade de erro e acerto dos nossos modelos, e agora você conhecerá algumas das mais utilizadas. No artigo, usarei a palavra métrica para me referir a essas funções.

Tão importante quanto saber escolher um bom modelo, é saber escolher a métrica correta para decidir qual é o melhor entre eles.

Existem métricas mais simples, outras mais complexas, algumas que funcionam melhor para datasets com determinadas características, ou outras personalizadas de acordo com o objetivo final do modelo.

Ao escolher uma métrica deve-se levar em consideração fatores como a proporção de dados de cada classe no dataset e o objetivo da previsão (probabilidade, binário, ranking, etc). Por isso é importante conhecer bem a métrica que será utilizada, já que isso pode fazer a diferença na prática.

Nenhuma destas funções é melhor do que as outras em todos os casos. É sempre importante levar em consideração a aplicação prática do modelo. O objetivo deste artigo não é ir a fundo em cada uma delas, mas apresentá-las para que você possa pesquisar mais sobre as que achar interessante.

Classificação

Estas métricas são utilizadas em tarefas de classificação, e a maioria delas pode ser adaptada tanto para classificação binária quanto de múltiplas classes. Nas tarefas de classificação buscamos prever qual é a categoria a que uma amostra pertence como, por exemplo, determinar se uma mensagem é spam.

Precisão Geral (Accuracy)

$Precis\tilde{a}o\ Geral = \frac{P}{P\ +\ N}$

Esta é a métrica mais simples. É basicamente o número de acertos (positivos) divido pelo número total de exemplos. Ela deve ser usada em datasets com a mesma proporção de exemplos para cada classe, e quando as penalidades de acerto e erro para cada classe forem as mesmas.

Em problemas com classes desproporcionais, ela causa uma falsa impressão de bom desempenho. Por exemplo, num dataset em que 80% dos exemplos pertençam a uma classe, só de classificar todos os exemplos naquela classe já se atinge uma precisão de 80%, mesmo que todos os exemplos da outra classe estejam classificados incorretamente.

F1 Score

$F1 = \frac{2\ *\ precis\tilde{a}o\ *\ recall}{precis\tilde{a}o\ +\ recall}$

O F1 Score é uma média harmônica entre precisão (que, apesar de ter o mesmo nome, não é a mesma citada acima) e recall. Veja abaixo as definições destes dois termos.

Ela é muito boa quando você possui um dataset com classes desproporcionais, e o seu modelo não emite probabilidades. Isso não significa que não possa ser usada com modelos que emitem probabilidades, tudo depende do objetivo de sua tarefa de machine learning.

Em geral, quanto maior o F1 score, melhor.

Precisão (Precision)

$Precis\tilde{a}o = \frac{PV}{PV\ +\ FP}$

Número de exemplos classificados como pertencentes a uma classe, que realmente são daquela classe (positivos verdadeiros), dividido pela soma entre este número, e o número de exemplos classificados nesta classe, mas que pertencem a outras (falsos positivos).

Recall

$Recall = \frac{PV}{P}$

Número de exemplos classificados como pertencentes a uma classe, que realmente são daquela classe, dividido pela quantidade total de exemplos que pertencem a esta classe, mesmo que sejam classificados em outra. No caso binário, positivos verdadeiros divididos por total de positivos.

AUC – Area Under the ROC Curve

auc — Imagem retirada da documentação do Scikit-Learn

Esta é uma métrica interessante para tarefas com classes desproporcionais. Nela, mede-se a área sob uma curva formada pelo gráfico entre a taxa de exemplos positivos, que realmente são positivos, e a taxa de falsos positivos.

Uma das vantagens em relação ao F1 Score, é que ela mede o desempenho do modelo em vários pontos de corte, não necessariamente atribuindo exemplos com probabilidade maior que 50% para a classe positiva, e menor, para a classe negativa.

Em sistemas que se interessam apenas pela classe, e não pela probabilidade, ela pode ser utilizada para definir o melhor ponto de corte para atribuir uma ou outra classe a um exemplo. Este ponto de corte normalmente é o ponto que se localiza mais à esquerda, e para o alto, no gráfico, mas depende bastante do custo do erro na previsão de uma determinada classe.

Log Loss

$Log\ Loss = -\frac{1}{N}\sum_{i=1}^N {(y_i\log(p_i) + (1 - y_i)\log(1 - p_i))}$

A fórmula do exemplo é para o caso binário, neste caso: p é a probabilidade do exemplo pertencer à classe 1, e y é o valor real da variável dependente.

Esta função pune previsões incorretas muito confiantes. Por exemplo, prever uma classe com uma probabilidade de 95%, e na realidade a correta ser outra. Ela pode ser utilizada para problemas binários ou com múltiplas classes, mas eu particularmente não gosto de usar em datasets com classes desproporcionais. O valor dela sempre terá a tendência de melhorar se o modelo estiver favorecendo a maior classe presente.

Tomando os cuidados acima, nas situações em que a probabilidade de um exemplo pertencer a uma classe for mais importante do que classificá-lo diretamente, esta função é preferível a usar simplesmente a precisão geral.

Regressão

Neste parte estão as funções mais comuns utilizadas para avaliar o desempenho de modelos de regressão. Na regressão buscamos prever um valor numérico, como, por exemplo, as vendas de uma empresa para o próximo mês. Nos exemplos abaixo:

$y_i = valor\ real\\\hat{y}_i = valor\ previsto$

Mean Squared Error – MSE

$MSE = \frac{1}{n} \sum_{i=1}^{n} (y_i - \hat{y}_i)^2$

Talvez seja a mais utilizada, esta função calcula a média dos erros do modelo ao quadrado. Ou seja, diferenças menores têm menos importância, enquanto diferenças maiores recebem mais peso.

Existe uma variação, que facilita a interpretação: o Root Mean Squared Error. Ele é simplesmente a raiz quadrada do primeiro. Neste caso, o erro volta a ter as unidades de medida originais da variável dependente.

Mean Absolute Error – MAE

$MAE = \frac{1}{n} \sum_{i=1}^{n} |y_i - \hat{y}_i|$

Bastante parecido com MSE, em vez de elevar a diferença entre a previsão do modelo, e o valor real, ao quadrado, ele toma o valor absoluto. Neste caso, em vez de atribuir um peso de acordo com a magnitude da diferença, ele atribui o mesmo peso a todas as diferenças, de maneira linear.

Se imaginarmos um exemplo simples, onde temos apenas a variável que estamos tentando prever, podemos ver um fato interessante que difere o MSE do MAE, e que devemos levar em conta ao decidir entre os dois: o valor que minimizaria o primeiro erro seria a média, já no segundo caso, a mediana.

Mean Absolute Percentage Error – MAPE

$MAPE = \frac{1}{n} \sum_{i=1}^{n} |\frac{y_i - \hat{y}_i}{y_i}|$

Este erro calcula a média percentual do desvio absoluto entre as previsões e a realidade. É utilizado para avaliar sistemas de previsões de vendas e outros sistemas nos quais a diferença percentual seja mais interpretável, ou mais importante, do que os valores absolutos.

Métricas específicas para tarefas

Em alguns casos, o ideal é usar uma métrica que tenha um significado específico para a tarefa em questão. Por exemplo, na segmentação de anúncios, a taxa de cliques; num sistema para comprar e vender ações, verificar o retorno médio. As métricas acima são importantes e podem ser utilizadas de maneira geral, mas se houver uma alternativa melhor, mais adequada ao contexto, ela deve ser utilizada.

MOOCs

Os Melhores Cursos Gratuitos Para Começar a Aprender Machine Learning e Data Science

maio 11, 2015 por Mario Filho

Data Science (ou ciência de dados) é uma área nova e muito promissora. Diversas fontes indicam que estamos precisando de mais data scientists do que conseguimos treinar. Neste artigo quero colocar cursos que você pode fazer para iniciar a sua jornada na área.

Todos eles são oferecidos por professores e universidades renomadas, e podem ser feitos gratuitamente pela internet. Além disso, possuem um conteúdo prático bastante forte, o que te dará condições para começar a aplicar as técnicas demonstradas em dados de seu interesse.

O foco destes cursos é em Machine Learning, que inclui boa parte das ferramentas mais utilizadas dentro da ciência de dados.

Machine Learning – Coursera

Este é o curso mais popular de introdução a machine learning online. Existem vários tipos de data scientists, mas boa parte do trabalho exige conhecimento de modelos que aprendem com dados. Neste curso é possível aprender superficialmente como funcionam diversos algoritmos, e como usá-los para resolver problemas reais.

O curso inicia explicando o que é Machine Learning e segue para a explicação de modelos simples, como a regressão linear, e vai construindo as bases para modelos mais complexos, e bastante utilizados, como Redes Neurais e SVMs. Além disso, a oportunidade de implementar partes destes algoritmos nos ajuda a entender melhor o funcionamento deles, e mesmo que você acabe nunca tendo que implementar um algoritmo destes do zero, isso vai te ajudar a conhecer suas características e limitações para usá-los em situações práticas.

Inicia com aprendizado supervisionado, e no fim faz uma breve excursão por tópicos de aprendizado não supervisionado e detecção de anomalias. Além disso, existem vídeos com técnicas e boas práticas para avaliar e otimizar um modelo, evitar erros como overfitting ou underfitting, e modificações para tornar possível o uso destes algoritmos em dados que não cabem na memória (o famoso big data).

O professor Andrew Ng é um dos fundadores do Coursera, além de professor da universidade de Stanford e cientista-chefe no Baidu.

É um curso bastante prático, e exige apenas familiaridade com álgebra linear. Dentro do curso existem vídeos que revisam conceitos que serão importantes. A linguagem utilizada para os trabalhos de programação é Octave. Existem vídeos dentro do curso que ensinam o básico sobre a linguagem, apenas o necessário para poder completar o curso.

Este curso é oferecido “on demand”, ou seja, pode ser iniciado a qualquer tempo.

Statistical Learning – StanfordX

Uma introdução ao Machine Learning do ponto de vista estatístico. Este curso oferecido por dois professores bastante respeitados na área de estatística, Rob Tibshrani e Trevor Hastie, é bastante prático e traz uma abordagem mais próxima dos conceitos estatísticos por trás de cada modelo, sem tanto foco na parte computacional.

Começa com uma visão geral da área de Statistical Learning, explica problemas de classificação e regressão e ferramentas básicas para modelagem linear. Depois nos leva por métodos de avaliação de modelos, e técnicas para otimizar o modelo levando em conta a generalização para novos dados. Por fim são apresentados algoritmos mais avançados, como SVM e Random Forests, além de uma breve passagem por métodos de aprendizado não supervisionado.

Este curso utiliza a linguagem R, que é largamente utilizada na área de estatística. Existem trabalhos de programação, mas são voltados para a utilização dos modelos, e não para a implementação.

Este curso normalmente é oferecido uma vez por ano, em meados de Janeiro.

Learning from Data – EdX

Este curso, oferecido pelo professor Yaser Abu-Mostafa aplica uma abordagem mais computacional mas, diferente do curso do Andrew Ng, existe uma boa dose de teoria, que acaba ajudando a entender o funcionamento dos modelos de forma mais completa.

Nas primeiras aulas o professor nos explica sobre o conceito de aprendizado de máquina, e a matemática que fundamenta a teoria que nos garante a possibilidade de usar um algoritmo para aprender uma tarefa através de dados. A teoria é apresentada de dois pontos de vista: VC Dimension, e Bias-Variance Tradeoff.

Depois disso, ele apresenta alguns modelos como regressão logística, redes neurais e SVMs. São demonstradas algumas técnicas para otimizar os modelos de maneira que sejam úteis numa aplicação prática.

Por fim, somos apresentados a Kernels, que são transformações de variáveis, bastante importantes, principalmente pelo sucesso dos SVMs, e uma visão geral sobre áreas de estudo em Machine Learning que podem ser seguidas após o curso.

Apesar de não haver uma confirmação sobre uma próxima sessão, vale a pena citá-lo e todos os materiais (vídeos e tarefas) estão disponíveis online. Ele não exige uma linguagem de programação específica, é possível completar as tarefas em qualquer linguagem.

Bônus – The Analytics Edge – EdX

Este curso, oferecido pelo MIT e pelo professor Dimitris Bertsimas e sua equipe através da plataforma EdX foca bastante na aplicação de métodos de machine learning usando R. Boa parte do curso é passada em exemplos de uso das técnicas, e as tarefas são extensivas, para que os alunos possam explorar todos os comandos ensinados. Toda semana temos um novo estudo de caso.

Além de ensinar métodos de aprendizado supervisionado, e não supervisionado, no fim eles falam sobre métodos de otimização que, além de interessantes por si só, são a base de muitos algoritmos de machine learning.

Um ponto que chama a atenção neste curso e que, no momento em que escrevo este artigo, não é oferecido nos outros, é que durante o curso uma das tarefas é participar de uma competição, oferecida apenas para alunos do curso, na plataforma Kaggle. Esta é uma oportunidade de usar as ferramentas num caso próximo da realidade, tendo que criar uma solução usando o conhecimento adquirido no curso.

Esta parte da competição é de extrema importância. Por experiência própria eu digo que nada ensina mais do que ver um dataset na sua frente e ter que decidir sozinho qual é a melhor direção a se tomar para fazer a análise.

Como escolher um curso?

Se você me perguntar qual deles deve fazer, eu vou responder “todos”. Apesar de boa parte do material ser o mesmo, cada um deles te dá uma visão diferente do Machine Learning. Com Andrew Ng você tem uma apresentação rápida e superficial, bastante prática, dos algoritmos. Já no Statistical Learning, apesar de bem prático, há uma preocupação maior com conceitos estatísticos clássicos, como p-value e intervalos de confiança. No curso da Caltech, Learning from Data, é possível entender a teoria que fundamenta o Machine Learning, a razão matemática para um algoritmo conseguir aprender através dos dados.

Se você estiver disposto a fazer todos, eu recomendo fazê-los na ordem que eles estão dispostos no artigo.

Big Data, Classificação, Online Learning, Python

Como Criar um SVM para Big Data em 10 Minutos Usando Python

abr 27, 2015 por Mario Filho

SVMs são bastante populares em aplicações de machine learning. Apesar de seguirem princípios simples, eles já demonstraram ser muito poderosos em diversas tarefas e datasets. Neste artigo quero demonstrar como implementar um SVM capaz de lidar com dados que chegam em tempo real, sem ter que armazená-los na memória.

Clique aqui para acessar o código e o dataset deste artigo. Eu recomendo rodá-lo usando Pypy.

Principais modos de treinar um modelo de machine learning

Existem três maneiras populares de treinar um modelo: batch learning, mini-batch learning, e stochastic learning.

Batch Learning: no primeiro modo, nós armazenamos todos os dados de treino numa matriz e alimentamos o algoritmo, reduzindo a loss function baseado em todos os exemplos de uma vez. Isso nem sempre é possível devido ao tamanho do dataset. Nestes casos temos que recorrer às duas outras maneiras.

Mini-Batch Learning: neste caso, selecionamos um número N de exemplos e dividimos o training set em blocos deste tamanho. Então treinamos o modelo em um bloco de cada vez.

Stochastic learning: é uma variação do mini-batch, mas com N = 1. Nele utilizamos apenas um exemplo de cada vez para treinar o modelo. Este é o modo preferido em soluções de big data.

Dentro deste modo temos duas divisões principais:

Offline: neste caso temos todos os dados armazenados num banco de dados, mas eles não cabem na memória, então precisamos treinar em um exemplo de cada vez.

Online: neste caso recebemos os exemplos em tempo real, treinamos e descartamos, sem a necessidade de armazená-los.

Além de não ter que carregar todos os dados na memória para treinar o modelo, em alguns casos, usar algoritmos que treinam em um exemplo de cada vez pode ser mais rápido do que o formato que utiliza todos os dados ao mesmo tempo.

Dataset

Neste artigo será usado o Skin Segmentation Data Set que pode ser baixado do banco de dados da UCI. Estes dados se referem a uma tarefa de classificação na qual são extraídos valores RGB de pontos aleatórios da foto do rosto de uma pessoa, e a tarefa é determinar, de acordo com estes valores, se aquele ponto corresponde à pele, ou outra parte da imagem.

Segundo os pesquisadores, foram coletados valores de imagens de diversas idades, raças e gêneros. Uma aplicação prática desta tarefa seria identificar imagens com conteúdo impróprio para menores na internet.

São cerca de 245 mil exemplos, e é um problema de classificação binária.

Algoritmo de treino

O algoritmo utilizado para treino será o PEGASOS. O trabalho acadêmico (em inglês) com os detalhes pode ser encontrado aqui.

Na sua versão mais simples, vamos minimizar a função primal do SVM utilizando um exemplo de cada vez. Como este método se aplica ao subgradiente, ele não garante a diminuição do erro em todas as iterações. Mesmo assim, existem garantias de convergência.

O único parâmetro a escolher é a constante que mede o grau de regularização. Normalmente ela é escolhida usando algum método de validação. Como o objetivo do artigo é demonstrar uma implementação de SVM que pode ser facilmente adaptada para dados que não cabem na memória, ou são recebidos em tempo real, vou fixa-la em 0,1.

Numa aplicação real, pode ser interessante ajustar esse parâmetro usando uma amostra menor dos dados, off-line, com validação cruzada.

Avaliação do resultado

Num ambiente real online, uma maneira de avaliar a performance do modelo seria calcular a loss antes de prever um novo exemplo e fazer a média numa janela de tempo. Este método é chamado de validação progressiva.

No nosso caso, estamos simulando um ambiente online, e para tornar mais clara a eficácia do modelo, vou usar um test set separado, com 20% dos dados.

No repositório deste artigo você encontra os arquivos train.csv e test.csv, já preparados. Cada um contém aproximadamente 25% de exemplos positivos, então além da loss, quero avaliar o número de acertos e erros em cada classe (positivos verdadeiros e negativos verdadeiros).

O training set possui 196.045 exemplos, e o test set, 49.012. Vamos utilizar cada exemplo do training set apenas uma vez, e depois avaliar a performance do modelo no test set.

Como nós não vamos utilizar este test set para ajustar hiperparâmetros, ele é uma estimativa legítima da performance do modelo em novos exemplos. Se você for usá-lo para validar parâmetros, precisará de outro test set.

Implementation

We start by creating a function to read each line of the CSV and transform it into a list with the values of the variables, and the bias (always 1). The function is divided between the training set and test set.

def data(self,test=False): if test: with open('test.csv','r') as f: samples = f.readlines() for t,row in enumerate(samples): row = row.replace('\n','') row = row.split(',') target = -1. if row[3] == '1': target = 1. del row[3] x = [1.] + [float(c) for c in row] # bias + inputs yield t, x,target else: with open('train.csv','r') as f: samples = f.readlines() shuffle(samples) for t,row in enumerate(samples): row = row.replace('\n','') row = row.split(',') target = -1. if row[3] == '1': target = 1. del row[3] x = [1.] + [float(c) for c in row] # bias + inputs yield t, x,target

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45

 def data(self,test=False):

 if test:
 with open('test.csv','r') as f:
 samples = f.readlines()

 for t,row in enumerate(samples):

 row = row.replace('\n','')
 row = row.split(',')

 target = -1.

 if row[3] == '1':
 target = 1.
 del row[3]


 x = [1.] + [float(c) for c in row] # bias + inputs

 yield t, x,target

 else:

 with open('train.csv','r') as f:
 samples = f.readlines()
 shuffle(samples)

 for t,row in enumerate(samples):


 row = row.replace('\n','')
 row = row.split(',')

 target = -1.

 if row[3] == '1':
 target = 1.
 del row[3]



 x = [1.] + [float(c) for c in row] # bias + inputs

 yield t, x,target

Depois vamos à implementação do algoritmo para treinar o SVM. Criamos algumas funções de apoio, uma para calcular a hinge loss, que é a loss que estamos minimizando, e outra para determinar o sinal da output do SVM. Os pesos iniciais do SVM podem ser 0, já que não temos problemas com simetria, como em redes neurais.

def init(self,lmbd,D): self.lmbd = lmbd # coeficiente de regularização self.D = D + 1 # dimensão dos atributos mais bias self.w = [0.] * self.D # pesos do modelo # Função que atribui a classe do exemplo... def sign(self, x): return -1. if x <= 0 else 1. #Hinge Loss, a função que será minimizada def hinge_loss(self,target,y): return max(0, 1 - targety)

1
2
3
4
5
6
7
8
9
10
11
12

 def init(self,lmbd,D):
 self.lmbd = lmbd # coeficiente de regularização
 self.D = D + 1 # dimensão dos atributos mais bias
 self.w = [0.] self.D # pesos do modelo

 # Função que atribui a classe do exemplo...
 def sign(self, x):
 return -1. if x <= 0 else 1.

 #Hinge Loss, a função que será minimizada
 def hinge_loss(self,target,y):
 return max(0, 1 - targety)

Depois criamos uma função para gerar as previsões. Ela simplesmente fará o produto escalar entre o vetor de pesos, e o vetor de atributos.

#Calcula o produto escalar entre os pesos e os atributos do exemplo def predict(self,x): wTx = 0. for i in xrange(len(x)): wTx += self.w[i]x[i] return wTx

1
2
3
4
5
6
7
8

 #Calcula o produto escalar entre os pesos e os atributos do exemplo
 def predict(self,x):
 wTx = 0.
 for i in xrange(len(x)):

 wTx += self.w[i]x[i]

 return wTx

Agora vamos à função para que o modelo aprenda os dados. Nela vamos atualizar os pesos com a learning rate (alpha) e o gradiente. Esta função treina em um exemplo de cada vez. A learning rate se ajusta, sendo menor a cada novo exemplo visto.

O gradiente é diferente para o caso em que a classificação atual está errada (o exemplo está do lado errado do hiperplano), e para os exemplos que estão classificados corretamente.

Seria possível implementar esta função usando menos linhas, mas neste caso optei pela clareza. Para verificar esta possibilidade, veja a fórmula na publicação original deste algoritmo.

# Treina em um exemplo def train(self,x,y,alpha): #se y estiver incorreto if yself.predict(x) < 1: for i in xrange(len(x)): self.w[i] = (1. - alphaself.lmbd)self.w[i] + alphayx[i] else: for i in xrange(len(x)): self.w[i] = (1. - alphaself.lmbd)self.w[i]

1
2
3
4
5
6
7
8
9
10
11
12

 # Treina em um exemplo
 def train(self,x,y,alpha):

 #se y estiver incorreto
 if yself.predict(x) < 1:

 for i in xrange(len(x)):
 self.w[i] = (1. - alphaself.lmbd)self.w[i] + alphay*x[i]

 else:
 for i in xrange(len(x)):
 self.w[i] = (1. - alphaself.lmbd)self.w[i]

Agora criamos a função para treinar e avaliar o modelo nos dados. No caso de online learning, receberíamos os dados através de uma stream, então eles seriam vistos apenas uma vez.

Para cada exemplo do training set invocamos a função ‘train’, que ajusta os pesos, de acordo com a learning rate, na direção do gradiente. Lembrando que, por ser um método que utiliza o subgradiente, ele não garante a diminuição da loss em todas as iterações.

def fit(self): ... #last = 0 for t, x,target in self.data(test=False): #if target == last: # continue alpha = 1./(self.lmbd(t+1.)) self.train(x,target,alpha) #last = target for t,x,target in self.data(test=True): pred = self.predict(x) loss += self.hinge_loss(target,pred) pred = self.sign(pred) ...

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22

 def fit(self):

 ...

 #last = 0
 for t, x,target in self.data(test=False):

 #if target == last:
 # continue

 alpha = 1./(self.lmbd(t+1.))
 self.train(x,target,alpha)
 #last = target

 for t,x,target in self.data(test=True):

 pred = self.predict(x)
 loss += self.hinge_loss(target,pred)

 pred = self.sign(pred)

 ...

Consideramos um exemplo positivo se a distância dele para o hiperplano for positiva, e o contrário para exemplos negativos. É extremamente raro que um exemplo seja absolutamente 0, mas vamos considerar que, se isso acontecer, ele será classificado como negativo.

Além disso, devido à natureza estocástica do algoritmo, vamos treiná-lo várias vezes, e fazer a média da performance modificando a ordem de chegada dos exemplos. Para que os resultados possam ser reproduzidos, determinamos uma seed para o RNG.

for i in range(100): seed(i) svm = SVMPegasos(1e-1,3) l,tp,tn = svm.fit()

1
2
3
4
5
6
7

for i in range(100):

 seed(i)

 svm = SVMPegasos(1e-1,3)

 l,tp,tn = svm.fit()

Influência do Training Set

Por ser um algoritmo que considera um exemplo de cada vez, a ordem na qual eles são apresentados muda o comportamento, fazendo com que ele busque caminhos diferentes para atingir a mínima. Todos os números reportados nesta parte são avaliados nos exemplos do arquivo “test.csv”.

Treinando em nosso training set original, sem embaralharmos os dados, a performance dele é a seguinte:

Loss 0.857527620114 
Negativos verdadeiros 95.43 % 
Positivos verdadeiros 6.66 % 
Precisão Geral 77.05 %

Loss 0.857527620114

Negativos verdadeiros 95.43 %

Positivos verdadeiros 6.66 %

Precisão Geral 77.05 %

Agora, para termos uma estimativa real da performance na prática, vamos testá-lo 100 vezes, modificando a ordem de chegada dos exemplos, a performance média neste caso é:

Loss Média 0.589555758345
Positivos verdadeiros 63.86 %
Negativos verdadeiros 92.20 %

Loss Média 0.589555758345

Positivos verdadeiros 63.86 %

Negativos verdadeiros 92.20 %

Um detalhe deste dataset é o fato que a proporção das classes é desequilibrada (75/25). Um remédio simples para isto, no caso de termos muitos dados, é diminuir a proporção da classe maior até atingir o equilíbrio (50/50). Além disso, alguns especialistas sugerem que garantir que o próximo exemplo visto pelo algoritmo seja diferente do anterior possa melhorar o desempenho. Vamos testar essas duas opções adicionando o seguinte código:

        <strong>last = 0</strong>
        for t, x,target in self.data(test=False):
             <strong>if target == last:
                  continue</strong>
            
            alpha = 1./(self.lmbd*(t+1.))
            self.train(x,target,alpha)
            <strong>last = target</strong>

last = 0

for t, x,target in self.data(test=False):

if target == last:

continue

alpha = 1./(self.lmbd*(t+1.))

self.train(x,target,alpha)

last = target

Ele fará com que o algoritmo sempre receba um exemplo diferente do anterior, e que a proporção das classes seja aproximadamente igual. Não vamos mexer no test set, ele continua com a proporção original (75/25).

Neste caso, usando apenas a ordem do training set original, a performance se torna:

Loss 70.5074554527 
Negativos verdadeiros 91.55 %
Positivos verdadeiros 95.87 % 
Precisão Geral 92.45 %

Loss 70.5074554527

Negativos verdadeiros 91.55 %

Positivos verdadeiros 95.87 %

Precisão Geral 92.45 %

E fazendo a média de 100 datasets embaralhados:

Loss Média 6.59906602428
Positivos verdadeiros 94.16 %
Negativos verdadeiros 90.48 %

Loss Média 6.59906602428

Positivos verdadeiros 94.16 %

Negativos verdadeiros 90.48 %

Vemos uma melhora significativa na quantidade de acertos dos positivos verdadeiros, apesar de um grande aumento na loss. Esse aumento se deve ao fato da loss ser calculada baseando-se em todos os exemplos. Como estamos “puxando” o hiperplano que separa as classes de maneira a aceitar mais erros na classe negativa, em troca de melhor performance na positiva, acabamos aumentando a distância de alguns pontos classificados de maneira incorreta.

Uma alternativa seria manter todos os dados e atribuir uma loss maior, proporcional, a erros em exemplos da classe menor.

Em casos que apresentam classes desproporcionais, como este, é interessante usar métricas que não sejam sensíveis à quantidade de exemplos. Por isso a escolha em reportar os positivos/negativos verdadeiros.

Resultado e Características

O primal do SVM é uma função convexa, então temos a possibilidade de atingir o ponto mínimo global. Entretanto, o fato de recebermos exemplos de maneira aleatória, treinando apenas uma vez em cada um, com uma learning rate diminuindo com o tempo, pode nos impedir de chegar a este ponto numa aplicação prática. Apesar disso, há garantias teóricas de convergência para a mínima global.

Em Machine Learning buscamos minimizar o erro de generalização, ou seja, fora dos exemplos utilizados para treino. Na maioria dos casos, a melhor performance de generalização ocorre antes de atingirmos a mínima no training set.