Churn-at-Risk: Aplicação de Survival Analysis no controle de churn de assinaturas em Telecom

Introdução

Um dos assuntos mais recorrentes em qualquer tipo de serviço de assinatura é como reduzir o Churn (saída de clientes), dado que conquistar novos clientes é bem mais difícil (e caro) do que manter os antigos.

Cerca de 70% das empresas sabem que é mais barato manter um cliente do que ter que ir atrás de um novo.

Fazendo uma analogia simples, o lucro dos serviços de assinatura são como uma espécie de sangue na corrente sanguínea de uma empresa e uma interrupção de qualquer natureza prejudica todo o negócio, dado que esse é um modelo de receita que se baseia na recorrência de tarifação e não no desenvolvimento, ou mesmo venda de outros produtos.

Em modelos de negócios baseados no volume de pessoas que estão dispostas a terem uma cobrança recorrente o negócio fica bem mais complicado, dado que diferentemente de produtos que tem uma elasticidade maior o fluxo de receita é extremamente sujeito aos sabores do mercado e dos clientes.

Dentro desse cenário, para todas as empresas que tem o seu fluxo de receita baseado nesse tipo de business, saber quando um cliente entrará em uma situação de saída através do cancelamento do serviço (Churn) é fundamental para criar mecanismos de retenção mais efetivos, ou mesmo criação de réguas de contato com os clientes para evitar ou minimizar a chance de um cliente sair da base de dados.

Sendo assim, qualquer mecanismo ou mesmo esforço para minimizar esse efeito é de grande valia. Nos baseamos na teoria estatística buscar respostas para as seguintes perguntas:

  • Como diminuir o Churn?
  • Como identificar um potencial cliente que irá entrar em uma situação de Churn? Quais estratégias seguir para minimizar esse Churn?
  • Quais réguas de comunicação com os clientes devemos ter para entender os motivos que estão fazendo um assinante cancelar o serviço e quais são as estratégias de customer winback possíveis nesse cenário?

E pra responder essa pergunta, fomos buscar as respostas na análise de sobrevivência dado que essa área da estatística é uma das que lidam melhor em termos de probabilidade de tempo de vida com dados censurados, seja de materiais (e.g. tempo de falha de algum sistema mecânico) ou no tempo de vida de pessoas propriamente ditas (e.g. dado uma determinada posologia qual é a estimativa de um paciente sobreviver a um câncer), e no nosso caso quanto tempo de vida um assinante tem até deixar cancelar a sua assinatura.

Análise de Sobrevivência

A análise de sobreviência é uma técnica estatístisca que foi desenvolvida na medicina e tem como principal finalidade estimar o tempo de sobrevivência ou tempo de morte de um determinado paciente dentro de um horizonte do tempo.

O estimador de Kaplan-Meier (1958) utiliza uma função de sobrevivência que leva em consideração uma divisão entre o número de observações que não falharam no tempo t pelo número total de observações no estudo em que cada intervalo de tempo tem-se o número de falhas/mortes/churn distintos bem como é calculado o risco de acordo com o número de indivíduos restantes no tempo subsequente.

Já o estimador Nelson-Aalen (1978) é um estimador que tem as mesmas características do Kaplan-Meier, com a diferença que esse estimador trabalha com uma função de sobrevivência que é a cumulative hazard rate function.

Os elementos fundamentais para caracterização de um estudo que envolve análise de sobrevivência são, o (a) tempo inicial, (b)escala de medida do intervalo de tempo e (c) se o evento de churn ocorreu.

Os principais artigos são de Aalen (1978), Kaplan-Meier (1958) e Cox (1972).

Esse post não tem como principal objetivo dar algum tipo de introdução sobre survival analysis, dado que tem muitas referências na internet sobre o assunto e não há nada a ser acrescentado nesse sentido por este pobre blogueiro.

Assim como a análise de cohort, a análise de sobrevivência tem como principal característica ser um estudo de natureza longitudinal, isto é, os seus resultados tem uma característica de temporalidade seja em aspectos de retrospecção, quanto em termos de perspectivas, isso é, tem uma resposta tipicamente temporal para um determinado evento de interesse.

O que vamos usar como forma de comparação amostral é o comportamento longitudinal__, de acordo com determinadas características de amostragens diferentes ao longo do tempo, e os fatores que influenciam no churn.

Devido a questões óbvias de NDA não vamos postar aqui características que possam indicar qualquer estratégia de negócios ou mesmo caracterização de alguma informação de qualquer natureza.

Podemos dizer que a análise de sobrevivência aplicada em um caso de telecom, pode ajudar ter uma estimativa em forma de probabilidade em relação ao tempo em que uma assinatura vai durar até o evento de churn (cancelamento) e dessa forma elaborar estratégias para evitar esse evento, dado que adquirir um novo cliente é mais caro do que manter um novo e entra totalmente dentro de uma estratégia de Customer Winback (Nota: Esse livro Customer Winback do Jill Griffin e do Michael Lowenstein é obrigatório para todos que trabalham com serviços de assinaturas ou negócios que dependam de uma recorrência muito grande como comércio).

No nosso caso o tempo de falha ou tempo de morte, como estamos falando de serviços de assinaturas, o nosso evento de interesse seria o churn, ou cancelamento da assinatura. Em outras palavras teríamos algo do tipo Time-to-Churn ou um Churn-at-Risk. Guardem esse termo.

Metodologia

Usamos dados de dois produtos antigos em que os dados foram anonimizados e aplicados um hash de embaralhamento uniforme (que obedece uma distribuição específica) nos atributos (por questões de privacidade) que são:

  • id = Identificador do registro;
  • product = produto;
  • channel = canal no qual o cliente entrou na base de dados;
  • free_user = flag que indica se o cliente entrou na base em gratuidade ou não;
  • user_plan = se o usuário é pré-pago ou pós-pago;
  • t = tempo que o assinante está na base de dados; e
  • c = informa se o evento de interesse (no caso o churn (cancelamento da assinatura) ocorreu ou não.

Eliminamos o efeito de censura à esquerda retirando os casos de reativações, dado que queríamos entender a jornada do assinante como um todo sem nenhum tipo de viés relativo a questões de customer winback. Em relação à censura à direita temos alguns casos bem específicos que já se passaram alguns meses desde que essa base de dados foi extraída.

Um aspecto técnico importante a ser considerado é que esses dois produtos estão em categorias de comparabilidade, dado que sem isso nenhum tipo de caractericação seria nula.

No fim dessa implementação teremos uma tabela de vida em relação a esses produtos.

Implementação

Primeiramente vamos importar as bibliotecas: Pandas (para manipulação de dados), matplotlib (para a geração de gráficos), e lifelines para aplicação da análise de sobrevivência:

%matplotlib inline import matplotlib.pyplot as plt import numpy as np import pandas as pd import lifelines

Após realizar a importação das bibliotecas, vamos ajustar o tamanho das imagens para uma melhor visualização:

%pylab inline pylab.rcParams[‘figure.figsize’] = (14, 9)

Vamos realizar o upload da nossa base de dados criando um objecto chamado df e usando a classe read_csv do Pandas:

df = pd.read_csv(‘https://raw.githubusercontent.com/fclesio/learning-space/master/Datasets/07%20-%20Survival/survival_data.csv’)

Vamos checar a nossa base de dados:

df.head()

id

product

channel

free_user

user_plan

t

c

0

3315

B

HH

1

0

22

0

1

2372

A

FF

1

1

16

0

2

1098

B

HH

1

1

22

0

3

2758

B

HH

1

1

4

1

4

2377

A

FF

1

1

29

0

Então como podemos ver temos as 7 variáveis na nossa base de dados.

Na sequência vamos importar a biblioteca do Lifelines, em especial o estimador de KaplanMaier:

from lifelines import KaplanMeierFitter kmf = KaplanMeierFitter()

Após realizar a importação da classe relativa ao estimador de Kaplan Meier no objeto kmf, vamos atribuir as nossas variáveis de tempo (T) e evento de interesse (C)

T = df[“t”] C = df[“c”]

O que foi feito anteriormente é que buscamos no dataframe df o array t e atribuímos no objeto T, e buscamos o array da coluna c no dataframe e atribuímos no objeto C. Agora vamos chamar o método fit usando esses dois objetos no snippet abaixo:

kmf.fit(T, event_observed=C )

Out[7]:

<lifelines.KaplanMeierFitter: fitted with 10000 observations, 6000 censored>

Objeto ajustado, vamos agora ver o gráfico relativo a esse objeto usando o estimador de Kaplan Meier.

kmf.survival_function_.plot() plt.title(‘Survival function of Service Valued Add Products’); plt.ylabel(‘Probability of Living (%)’) plt.xlabel(‘Lifespan of the subscription (in days)’)

Out[8]:

<matplotlib.text.Text at 0x101b24a90>

 1

Como podemos ver no gráfico, temos algumas observações pertinentes, quando tratamos a probabilidade de sobrevivência desses dois produtos no agregado que são:

  • Logo no primeiro dia há uma redução substancial do tempo de sobrevivência da assinatura em aproximadamente 22%;
  • Há um decaimento quase que linear depois do quinto dia de assinatura; e
  • Depois do dia número 30, a probabilidade de sobrevivência de uma assinatura é de aproximadamente de 50%. Em outras palavras: depois de 30 dias, metade dos novos assinantes já estarão fora da base de assinantes.

No entanto, vamos plotar a mesma função de sobrevivência considerando os intervalos de confiança estatística.

kmf.plot() plt.title(‘Survival function of Service Valued Add Products - Confidence Interval in 85-95%’); plt.ylabel(‘Probability of Living (%)’) plt.xlabel(‘Lifespan of the subscription’)

Out[9]:

<matplotlib.text.Text at 0x10ad8e0f0>

 2

Contudo nesse modelo inicial temos duas limitações claras que são:

  • Os dados no agregado não dizem muito em relação à dinâmicas que podem estar na especificidade de alguns atributos/dimensões;
  • Não são exploradas as dimensões (ou quebras) de acordo com os atributos que vieram na base de dados; e
  • Não há a divisão por produto.

Para isso, vamos começar a entrar no detalhe em relação a cada uma das dimensões e ver o que cada uma tem de influência em relação à função de sobrevivência. Vamos começar realizando a quebra pela dimensão que determina se o cliente entrou via gratuidade ou não (free_user).

ax = plt.subplot(111) free = (df[“free_user”] == 1) kmf.fit(T[free], event_observed=C[free], label=“Free Users”) kmf.plot(ax=ax, ci_force_lines=True) kmf.fit(T[~free], event_observed=C[~free], label=“Non-Free Users”) kmf.plot(ax=ax, ci_force_lines=True) plt.ylim(0,1); plt.title(“Lifespans of different subscription types”); plt.ylabel(‘Probability of Living (%)’) plt.xlabel(‘Lifespan’)

Out[10]:

<matplotlib.text.Text at 0x10ad8e908>

 3

Este gráfico apresenta algumas informações importantes para os primeiros insights em relação a cada uma das curvas de sobrevivência em relação ao tipo de gratuidade oferecida como fator de influência para o churn que são:

  • Os assinantes que entram como não gratuitos (i.e. não tem nenhum tipo de gratuidade inicial) após o 15o dia apresenta um decaimento brutal de mais de 40% da chance de sobrevivência (tratando-se do intervalo de confiança);
  • Após o 15o dia os assinantes que não desfrutam de gratuidade tem a sua curva de sobrevivência em uma relativa estabilidade em torno de 60% na probabilidade de sobrevivência até o período censurado;
  • Ainda nos usuários sem gratuidade, dado o grau de variabilidade do intervalo de confiança podemos tirar como conclusão que muitos cancelamentos estão ocorrendo de forma muito acelerada, o que deve ser investigado com mais calma pelo time de produtos; e
  • Já os usuários que entram via gratuidade (i.e. ganham alguns dias grátis antes de serem tarifados) apresenta um nível de decaimento do nível de sobrevivência maior seja no período inicial, quando ao longo do tempo, contudo uma estabilidade é encontrada ao longo de toda a série sem maiores sobressaltos.

Dado essa análise inicial das curvas de sobrevivência, vamos avaliar agora as probabilidades de sobrevivência de acordo com o produto.

ax = plt.subplot(111) product = (df[“product”] == “A”) kmf.fit(T[product], event_observed=C[product], label=“Product A”) kmf.plot(ax=ax, ci_force_lines=True) kmf.fit(T[~product], event_observed=C[~product], label=“Product B”) kmf.plot(ax=ax, ci_force_lines=True) plt.ylim(0,1); plt.title(“Survival Curves of different Products”); plt.ylabel(‘Probability of Living (%)’) plt.xlabel(‘Lifespan’)

Out[11]:

<matplotlib.text.Text at 0x10aeaabe0>

 4

Este gráfico apresenta a primeira distinção entre os dois produtos de uma forma mais clara. Mesmo com os intervalos de confiança com uma variação de 5%, podemos ver que o produto A (linha azul) tem uma maior probabilidade de sobrevivência com uma diferença percentual de mais de 15%; diferença essa amplificada depois do vigésimo dia. Em outras palavras: Dado um determinada safra de usuários, caso o usuário entre no produto A o mesmo tem uma probabilidade de retenção de cerca de 15% em relação a um usuário que por ventura entre no produto B, ou o produto A apresenta uma cauda de retenção superior ao produto B. Empiricamente é sabido que um dos principais fatores de influência de produtos SVA são os canais de mídia os quais esses produtos são oferecidos. O canal de mídia é o termômetro em que podemos saber se estamos oferencendo os nossos produtos para o público alvo correto. No entanto para um melhor entendimento, vamos analisar os canais nos quais as assinaturas são originadas. A priori vamos normalizar a variável channel para realizar a segmentação dos canais de acordo com o conjunto de dados.

df[‘channel’] = df[‘channel’].astype(‘category’); channels = df[‘channel’].unique()

Após normalização e transformação da variável para o tipo categórico, vamos ver como está o array.

channels

Out[13]:

[HH, FF, CC, AA, GG, …, BB, EE, DD, JJ, ZZ] Length: 11 Categories (11, object): [HH, FF, CC, AA, …, EE, DD, JJ, ZZ]

Aqui temos a representação de 11 canais de mídia os quais os clientes entraram no serviço. Com esses canais, vamos identificar a probabilidade de sobrevivência de acordo com o canal.

for i,channel_type in enumerate(channels): ax = plt.subplot(3,4,i+1) ix = df[‘channel’] == channel_type kmf.fit( T[ix], C[ix], label=channel_type ) kmf.plot(ax=ax, legend=True) plt.title(channel_type) plt.xlim(0,40) if i==0: plt.ylabel(‘Probability of Survival by Channel (%)’) plt.tight_layout()

5

Fazendo uma análise sobre cada um desses gráficos temos algumas considerações sobre cada um dos canais:

  • HH, DD: Uma alta taxa de mortalidade (churn) logo antes dos primeiros 5 dias, o que indica uma característica de efemeridade ou atratividade no produto para o público desse canal de mídia.

  • FF: Apresenta menos de 10% de taxa de mortalidade nos primeiros 20 dias, e tem um padrão muito particular depois do 25o dia em que praticamente não tem uma mortalidade tão alta. Contém um intervalo de confiança com uma oscilação muito forte.

  • CC: Junto com o HH apesar de ter uma taxa de mortalidade alta antes do 10o dia, apresenta um grau de previsibilidade muito bom, o que pode ser utilizado em estratégias de incentivos de mídia que tenham que ter uma segurança maior em termos de retenção a médio prazo.

  • GG, BB: Apresentam uma boa taxa de sobrevivência no inicio do período, contudo possuem oscilações severas em seus respectivos intervalos de confiança. Essa variável deve ser considerada no momento de elaboração de uma estratégia de investimento nesses canais.

  • JJ: Se houvesse uma definição de incerteza em termos de sobrevivência, esse canal seria o seu melhor representante. Com os seus intervalos de confiança oscilando em mais de 40% em relação ao limite inferior e superior, esse canal de mídia mostra-se extremamente arriscado para os investimentos, dado que não há nenhum tipo de regularidade/previsibilidade de acordo com esses dados.

  • II: Apesar de ter um bom grau de previsibilidade em relação à taxa de sobrevivência nos primeiros 10 dias, após esse período tem uma curva de hazard muito severa, o que indica que esse tipo de canal pode ser usado em uma estratégia de curto prazo.

  • AA, EE, ZZ: Por haver alguma forma de censura nos dados, necessitam de mais análise nesse primeiro momento. (Entrar no detalhe dos dados e ver se é censura à direita ou algum tipo de truncamento).

Agora que já sabemos um pouco da dinâmica de cada canal, vamos criar uma tabela de vida para esses dados. A tabela de vida nada mais é do que uma representação da função de sobrevivência de forma tabular em relação aos dias de sobrevivência. Para isso vamos usar a biblioteca utils do lifelines para chegarmos nesse valor.

from lifelines.utils import survival_table_from_events

Biblioteca importada, vamos usar agora as nossas variáveis T e C novamente para realizar o ajuste da tabela de vida.

lifetable = survival_table_from_events(T, C)

Tabela importada, vamos dar uma olhada no conjunto de dados.

print (lifetable)

      removed  observed  censored  entrance  at\_risk event\_at 0            2250      2247         3     10000    10000 1             676       531       145         0     7750 2             482       337       145         0     7074 3             185       129        56         0     6592 4             232        94       138         0     6407 5             299        85       214         0     6175 6             191        73       118         0     5876 7             127        76        51         0     5685 8             211        75       136         0     5558 9            2924        21      2903         0     5347 10            121        27        94         0     2423 11             46        27        19         0     2302 12             78        26        52         0     2256 13            111        16        95         0     2178 14             55        35        20         0     2067 15            107        29        78         0     2012 16            286        30       256         0     1905 17            156        23       133         0     1619 18            108        18        90         0     1463 19             49        11        38         0     1355 20             50        17        33         0     1306 21             61        13        48         0     1256 22            236        23       213         0     1195 23             99         6        93         0      959 24            168         9       159         0      860 25            171         7       164         0      692 26             58         6        52         0      521 27             77         2        75         0      463 28             29         6        23         0      386 29            105         1       104         0      357 30             69         0        69         0      252 31            183         0       183         0      183

Diferentemente do R que possuí a tabela de vida com a porcentagem relativa à probabilidade de sobrevivência, nesse caso vamos ter que fazer um pequeno ajuste para obter a porcentagem de acordo com o atributo entrance e at_risk. O ajuste se dará da seguinte forma:

survivaltable = lifetable.at_risk/np.amax(lifetable.entrance)

Ajustes efetuados, vamos ver como está a nossa tabela de vida.

survivaltable

Out[19]:

event_at 0 1.0000 1 0.7750 2 0.7074 3 0.6592 4 0.6407 5 0.6175 6 0.5876 7 0.5685 8 0.5558 9 0.5347 10 0.2423 11 0.2302 12 0.2256 13 0.2178 14 0.2067 15 0.2012 16 0.1905 17 0.1619 18 0.1463 19 0.1355 20 0.1306 21 0.1256 22 0.1195 23 0.0959 24 0.0860 25 0.0692 26 0.0521 27 0.0463 28 0.0386 29 0.0357 30 0.0252 31 0.0183 Name: at_risk, dtype: float64

Vamos transformar a nossa tabela de vida em um objeto do pandas para melhor manipulação do conjunto de dados.

survtable = pd.DataFrame(survivaltable)

Para casos de atualização de Churn-at-Risk podemos definir uma função que já terá a tabela de vida e poderá fazer a atribuição da probabilidade de sobrevivência de acordo com os dias de sobrevivência. Para isso vamos fazer uma função simples usando o próprio python.

def survival_probability( int ): survtable[“at_risk”].iloc[int] print (“The probability of Survival after”, int, “days is”, survtable[“at_risk”].iloc[int]*100, “%”) return;

Nesse caso vamos ver a chance de sobrevivência usando o nosso modelo Kaplan-Meier já ajustado para uma assinatura que tenha 22 dias de vida.

In [22]:

survival_probability(22)

The probability of Survival after 22 days is 11.95 %

Ou seja, essa assinatura tem apenas 11.95% de probabilidade de estar ativa, o que significa que em algum momento muito próximo ela pode vir a ser cancelada.

Conclusão

Como podemos ver acima, usando análise de sobrevivência podemos tirar insights interessantes em relação ao nosso conjunto de dados, em especial para descobrirmos a duração das assinaturas em nossa base de dados, e estimar um tempo até o evento de churn.

Os dados utilizados refletem o comportamento de dois produtos reais, porém, que foram anonimizados por questões óbvias de NDA. Contudo nada impede a utilização e a adaptação desse código para outros experimentos. Um ponto importante em relação a essa base de dados é que como pode ser observado temos uma censura à direita muito acentuada o que limita um pouco a visão dos dados a longo prazo, principalmente se houver algum tipo de cauda longa no evento de churn.

Como coloquei no São Paulo Big Data Meetup de Março há uma série de arquiteturas que podem ser combinadas com esse tipo de análise, em especial métodos de Deep Learning que podem ser um endpoint de um pipeline de predição.

Espero que tenham gostado e quaisquer dúvidas mandem uma mensagem para flavioclesio at gmail.com

PS: Agradecimentos especiais aos meus colegas e revisores Eiti Kimura, Gabriel Franco e Fernanda Eleuterio.