RFM com K=4: Quatro Clusters, mas Só Dois Contam Alguma Coisa

O capítulo 4 treinou os modelos e comparou métrica. Agora eu abro a caixa-preta de verdade: SHAP pros cenários de atraso e nota, e a decisão final de quantos grupos de cliente existem no RFM. Notebook completo, executado, público no Colab. Esse notebook também parte do zero (recarrega os 9 CSVs, remonta o dataframe mestre) e reconstrói as features do capítulo 3:

geo_media = geolocation.groupby('geolocation_zip_code_prefix')[['geolocation_lat', 'geolocation_lng']].mean().reset_index()

def haversine(lat1, lon1, lat2, lon2):
    R = 6371
    phi1, phi2 = np.radians(lat1), np.radians(lat2)
    dphi = np.radians(lat2 - lat1)
    dlambda = np.radians(lon2 - lon1)
    a = np.sin(dphi / 2) ** 2 + np.cos(phi1) * np.cos(phi2) * np.sin(dlambda / 2) ** 2
    return 2 * R * np.arcsin(np.sqrt(a))

m = mestre.copy()
m['customer_zip_code_prefix'] = m['customer_zip_code_prefix'].astype('int64')
m['seller_zip_code_prefix'] = m['seller_zip_code_prefix'].astype('int64')
geo_media['geolocation_zip_code_prefix'] = geo_media['geolocation_zip_code_prefix'].astype('int64')
m = m.merge(geo_media.rename(columns={'geolocation_zip_code_prefix': 'customer_zip_code_prefix', 'geolocation_lat': 'cust_lat', 'geolocation_lng': 'cust_lng'}), on='customer_zip_code_prefix', how='left')
m = m.merge(geo_media.rename(columns={'geolocation_zip_code_prefix': 'seller_zip_code_prefix', 'geolocation_lat': 'sell_lat', 'geolocation_lng': 'sell_lng'}), on='seller_zip_code_prefix', how='left')
m['distance_km'] = haversine(m['cust_lat'], m['cust_lng'], m['sell_lat'], m['sell_lng'])
m['approval_hours'] = (m['order_approved_at'] - m['order_purchase_timestamp']).dt.total_seconds() / 3600
m['month'] = m['order_purchase_timestamp'].dt.month
m['atraso_dias'] = (m['order_delivered_customer_date'] - m['order_estimated_delivery_date']).dt.days
m['delivery_days'] = (m['order_delivered_customer_date'] - m['order_purchase_timestamp']).dt.days

Mesmo código de sempre: distance_km (Haversine, já com o fix do int64 do capítulo 3), approval_hours, month, atraso_dias e delivery_days, tudo já coberto com calma lá atrás.

SHAP: por que o modelo pensa que vai atrasar

feature_importances_, que eu já usei no capítulo 3, é uma média global calculada durante o treino, mede o quanto cada feature ajudou a separar as classes na árvore inteira. SHAP (SHapley Additive exPlanations) é mais rigoroso: pra cada previsão individual, ele calcula o quanto cada feature empurrou aquela previsão específica pra cima ou pra baixo, partindo de uma "previsão média" de referência. A soma dessas contribuições explica exatamente o valor final. Fazendo essa conta pra cada pedido do conjunto de teste e tirando a média do valor absoluto, chego numa importância mais confiável que a do capítulo 3.

Usei o XGBoost aqui (não o Random Forest), porque foi o único que realmente "aposta" no positivo no capítulo 4, faz mais sentido abrir a caixa-preta de um modelo que decide alguma coisa em vez de um que nunca cruza o threshold.

delay_df = (
    m.drop_duplicates(subset='order_id')
    [['order_id', 'distance_km', 'product_weight_g', 'price', 'freight_value', 'approval_hours', 'month', 'atraso_dias']]
    .dropna()
)
delay_df['atrasado'] = (delay_df['atraso_dias'] > 0).astype(int)

features_1 = ['distance_km', 'product_weight_g', 'price', 'freight_value', 'approval_hours', 'month']
X1, y1 = delay_df[features_1], delay_df['atrasado']
X1tr, X1te, y1tr, y1te = train_test_split(X1, y1, test_size=0.2, random_state=42, stratify=y1)

modelo_delay = xgb.XGBClassifier(n_estimators=100, max_depth=6, random_state=42, eval_metric='logloss')
modelo_delay.fit(X1tr, y1tr)

explainer_delay = shap.TreeExplainer(modelo_delay)
shap_values_delay = explainer_delay.shap_values(X1te)
importancia_delay = np.abs(shap_values_delay).mean(axis=0)

Retreino o XGBoost com os mesmos hiperparâmetros do capítulo 4 e chamo shap.TreeExplainer(modelo_delay), um explicador desenhado especificamente pra modelos de árvore, que calcula a contribuição exata de cada feature sem precisar de aproximação por amostragem. explainer_delay.shap_values(X1te) devolve uma matriz do tamanho do conjunto de teste inteiro (uma linha por pedido, uma coluna por feature), cada célula é o quanto aquela feature específica empurrou aquela previsão específica pra cima ou pra baixo. np.abs(...).mean(axis=0) tira o valor absoluto (não importa a direção, só o quanto empurrou) e faz a média por coluna, isso é o número que vira o gráfico abaixo.

Carregando dados reais...

Mesma ordem do capítulo 3: month disparado na frente, distance_km em segundo. Dois métodos diferentes (importância de árvore e SHAP) concordando na mesma ordem é um sinal bom, o achado não é artefato de um jeito específico de medir, é robusto ao método.

SHAP: o que pesa na nota

Mesma lógica pro XGBoost que prevê review_score:

s2 = (
    m.drop_duplicates(subset='order_id')
    [['order_id', 'atraso_dias', 'delivery_days', 'price', 'freight_value', 'payment_installments', 'review_score']]
    .dropna()
)
features_2 = ['atraso_dias', 'delivery_days', 'price', 'freight_value', 'payment_installments']
X2, y2 = s2[features_2], s2['review_score']
X2tr, X2te, y2tr, y2te = train_test_split(X2, y2, test_size=0.2, random_state=42)

modelo_review = xgb.XGBRegressor(n_estimators=100, max_depth=6, random_state=42)
modelo_review.fit(X2tr, y2tr)

explainer_review = shap.TreeExplainer(modelo_review)
shap_values_review = explainer_review.shap_values(X2te)
importancia_review = np.abs(shap_values_review).mean(axis=0)

Mesma receita da seção anterior (treina o XGBoost, TreeExplainer, shap_values, média do valor absoluto), só troca o que sai do outro lado porque aqui o alvo é contínuo (review_score), não uma classe: é um regressor, não um classificador.

Carregando dados reais...

De novo a mesma ordem da informação mútua do capítulo 3: atraso e tempo de entrega dominam, o resto é ruído comparado a esses dois. Repara numa coisa: o XGBoost apareceu nos dois cenários que eu escolhi pra abrir com SHAP, mas isso não quer dizer que ele é sempre o melhor modelo. No sweep de threshold que já deixei pronto pro capítulo 0, no ponto ótimo (0,14, não o 0,5 padrão que o capítulo 4 usou) o Random Forest ultrapassa o XGBoost em F1. Guarda essa, o capítulo 0 deixa você brincar com esse limiar ao vivo.

RFM: a segmentação final

O capítulo 4 varreu K de 2 a 8 por silhueta. K=2 vencia disparado (0,743), mas era só o corte trivial que o capítulo 3 já tinha achado sem clustering nenhum: comprou mais de uma vez ou não. Tirando o K=2, o próximo melhor por silhueta era K=4 (0,49, maior que 5, 6, 7 e 8), então é esse que uso pra segmentação de verdade, nos 95.560 clientes únicos do capítulo 3.

data_max = orders['order_purchase_timestamp'].max()
orders_validos = orders[orders['order_status'] != 'canceled']
oc = orders_validos.merge(customers, on='customer_id', how='left')
pagamento_por_pedido = payments.groupby('order_id')['payment_value'].sum().reset_index()
oc = oc.merge(pagamento_por_pedido, on='order_id', how='left')

rfm = oc.groupby('customer_unique_id').agg(
    recencia_dias=('order_purchase_timestamp', lambda x: (data_max - x.max()).days),
    frequencia=('order_id', 'nunique'),
    valor_monetario=('payment_value', 'sum'),
).reset_index().dropna()

scaler = StandardScaler()
rfm_escalado = scaler.fit_transform(rfm[['recencia_dias', 'frequencia', 'valor_monetario']])

kmeans_final = KMeans(n_clusters=4, random_state=42, n_init=10)
rfm['cluster'] = kmeans_final.fit_predict(rfm_escalado)

perfil = rfm.groupby('cluster').agg(
    recencia_media=('recencia_dias', 'mean'),
    frequencia_media=('frequencia', 'mean'),
    valor_medio=('valor_monetario', 'mean'),
    n_clientes=('customer_unique_id', 'count'),
).reset_index()

Reconstruo o RFM idêntico ao capítulo 4, mas agora com K fixo em 4 em vez do sweep de 2 a 8: kmeans_final.fit_predict(rfm_escalado) treina e já devolve o cluster de cada cliente numa passada só, vira uma coluna nova (rfm['cluster']). O .groupby('cluster').agg(...) seguinte é o que constrói a tabela de perfil abaixo, a média de cada variável RFM dentro de cada cluster mais a contagem de clientes.

Batizando cada um pelo que os números realmente mostram: cluster 0 é o cliente recente de entrada (comprou uma vez, faz pouco tempo, ticket baixo), cluster 2 é o cliente dormente (comprou uma vez, faz quase um ano e dois meses, mesmo ticket baixo do cluster 0), cluster 1 é o grande comprador esporádico (ticket quase 9 vezes maior que os dois anteriores, mas praticamente nunca volta), e cluster 3 é o cliente recorrente de verdade, o único com frequência média acima de 1 (2,114 pedidos), batendo com os 3,06% de recorrentes que o capítulo 3 já tinha achado sem clustering nenhum.

Aqui mora o achado real do capítulo: clusters 0 e 2 juntos são 94,3% da base inteira, e a única coisa que os separa é recência, o ticket médio e a frequência são praticamente idênticos entre eles (R$ 134,58 contra R$ 133,72, ambos com frequência 1,000). K=4 é estatisticamente o melhor não-trivial, mas na prática só entrega duas segmentações que importam pro negócio (o grande comprador esporádico e o recorrente), o resto é uma forma mais chique de dizer "comprou recentemente" versus "comprou há muito tempo". K estatisticamente ótimo não é sinônimo de segmentação rica, é preciso olhar o centróide de cada cluster e perguntar se a diferença importa de verdade pra alguma decisão de negócio.

Preparando a vitrine interativa

Três exportações a mais que não geram gráfico aqui no capítulo, mas alimentam os widgets interativos do capítulo 0: uma amostra do RFM com o cluster de cada cliente, o sweep de threshold nos dois modelos que realmente decidem alguma coisa, e o fluxo região → entrega → nota.

amostra = rfm.sample(n=min(500, len(rfm)), random_state=42)[['recencia_dias', 'frequencia', 'valor_monetario', 'cluster']]

Só uma amostra de 500 clientes, rodar o scatter 3D com quase 95 mil pontos travaria o navegador sem ganhar nada em clareza visual, guardando as três variáveis RFM originais mais o cluster de cada um.

modelos_sweep = {
    'Random Forest': RandomForestClassifier(n_estimators=100, max_depth=8, random_state=42),
    'XGBoost': xgb.XGBClassifier(n_estimators=100, max_depth=6, random_state=42, eval_metric='logloss'),
}

sweep = []
for nome, modelo in modelos_sweep.items():
    modelo.fit(X1tr, y1tr)
    proba = modelo.predict_proba(X1te)[:, 1]
    for limiar in np.arange(0.0, 1.01, 0.02):
        pred = (proba >= limiar).astype(int)
        sweep.append({
            'model': nome,
            'threshold': round(float(limiar), 2),
            'precision': round(float(precision_score(y1te, pred, zero_division=0)), 4),
            'recall': round(float(recall_score(y1te, pred, zero_division=0)), 4),
            'f1': round(float(f1_score(y1te, pred, zero_division=0)), 4),
        })

Retreino Random Forest e XGBoost, os dois que realmente apostam em alguma previsão, e varro np.arange(0.0, 1.01, 0.02), 51 valores de threshold igualmente espaçados. Pra cada um, (proba >= limiar).astype(int) recalcula a decisão binária naquele corte específico e recalcula as três métricas: é o dado bruto por trás do slider que você arrasta no capítulo 0, cada posição do slider é, literalmente, uma linha diferente dessa lista sweep.

regiao_por_uf = {
    'AC': 'Norte', 'AP': 'Norte', 'AM': 'Norte', 'PA': 'Norte', 'RO': 'Norte', 'RR': 'Norte', 'TO': 'Norte',
    'AL': 'Nordeste', 'BA': 'Nordeste', 'CE': 'Nordeste', 'MA': 'Nordeste', 'PB': 'Nordeste', 'PE': 'Nordeste', 'PI': 'Nordeste', 'RN': 'Nordeste', 'SE': 'Nordeste',
    'DF': 'Centro-Oeste', 'GO': 'Centro-Oeste', 'MT': 'Centro-Oeste', 'MS': 'Centro-Oeste',
    'ES': 'Sudeste', 'MG': 'Sudeste', 'RJ': 'Sudeste', 'SP': 'Sudeste',
    'PR': 'Sul', 'RS': 'Sul', 'SC': 'Sul',
}

sk = m.drop_duplicates(subset='order_id')[['order_id', 'customer_state', 'atraso_dias', 'review_score']].dropna()
sk['regiao'] = sk['customer_state'].map(regiao_por_uf)
sk['status_entrega'] = np.where(sk['atraso_dias'] > 0, 'Atrasado', 'No prazo')
sk['faixa_nota'] = pd.cut(sk['review_score'], bins=[0, 2, 3, 5], labels=['Nota baixa (1-2)', 'Nota média (3)', 'Nota alta (4-5)'])

fluxo1 = sk.groupby(['regiao', 'status_entrega']).size().reset_index(name='value')
fluxo2 = sk.groupby(['status_entrega', 'faixa_nota'], observed=True).size().reset_index(name='value')

E o fluxo do sankey em duas etapas: primeiro mapeio cada estado pra uma das cinco regiões do Brasil (regiao_por_uf), depois construo três colunas categóricas (região, se atrasou ou não, faixa de nota) e conto quantos pedidos caem em cada combinação região→status e status→faixa. Duas contagens separadas (fluxo1, fluxo2) que, juntas no gráfico, formam o funil região → entrega → nota inteiro que você vê no capítulo 0.

Fechando o projeto

Os quatro cenários, resumidos pelo que realmente entregam: atraso tem AUC honesta (0,74) mas nenhum modelo pronto pra decisão automática sem ajustar o threshold primeiro. Nota esbarra num teto baixo porque as únicas features realmente informativas (atraso, tempo de entrega) já foram esgotadas no capítulo 4, sobra pouco sinal pra qualquer modelo explorar. Frete é o cenário mais "resolvido" de todos, XGBoost com R² 0,64 já é utilizável pra estimativa de verdade. RFM tem quatro clusters no papel, mas só dois de fato diferenciam alguém (juntos, apenas 5,6% da base), o resto é ruído de recência disfarçado de segmentação.

Limitações que valem registrar: o dataset é uma foto parada de 2016 a 2018, sem atualização contínua. Não tem dado de churn de verdade, um cliente "dormente" pode ter simplesmente parado de comprar, ou migrado pra outro marketplace, o dataset não distingue os dois casos. SHAP só nos dois cenários mais ricos narrativamente (atraso e nota), escolha editorial, não limitação técnica. Nenhum modelo passou por busca de hiperparâmetro, os números de n_estimators/max_depth são escolhas razoáveis, não um ótimo.

Se quiser reviver tudo isso sem ler mais nenhuma linha de prosa, o capítulo 0 é a vitrine: o mesmo threshold que só apareceu em texto aqui vira um slider que você arrasta com a mão, o RFM vira um scatter 3D que gira, e mais um punhado de jeitos diferentes de olhar pros mesmos dados que atravessaram esses cinco capítulos. E se quiser mais um ângulo técnico, o capítulo 6 fecha o projeto com a única pergunta que os quatro cenários originais nunca tocaram: quantos pedidos vêm por aí.