MLOps

Um modelo num notebook não gera valor nenhum; o que gera valor são predições entregues de forma confiável a quem age com base nelas. Até aqui, toda aula terminava quando o escore de teste era impresso. Esta começa aí. MLOps (machine-learning operations) é o conjunto de hábitos de engenharia que levam um pipeline ajustado para produção e o mantêm confiável depois que ele chega lá: empacotamento, serving, versionamento, testes, monitoramento e retreino. O diagrama do fluxo de trabalho das primeiras aulas prometia que "implantar e monitorar" volta para os dados; esta aula é esse laço.

Cada seção abaixo torna uma ideia concreta com um pequeno experimento medido: o que um arquivo de modelo contém de verdade, o que um pipeline faz com uma requisição malformada, quanto uma semente aleatória muda e como uma estatística como o PSI detecta que o mundo se moveu. O lab do final leva um modelo real pelo laço inteiro, com drift real: o modelo de telemarketing de um banco é treinado com ligações de 2008–2009 e depois encontra as ligações de 2009–2010, quando a fração de clientes que disseram sim passou de cerca de 7% para mais de 50%.

O laço desenhado na primeira seção é o sistema; a aula o percorre em cinco etapas, desenhadas aqui: 1. empacotar o modelo, o que exatamente vai no arquivo; 2. servir predições, como as requisições chegam ao modelo e o que ele faz com as malformadas; 3. reproduzir e versionar, para que o modelo em produção possa ser reconstruído; 4. monitorar o drift, perceber que o mundo mudou; 5. mudar com segurança, o desalinhamento treino–serving, os laços de realimentação e a troca de um modelo sem quebrar nada. Cada etapa leva à sua seção; nos slides, uma faixa embaixo mostra em que etapa você está.

Slides

Clique nos slides e use as setas ← → para navegar; F abre em tela cheia e Esc mostra todos os slides. Abrir os slides numa aba nova

O modelo é a caixa pequena

Sculley et al. (2015), em Hidden Technical Debt in Machine Learning Systems, fizeram a observação que definiu a agenda desta área: num sistema de ML em produção, o código que aprende é uma caixinha cercada de infraestrutura (coleta de dados, validação, extração de atributos, serving, monitoramento, configuração); a maioria das falhas acontece na infraestrutura. O próprio lab desta aula ilustra a proporção: das suas 212 linhas de código (sem contar linhas em branco e comentários), 4 chamam fit e outras 5 definem o modelo. O resto carrega, verifica, empacota, serve, monitora e compara.

Cada seta é um lugar onde algo pode quebrar sem levantar erro nenhum:

Empacotando um modelo

O artefato é o pipeline inteiro

Etapa 1. A aula de pipelines terminou com uma linha, joblib.dump(pipe, ...). Aqui ela está em contexto. Os dados são a tabela de crédito alemã (1.000 empréstimos passados; o alvo diz se o empréstimo deu errado), com três colunas numéricas e cinco categóricas:

import os
import tempfile

import joblib
import numpy as np
import pandas as pd
from sklearn.compose import make_column_transformer
from sklearn.datasets import fetch_openml
from sklearn.linear_model import LogisticRegression
from sklearn.metrics import roc_auc_score
from sklearn.model_selection import train_test_split
from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import OneHotEncoder, StandardScaler

credit = fetch_openml(data_id=31, as_frame=True)          # crédito alemão: 1.000 empréstimos passados
num = ['duration', 'credit_amount', 'age']
cat = ['checking_status', 'savings_status', 'purpose', 'employment', 'housing']
X = credit.data[num + cat]
y = (credit.target == 'bad').astype(int)                  # 1 = o empréstimo deu errado
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.25, stratify=y,
                                                    random_state=0)

pipe = make_pipeline(
    make_column_transformer((StandardScaler(), num),
                            (OneHotEncoder(handle_unknown='ignore'), cat)),
    LogisticRegression(max_iter=1000))
pipe.fit(X_train, y_train)
proba = pipe.predict_proba(X_test)[:, 1]
print(f"test ROC AUC: {roc_auc_score(y_test, proba):.3f}")

folder = tempfile.mkdtemp()                               # uma pasta de rascunho para esta aula
path = os.path.join(folder, 'credit-pipeline.joblib')
joblib.dump(pipe, path)                                   # o pipeline INTEIRO, um arquivo
loaded = joblib.load(path)
print(f"{os.path.getsize(path):,} bytes on disk")
print("identical probabilities after reload:",
      np.array_equal(proba, loaded.predict_proba(X_test)[:, 1]))

test ROC AUC: 0.783
4,721 bytes on disk
identical probabilities after reload: True

A primeira verificação que qualquer empacotamento precisa passar é a ida e volta: salvar, carregar e comparar as predições bit a bit nas mesmas linhas. np.array_equal é o teste certo aqui, não "perto o bastante": nada em salvar e carregar deveria mudar um único dígito.

Por que o pipeline inteiro e não só o modelo? Tente servir o modelo sozinho:

model_only = pipe[-1]                     # a LogisticRegression, sem o pré-processamento
try:
    model_only.predict_proba(X_test)
except ValueError as e:
    print('ValueError:', e)
print('the model expects', model_only.n_features_in_, 'columns; the raw table has', X_test.shape[1])

ValueError: X has 8 features, but LogisticRegression is expecting 30 features as input.
the model expects 30 columns; the raw table has 8

A regressão logística nunca viu a tabela crua. Ela viu as 30 colunas produzidas pelo column transformer: 3 números padronizados e 27 indicadores one-hot. Um programa de serving que distribui só o modelo precisa reconstruir essas 30 colunas à mão, com as mesmas médias e desvios padrão, a mesma ordem de categorias e o mesmo tratamento de valores desconhecidos. Qualquer diferença produz predições erradas sem erro nenhum. Distribuir o pipeline torna o pré-processamento parte do artefato.

Um pickle é um programa

Arquivos do joblib são pickles do Python. Um pickle não é uma tabela de números; é uma sequência de instruções para reconstruir objetos, instruções que podem incluir "chame esta função com estes argumentos". Qualquer classe pode escolher a função por meio de __reduce__:

import pickle

class Innocent:
    def __reduce__(self):                 # diz ao pickle como reconstruir este objeto...
        return (print, ('this line ran while the file was being LOADED',))

blob = pickle.dumps(Innocent())
obj = pickle.loads(blob)                  # ...e o carregamento obedece: chama print(...)

this line ran while the file was being LOADED

Nada foi chamado explicitamente; pickle.loads executou print porque o arquivo mandou. Troque print por os.system e o arquivo roda um comando de shell na máquina que o carrega. A regra sai direto daí: carregue arquivos de modelo apenas de fontes em que você confiaria para rodar código na sua máquina, o que, para um registro de modelos, significa controlar quem pode escrever nele.

Fixe as versões

Um pickle também registra a versão da biblioteca que o escreveu. O scikit-learn carimba sua versão em todo estimador que salva e compara o carimbo ao carregar. Para ver o que acontece com um arquivo de outra versão, forjamos um editando um único dígito do carimbo:

import warnings

import sklearn
from sklearn.exceptions import InconsistentVersionWarning

raw = open(path, 'rb').read()
v = sklearn.__version__.encode()
print(f"the bytes {v!r} appear {raw.count(v)} time(s) inside the file")

# forja um arquivo escrito por outra versão: mesmos bytes, um dígito do carimbo trocado
fake = os.path.join(folder, 'old-pipeline.joblib')
open(fake, 'wb').write(raw.replace(v, b'0' + v[1:]))
with warnings.catch_warnings(record=True) as caught:
    warnings.simplefilter('always')
    joblib.load(fake)
stale = [w.message for w in caught if isinstance(w.message, InconsistentVersionWarning)]
print('InconsistentVersionWarning for:', list(dict.fromkeys(m.estimator_name for m in stale)))
print(str(stale[0])[:100] + '...')

the bytes b'1.9.0' appear 1 time(s) inside the file
InconsistentVersionWarning for: ['StandardScaler', 'OneHotEncoder', 'ColumnTransformer', 'LogisticRegression', 'Pipeline']
Trying to unpickle estimator StandardScaler from version 0.9.0 when using version 1.9.0. This might ...

Todo estimador dentro do pipeline reclama; ainda assim, o arquivo carrega. É esse o perigo: os atributos internos de um estimador podem mudar entre versões, então um modelo pode carregar e depois predizer diferente ou falhar mais tarde num caminho de código que o teste não percorreu. O aviso passa despercebido com facilidade no log de um servidor. Três hábitos resolvem:

• guardar as versões (Python, scikit-learn, NumPy) ao lado do modelo, num arquivo de metadados (o lab faz isso);
• fixá-las no ambiente de serving (scikit-learn==1.9.0 no requirements.txt ou uma imagem de contêiner);
• ao atualizar a biblioteca, retreinar e revalidar em vez de recarregar o arquivo antigo. Numa suíte de testes, warnings.simplefilter('error', InconsistentVersionWarning) transforma o aviso numa falha.

Um formato mais seguro: skops

A biblioteca skops salva objetos do scikit-learn num formato que não executa código arbitrário: ao carregar, ela só reconstrói tipos de uma lista permitida; qualquer outro tipo precisa ser declarado confiável explicitamente.

import skops.io as sio                    # pip install skops
from sklearn.preprocessing import FunctionTransformer

sk_path = os.path.join(folder, 'credit-pipeline.skops')
sio.dump(pipe, sk_path)
print('untrusted types in the plain pipeline:', sio.get_untrusted_types(file=sk_path))
same = np.array_equal(proba, sio.load(sk_path).predict_proba(X_test)[:, 1])
print('loaded without trusting anything, identical probabilities:', same)

def cap_amount(X):                        # uma função própria dentro do pipeline
    return np.minimum(X, 15000)

custom = make_pipeline(FunctionTransformer(cap_amount), LogisticRegression())
custom.fit(X_train[['credit_amount']], y_train)
custom_path = os.path.join(folder, 'custom.skops')
sio.dump(custom, custom_path)
print('untrusted types in the custom pipeline:', sio.get_untrusted_types(file=custom_path))
try:
    sio.load(custom_path)
except Exception as e:
    print('loading it without trust raises', type(e).__name__)

untrusted types in the plain pipeline: []
loaded without trusting anything, identical probabilities: True
untrusted types in the custom pipeline: ['__main__.cap_amount']
loading it without trust raises UntrustedTypesFoundException

O pipeline simples só contém tipos do scikit-learn e do NumPy, então carrega sem perguntas e prediz de forma idêntica. A função própria cap_amount é código que o skops não conhece, então o carregamento é interrompido até que alguém leia a lista e passe trusted=[...]. Para servir fora do Python, o formato ONNX (via skl2onnx) converte um pipeline ajustado num grafo portátil que runtimes em C++, Java ou JavaScript conseguem executar.

Servindo predições

Batch, online, streaming, edge

Etapa 2. A forma como o modelo é chamado define todo o resto da arquitetura:

Prefira batch sempre que a decisão puder esperar: é mais barato (medido abaixo), mais fácil de monitorar e mais fácil de rodar de novo. O resto desta seção trata do caso online, em que a requisição de um desconhecido chega ao modelo.

O pipeline aceita lixo

Um pipeline verifica muito pouco. Cinco requisições, uma boa e quatro quebradas de jeitos que acontecem em sistemas reais:

good = {'duration': 24, 'credit_amount': 4000, 'age': 35, 'checking_status': '<0',
        'savings_status': '<100', 'purpose': 'new car', 'employment': '1<=X<4',
        'housing': 'own'}
requests = {
    'good request':        good,
    'age = -3':            {**good, 'age': -3},
    'amount in thousands': {**good, 'credit_amount': 4},
    'unseen purpose':      {**good, 'purpose': 'crypto'},
    'housing missing':     {k: v for k, v in good.items() if k != 'housing'},
}
for name, payload in requests.items():
    try:
        p = loaded.predict_proba(pd.DataFrame([payload]))[0, 1]
        print(f"{name:20s} P(bad) = {p:.3f}")
    except Exception as e:
        print(f"{name:20s} {type(e).__name__}: {str(e)[:45]}")

good request         P(bad) = 0.626
age = -3             P(bad) = 0.784
amount in thousands  P(bad) = 0.574
unseen purpose       P(bad) = 0.519
housing missing      ValueError: columns are missing: {'housing'}

Só a coluna ausente levanta erro. As outras três requisições quebradas recebem uma probabilidade que parece tão confiável quanto a da boa. Uma idade negativa é padronizada como qualquer número. Um empréstimo de 4 marcos alemães é só um empréstimo muito pequeno. O propósito desconhecido 'crypto' é codificado como tudo zero por handle_unknown='ignore', o que quer dizer "nenhum dos propósitos que conheço". Com o padrão handle_unknown='error', ele teria levantado erro: uma categoria nova surgida na origem dos dados derrubaria o serviço. Nenhuma das duas escolhas avisa ninguém de que algo está errado.

Um schema e uma função que o verifica

A correção é registrar num schema o que o modelo viu no treino e conferir cada requisição com ele antes de predizer:

def make_schema(X):
    """O que o modelo viu no treino: uma faixa por número, uma lista por categoria."""
    schema = {}
    for col in X.columns:
        if pd.api.types.is_numeric_dtype(X[col]):
            schema[col] = {'type': 'number', 'min': float(X[col].min()), 'max': float(X[col].max())}
        else:
            schema[col] = {'type': 'category', 'values': sorted(X[col].astype(str).unique())}
    return schema

def validate(payload, schema):
    """Devolve uma lista de problemas; lista vazia quer dizer que o payload é aceitável."""
    errors = [f"unexpected field {k!r}" for k in payload if k not in schema]
    for col, rule in schema.items():
        v = payload.get(col)
        if v is None:
            errors.append(f"{col}: missing")
        elif rule['type'] == 'number':
            if isinstance(v, bool) or not isinstance(v, (int, float)):
                errors.append(f"{col}: expected a number, got {v!r}")
            elif not rule['min'] <= v <= rule['max']:
                errors.append(f"{col}: {v} outside training range [{rule['min']:g}, {rule['max']:g}]")
        elif v not in rule['values']:
            errors.append(f"{col}: unknown category {v!r}")
    return errors

schema = make_schema(X_train)
print(schema['age'])
print(schema['housing'])

{'type': 'number', 'min': 19.0, 'max': 75.0}
{'type': 'category', 'values': ['for free', 'own', 'rent']}

A função de serving fica então com três linhas de lógica: validar, predizer e responder com a versão do modelo, para que toda resposta registrada possa ser rastreada até o artefato que a produziu.

def predict_one(payload, model=loaded, schema=schema, version='credit-2026.09.1'):
    """O contrato de serving inteiro: validar, predizer, responder com a versão do modelo."""
    errors = validate(payload, schema)
    if errors:
        return {'status': 'rejected', 'errors': errors, 'model_version': version}
    p = float(model.predict_proba(pd.DataFrame([payload]))[0, 1])
    return {'status': 'ok', 'p_bad': round(p, 3), 'model_version': version}

for name, payload in requests.items():
    answer = predict_one(payload)
    print(f"{name:20s} {answer['status']:8s} {answer.get('p_bad', '')} {answer.get('errors', '')}")

good request         ok       0.626
age = -3             rejected  ['age: -3 outside training range [19, 75]']
amount in thousands  rejected  ['credit_amount: 4 outside training range [250, 15945]']
unseen purpose       rejected  ["purpose: unknown category 'crypto'"]
housing missing      rejected  ['housing: missing']

A requisição boa recebe o mesmo 0,626 de antes; cada uma das quebradas é rejeitada com uma mensagem clara o bastante para que alguém possa agir. Se "fora da faixa de treino" deve rejeitar ou apenas avisar é uma decisão de negócio. Um empréstimo de 4 marcos alemães é certamente um bug de unidade, mas um cliente mais velho do que qualquer um nos dados de treino é uma pessoa real que o modelo simplesmente nunca viu. O lab usa dois níveis: rejeitar valores impossíveis e avisar (predizer, mas registrar o aviso) para valores que são apenas novos. Bibliotecas como pydantic (embutida no FastAPI), pandera e Great Expectations expressam as mesmas verificações de forma declarativa.

A mesma função atrás de um endpoint HTTP

O bloco abaixo não é executado nesta página: FastAPI é um framework web; um servidor não tem lugar num notebook. Salve make_schema, validate e predict_one num arquivo service.py (carregando o modelo com joblib.load); o serviço web inteiro é este:

# app.py  (ilustrativo: rode num terminal com uvicorn, não num notebook)
from fastapi import FastAPI

from service import predict_one           # a função testada nesta seção

app = FastAPI()

@app.post('/predict')
def predict(payload: dict):
    return predict_one(payload)

@app.get('/health')                        # para o balanceador de carga: "o serviço está vivo?"
def health():
    return {'status': 'ok'}

$ pip install fastapi uvicorn
$ uvicorn app:app --port 8000
$ curl -X POST localhost:8000/predict -H 'Content-Type: application/json' \
       -d '{"duration": 24, "credit_amount": 4000, "age": 35, "checking_status": "<0",
            "savings_status": "<100", "purpose": "new car", "employment": "1<=X<4",
            "housing": "own"}'
{"status":"ok","p_bad":0.626,"model_version":"credit-2026.09.1"}

Repare onde a lógica mora. A camada web tem um punhado de linhas sem nenhuma decisão; tudo o que pode estar errado está em predict_one, que é uma função comum que você testa num notebook, como acima. Deixar isso pronto para produção acrescenta um contêiner (Docker), para que o ambiente fixado seja o que de fato roda; um log de toda requisição e resposta (a matéria-prima do monitoramento); e uma versão na URL (/v1/predict), para que os clientes não quebrem com o próximo modelo.

Latência: paga-se por chamada, não por linha

Quanto tempo leva uma predição? Depende muito mais do número de chamadas do que do número de linhas:

import time

def per_row_ms(model, rows, repeats=5):
    """Melhor de `repeats` medições, em milissegundos por linha."""
    best = np.inf
    for _ in range(repeats):
        start = time.perf_counter()
        model.predict_proba(rows)
        best = min(best, time.perf_counter() - start)
    return 1000 * best / len(rows)

big = X_test.sample(10_000, replace=True, random_state=0)
one_at_a_time = np.mean([per_row_ms(loaded, big.iloc[[i]]) for i in range(200)])
in_one_batch = per_row_ms(loaded, big)
print(f"one row per call: {one_at_a_time:.3f} ms/row · 10,000 rows per call: {in_one_batch:.4f} ms/row"
      f" · ratio {one_at_a_time / in_one_batch:.0f}x")

Tempos dependem da máquina, então a saída deste bloco não é fixa. Na máquina que gerou esta página, uma execução imprimiu:

one row per call: 1.724 ms/row · 10,000 rows per call: 0.0006 ms/row · ratio 2688x

(script da figura)

Uma chamada custa cerca de 1,7 ms quer carregue uma linha, quer mil: esse tempo vai para trabalho fixo (montar o DataFrame, conferir nomes de colunas, despachar por três estimadores), não para aritmética. Pontuar 10.000 linhas numa só chamada custa uns 6 ms no total, entre 2.500 e 3.000 vezes menos por linha do que uma chamada por linha nesta máquina (a figura, que usa medianas de várias chamadas, fica perto do topo dessa faixa). Uma requisição HTTP soma seu próprio custo fixo (rede, parsing de JSON). Duas consequências: pontuar em batch sai quase de graça; um serviço online deve manter o pipeline leve, já que o custo fixo é pago a cada requisição.

Reprodutibilidade e versionamento

Etapa 3. "Qual modelo está rodando? Conseguiríamos reconstruí-lo exatamente?" Essas perguntas só têm resposta se cada ingrediente do modelo tiver uma versão:

Sementes

Muitos estimadores usam aleatoriedade: uma random forest sorteia amostras bootstrap e subconjuntos aleatórios de atributos. Duas florestas que diferem só na semente:

from sklearn.ensemble import RandomForestClassifier

def forest(seed):
    return make_pipeline(
        make_column_transformer((OneHotEncoder(handle_unknown='ignore'), cat),
                                remainder='passthrough'),
        RandomForestClassifier(n_estimators=200, random_state=seed)).fit(X_train, y_train)

p0, p1 = [forest(s).predict_proba(X_test)[:, 1] for s in (0, 1)]
p0_again = forest(0).predict_proba(X_test)[:, 1]
print(f"seed 0: AUC {roc_auc_score(y_test, p0):.3f} · seed 1: AUC {roc_auc_score(y_test, p1):.3f}")
print(f"seed 0 vs seed 1: {np.mean(p0 != p1):.0%} of probabilities differ, "
      f"{np.sum((p0 > 0.5) != (p1 > 0.5))} of {len(p0)} decisions flip")
print("seed 0 vs seed 0 again: identical =", np.array_equal(p0, p0_again))

seed 0: AUC 0.763 · seed 1: AUC 0.766
seed 0 vs seed 1: 96% of probabilities differ, 18 of 250 decisions flip
seed 0 vs seed 0 again: identical = True

Os dois AUCs são quase iguais, então no placar as florestas parecem o mesmo modelo. Para os solicitantes, não são: 18 dos 250 solicitantes de teste (7%) recebem uma decisão diferente dependendo da semente. Sem random_state, cada retreino sorteia uma semente nova. Um cliente que pergunta "por que fui recusado ontem e aprovado hoje?" merece uma resposta melhor do que "rodamos o script de novo". Com a semente fixa, a nova execução é idêntica. Não vá no sentido contrário, testando sementes até o escore de teste ficar o melhor possível: isso é escolher no conjunto de teste, sujeito à maldição do vencedor.

A impressão digital dos dados

Um modelo é função dos seus dados de treino, então os dados também precisam de uma identidade. Um hash de cada linha e valor fornece uma:

import hashlib

def data_hash(df):
    """Uma impressão digital das linhas, valores e ordem exatos usados no treino."""
    row_hashes = pd.util.hash_pandas_object(df, index=True).to_numpy()
    return hashlib.sha256(row_hashes.tobytes()).hexdigest()[:12]

changed = X_train.copy()
changed.iloc[0, changed.columns.get_loc('age')] += 1          # uma pessoa um ano mais velha
print('training data    ', data_hash(X_train))
print('one value changed', data_hash(changed))
print('rows shuffled    ', data_hash(X_train.sample(frac=1, random_state=0)))
print('same data again  ', data_hash(X_train.copy()))

training data     68ad024c9661
one value changed 7eb64a1901f0
rows shuffled     a662547f0332
same data again   68ad024c9661

Um único valor alterado entre 6.000 dá uma impressão digital completamente diferente; o mesmo acontece com os mesmos dados em outra ordem (a ordem importa: os folds da validação cruzada e alguns solvers dependem dela). Os mesmos dados dão sempre a mesma impressão digital. Guarde-a com o modelo; quando alguém perguntar se um modelo foi treinado com a extração deste mês ou do mês passado, a resposta é uma comparação de strings. Um hash apenas identifica os dados. Para reconstruir o modelo você também precisa das próprias linhas, que é o que ferramentas de versionamento de dados (DVC, lakeFS) ou simples snapshots datados e somente leitura oferecem.

Um registro de experimentos

Junte os ingredientes e cada execução de treino vira um registro. O rastreador honesto mais simples é um arquivo em que cada execução acrescenta uma linha e nada nunca é editado:

import json
import platform

from sklearn.model_selection import StratifiedKFold, cross_val_score

log_file = os.path.join(folder, 'runs.jsonl')

def log_run(params, metrics):
    """Acrescenta um experimento a um arquivo JSON-lines: uma linha por execução, nunca editada."""
    record = {'time': time.strftime('%Y-%m-%d %H:%M:%S'),
              'data': data_hash(X_train), 'python': platform.python_version(),
              'sklearn': sklearn.__version__, 'params': params, 'metrics': metrics}
    with open(log_file, 'a') as f:
        f.write(json.dumps(record) + '\n')

cv = StratifiedKFold(n_splits=5, shuffle=True, random_state=0)
for C in [0.01, 0.1, 1.0, 10.0]:
    candidate = pipe.set_params(logisticregression__C=C)
    auc = cross_val_score(candidate, X_train, y_train, cv=cv, scoring='roc_auc')
    log_run({'model': 'logistic', 'C': C}, {'cv_auc': round(auc.mean(), 4), 'cv_std': round(auc.std(), 4)})

runs = pd.json_normalize([json.loads(line) for line in open(log_file)])
print(runs[['data', 'params.C', 'metrics.cv_auc', 'metrics.cv_std']].to_string(index=False))

        data  params.C  metrics.cv_auc  metrics.cv_std
68ad024c9661      0.01          0.7495          0.0292
68ad024c9661      0.10          0.7599          0.0224
68ad024c9661      1.00          0.7522          0.0244
68ad024c9661     10.00          0.7442          0.0230

O registro responde a perguntas que a memória não responde: que configurações foram tentadas, com quais dados, com qual biblioteca. C = 0.1 tem a melhor média (0,7599), mas sua vantagem sobre C = 1 (0,0077) fica muito abaixo dos desvios padrão entre folds (cerca de 0,02): sem uma comparação fold a fold, nada diz que ela é real, então trate as configurações como empatadas. Guarde também as execuções perdedoras: elas são o registro de quantas coisas foram tentadas, que é justamente aquilo de que depende a maldição do vencedor.

O MLflow faz o mesmo, com interface web, armazenamento de artefatos e um registro de modelos que promove versões de "staging" para "production". O bloco abaixo é ilustrativo e não é executado aqui (o MLflow é uma instalação grande):

import mlflow                                   # pip install mlflow

with mlflow.start_run():
    mlflow.log_params({'model': 'logistic', 'C': 0.1, 'data': data_hash(X_train)})
    mlflow.log_metric('cv_auc', 0.7599)
    mlflow.sklearn.log_model(pipe, name='model', registered_model_name='credit')

Monitoramento e drift

Etapa 4. Software costuma quebrar fazendo barulho; modelos se degradam em silêncio. O serviço continua respondendo 200 OK enquanto suas predições pioram. O monitoramento tem três camadas:

1. Sistema: latência, vazão, erros, as preocupações usuais de qualquer serviço web;
2. Qualidade dos dados: violações do schema, picos de valores ausentes, valores fora da faixa de treino. Muitos incidentes atribuídos "ao modelo" são incidentes de dados em etapas anteriores ao modelo;
3. Drift estatístico: o mundo em que o modelo foi treinado já não é o mundo em que ele é usado.

Três maneiras de o mundo mudar

Um modelo é treinado com uma amostra de uma distribuição conjunta \(P(x, y)\), que pode mudar de três maneiras básicas. Elas são definidas pelo que fica fixo (Moreno-Torres et al., 2012):

Antes de medir o que cada uma faz com um modelo, precisamos de um jeito de medir uma mudança numa distribuição. Há dois padrões.

Medindo uma mudança: PSI

O índice de estabilidade populacional (population stability index) compara a fração de linhas que cai em cada faixa, com as faixas cortadas nos dados de referência (treino):

em que \(e_i\) é a fração esperada (referência) e \(a_i\) a fração observada (dados atuais) na faixa \(i\). Um exemplo resolvido com 4 faixas, cortadas nos quartis dos dados de treino, de modo que cada uma contém 25% das linhas de referência:

Todo termo é positivo, porque \(a_i - e_i\) e \(\ln(a_i/e_i)\) sempre têm o mesmo sinal. O PSI só é zero quando cada faixa mantém sua fração. (É a divergência de Kullback–Leibler tomada nos dois sentidos, \(\text{KL}(a\|e) + \text{KL}(e\|a)\).) O NumPy confirma a tabela:

e = np.array([0.25, 0.25, 0.25, 0.25])          # referência: os quartis dos dados de treino
a = np.array([0.10, 0.20, 0.30, 0.40])          # a fração das linhas desta semana em cada quartil
terms = (a - e) * np.log(a / e)
print(terms.round(4), '-> PSI =', terms.sum().round(3))

[0.1374 0.0112 0.0091 0.0705] -> PSI = 0.228

A leitura usual, uma convenção vinda do escore de crédito (Siddiqi, 2006), é: abaixo de 0,1, nenhuma mudança relevante; de 0,1 a 0,25, uma mudança moderada que merece atenção; acima de 0,25, uma mudança significativa. Para dar sentido a esses números, eis o PSI de um atributo normal cuja média se desloca por uma quantidade conhecida, com 10 faixas nos decis da referência:

def psi(ref, cur, bins=10, eps=1e-4):
    """Índice de estabilidade populacional de `cur` contra `ref`, nas faixas de quantis da referência."""
    edges = np.unique(np.quantile(ref, np.linspace(0, 1, bins + 1)))
    edges[0], edges[-1] = -np.inf, np.inf                    # faixas externas abertas
    e = np.histogram(ref, edges)[0] / len(ref)
    a = np.histogram(cur, edges)[0] / len(cur)
    e, a = np.clip(e, eps, None), np.clip(a, eps, None)      # uma faixa vazia daria log(0)
    return float(np.sum((a - e) * np.log(a / e)))

rng = np.random.default_rng(0)
ref = rng.normal(0, 1, 5000)
for shift in [0.0, 0.1, 0.25, 0.5, 1.0]:
    print(f"mean moved by {shift:4.2f} sd -> PSI {psi(ref, rng.normal(shift, 1, 5000)):.3f}")

mean moved by 0.00 sd -> PSI 0.005
mean moved by 0.10 sd -> PSI 0.015
mean moved by 0.25 sd -> PSI 0.064
mean moved by 0.50 sd -> PSI 0.250
mean moved by 1.00 sd -> PSI 0.929

O nível de alarme 0,25 corresponde a um deslocamento da média de meio desvio padrão; um desvio padrão inteiro dá cerca de 0,93. Sem deslocamento nenhum, o PSI não é exatamente zero (0,005): duas amostras de 5.000 linhas nunca têm histogramas idênticos. O corte eps importa na prática: uma faixa vazia nos dados atuais daria \(\ln 0\).

Desloque você mesmo a distribuição de produção. Encontre o deslocamento em que o PSI cruza 0,25; depois volte a média para 0, mude só a dispersão e, por fim, ligue o concept drift:

Com as frações exatas, meio desvio padrão dá PSI 0,24; troque para amostras de 5.000 linhas (o cenário do código acima) e o ruído de amostragem o move alguns centésimos para cima ou para baixo; é por isso que a tabela mostra 0,250. Uma dispersão maior ou menor com a mesma média também move o PSI, pelas faixas das pontas. E o concept drift, apresentado na tabela acima e medido mais abaixo, muda \(P(y \mid x)\) sem mover uma única linha entre as faixas: o PSI fica em 0 enquanto a acurácia cai.

Medindo uma mudança: o teste KS

A estatística de Kolmogorov–Smirnov compara as duas funções de distribuição acumulada empíricas e pega o maior salto vertical entre elas:

À mão, com cinco valores de cada lado:

from scipy.stats import ks_2samp

ref_s, cur_s = np.array([1, 2, 3, 4, 5]), np.array([3, 4, 5, 6, 7])
grid = np.union1d(ref_s, cur_s)
F_ref = np.searchsorted(np.sort(ref_s), grid, side='right') / len(ref_s)
F_cur = np.searchsorted(np.sort(cur_s), grid, side='right') / len(cur_s)
print(pd.DataFrame({'x': grid, 'F_ref': F_ref, 'F_cur': F_cur, 'gap': np.abs(F_ref - F_cur)})
      .T.to_string(header=False))
res = ks_2samp(ref_s, cur_s)
print(f"by hand D = {np.abs(F_ref - F_cur).max():.1f} · scipy D = {res.statistic:.1f}, p-value = {res.pvalue:.3f}")

x      1.0  2.0  3.0  4.0  5.0  6.0  7.0
F_ref  0.2  0.4  0.6  0.8  1.0  1.0  1.0
F_cur  0.0  0.0  0.2  0.4  0.6  0.8  1.0
gap    0.2  0.4  0.4  0.4  0.4  0.2  0.0
by hand D = 0.4 · scipy D = 0.4, p-value = 0.873

Em \(x = 2\), 40% da referência está nele ou abaixo dele; da amostra atual, nada: \(D = 0{,}4\). O teste pergunta então quão surpreendente seria um salto de 0,4 se as duas amostras viessem da mesma distribuição. Com cinco pontos de cada lado, nada surpreendente: \(p = 0{,}873\). A figura mostra as duas estatísticas num exemplo maior, 5.000 linhas cuja média se deslocou meio desvio padrão (script da figura):

As duas estatísticas respondem a perguntas diferentes, o que importa com dados grandes. Um deslocamento de 0,02 desvio padrão é irrelevante para qualquer modelo; veja o que cada estatística diz sobre ele conforme o número de linhas cresce:

rng = np.random.default_rng(1)
for n in [1_000, 10_000, 100_000, 1_000_000]:
    ref_n, cur_n = rng.normal(0, 1, n), rng.normal(0.02, 1, n)    # um deslocamento de 0,02 dp
    print(f"n = {n:>9,}: KS p-value {ks_2samp(ref_n, cur_n).pvalue:.2g} · PSI {psi(ref_n, cur_n):.4f}")

n =     1,000: KS p-value 0.097 · PSI 0.0344
n =    10,000: KS p-value 0.14 · PSI 0.0018
n =   100,000: KS p-value 0.00049 · PSI 0.0004
n = 1,000,000: KS p-value 1.6e-41 · PSI 0.0005

O p-valor do KS responde "existe alguma diferença?"; com linhas suficientes, a resposta é sempre sim: com um milhão de linhas, o mesmo deslocamento minúsculo dá \(p = 1{,}6 \times 10^{-41}\). O PSI mede quão grande é a diferença e fica perto de zero. Na outra ponta, com 1.000 linhas, o PSI é mais ruidoso (0,034, quase tudo ruído de amostragem). Daí saem duas regras práticas: com lotes grandes, dispare alertas por magnitude (PSI ou a estatística \(D\) do KS) e não por p-valor; e quando muitos atributos são testados de uma vez, espere alarmes falsos (50 atributos testados ao nível de 5% dão uns 2,5 alarmes por lote sem drift nenhum), então corrija o nível (Bonferroni: divida-o pelo número de testes), como fazem Rabanser et al. (2019) em seu benchmark de detectores de mudança.

As três mudanças, medidas

Agora os três tipos de mudança, simulados para sabermos exatamente qual aconteceu. Um atributo \(x\); a probabilidade da classe positiva sobe com \(x\), atinge o pico em \(x = 1{,}5\) e cai de novo. O modelo é um pequeno classificador de gradient boosting treinado com \(x \sim N(0, 1)\). Cada cenário sorteia 5.000 linhas novas:

from sklearn.ensemble import HistGradientBoostingClassifier

def true_p(x, move=0.0):
    """P(y = 1 | x): sobe com x, atinge o pico em x = 1.5 + move, depois cai."""
    z = x - move
    return 1 / (1 + np.exp(-(3 * z - z ** 2)))

def draw(n, mean, rng, move=0.0):
    x = rng.normal(mean, 1, n)
    return x, (rng.random(n) < true_p(x, move)).astype(int)

rng = np.random.default_rng(2)
x_ref, y_ref = draw(5000, 0, rng)
model = HistGradientBoostingClassifier(max_depth=3, random_state=0).fit(x_ref.reshape(-1, 1), y_ref)

pool_x, pool_y = draw(100_000, 0, rng)                           # para o prior shift:
pos, neg = pool_x[pool_y == 1][:1000], pool_x[pool_y == 0][:4000]  # mesma P(x | y), 20% positivos
scenarios = {'no change':       draw(5000, 0, rng),
             'covariate shift': draw(5000, 2.0, rng),            # P(x) muda, P(y | x) fixa
             'concept drift':   draw(5000, 0, rng, move=1.0),    # P(x) fixa, P(y | x) muda
             'prior shift':     (np.r_[pos, neg], np.r_[np.ones(1000), np.zeros(4000)].astype(int))}
rows = {}
for name, (x, y) in scenarios.items():
    p = model.predict_proba(x.reshape(-1, 1))[:, 1]
    rows[name] = {'PSI(x)': psi(x_ref, x), 'KS p': ks_2samp(x_ref, x).pvalue,
                  'accuracy': np.mean((p > 0.5) == y), 'precision': y[p > 0.5].mean(),
                  'actual rate': y.mean(), 'mean predicted': p.mean()}
print(pd.DataFrame(rows).T.round(3).to_string())

                 PSI(x)   KS p  accuracy  precision  actual rate  mean predicted
no change         0.004  0.695     0.816      0.769        0.462           0.465
covariate shift   3.505  0.000     0.728      0.729        0.720           0.780
concept drift     0.003  0.340     0.660      0.329        0.167           0.459
prior shift       0.152  0.000     0.794      0.491        0.200           0.339

(script da figura)

Leia a tabela uma linha de cada vez:

• Covariate shift: a população se moveu para \(x \approx 2\), PSI 3,5, p-valor do KS zero; todo monitor de entrada dispara. A acurácia cai porque o modelo agora é usado onde tinha poucos dados: depois do pico, onde a probabilidade verdadeira começa a cair, as árvores prolongam os últimos valores que aprenderam (a curva azul tracejada e plana à direita do segundo painel).
• Concept drift: as entradas não mudaram, PSI 0,003, p-valor do KS 0,34; nenhum monitor de entrada consegue enxergá-lo. O mesmo \(x\) agora merece outra resposta; a acurácia do modelo caiu. Só os rótulos revelam concept drift.
• Prior shift: os positivos foram de 46% para 20%. As entradas se mexem um pouco (PSI 0,15), a acurácia quase não muda, mas a precisão cai e as probabilidades agora estão altas demais: o modelo prevê 34% de positivos em média onde 20% são positivos. É o efeito da prevalência sobre a precisão, agora acontecendo em produção.

Só o concept drift é invisível sem rótulos, mas cada mudança quebra um número diferente, então um monitor precisa dos dois tipos de verificação: entradas (disponíveis na hora) e saídas comparadas com os desfechos (disponíveis depois).

Os rótulos chegam atrasados

A verdade costuma chegar com atraso: se um cliente cancela, sabe-se depois de 60 dias; se um empréstimo dá calote, depois de meses ou anos. O monitoramento funciona, portanto, em duas velocidades:

• imediatamente, sobre o que o serviço vê: verificações de qualidade dos dados, PSI ou KS de cada entrada em relação à referência de treino e a distribuição dos próprios escores do modelo (a fração sinalizada, a probabilidade média). Rabanser et al. (2019) mostraram que testar a distribuição das saídas do modelo está entre os detectores de mudança mais eficazes e custa uma linha;
• depois, quando os desfechos chegam: as métricas de verdade (AUC, precisão, calibração), calculadas em cada lote e comparadas com a estimativa de validação.

Quando os rótulos chegam muito tarde ou só para os casos que o modelo aprovou (um empréstimo recusado nunca ganha histórico de pagamento), uma pequena amostra auditada de casos rotulados é a única proteção contra concept drift.

Degradar, detectar, retreinar

Um drift gradual é o caso comum em produção. A mesma simulação, agora como um ano de lotes mensais: a população se move 0,3 desvio padrão por mês. O modelo estático nunca é tocado. O modelo monitorado é retreinado com os dados rotulados do mês sempre que o PSI da sua entrada, medido contra os dados com que ele foi treinado, passa de 0,25:

rng = np.random.default_rng(3)
new_model = lambda: HistGradientBoostingClassifier(max_depth=3, random_state=0)
static = new_model().fit(x_ref.reshape(-1, 1), y_ref)
live, live_ref = static, x_ref                 # o modelo servido e sua amostra de referência
print('month  mean   PSI  alert  static acc  live acc  best possible')
for month in range(13):
    x, y = draw(1000, 0.3 * month, rng)        # a população se move 0,3 dp por mês
    X_m = x.reshape(-1, 1)
    drift = psi(live_ref, x)
    alert = drift > 0.25
    print(f"{month:5d} {0.3 * month:5.1f} {drift:5.2f}  {'ALERT' if alert else '     '} "
          f"{np.mean(static.predict(X_m) == y):10.3f} {np.mean(live.predict(X_m) == y):9.3f}"
          f" {np.mean((true_p(x) > 0.5) == y):14.3f}")
    if alert:                                  # chegam os rótulos do mês: retreinar com eles
        live, live_ref = new_model().fit(X_m, y), x

month  mean   PSI  alert  static acc  live acc  best possible
    0   0.0  0.01             0.801     0.801          0.802
    1   0.3  0.09             0.809     0.809          0.815
    2   0.6  0.33  ALERT      0.780     0.780          0.785
    3   0.9  0.10             0.811     0.813          0.824
    4   1.2  0.33  ALERT      0.812     0.807          0.819
    5   1.5  0.14             0.783     0.805          0.814
    6   1.8  0.41  ALERT      0.759     0.814          0.830
    7   2.1  0.09             0.729     0.820          0.821
    8   2.4  0.30  ALERT      0.654     0.822          0.816
    9   2.7  0.18             0.507     0.807          0.805
   10   3.0  0.46  ALERT      0.417     0.803          0.802
   11   3.3  0.08             0.374     0.844          0.842
   12   3.6  0.34  ALERT      0.273     0.849          0.850

(script da figura)

O modelo estático se degrada sem uma única mensagem de erro: acurácia de 0,801 no mês 0, 0,654 no mês 8, 0,273 no mês 12. Isso é muito pior do que jogar uma moeda, porque a região para onde a população se moveu é uma em que as respostas antigas estão sistematicamente invertidas. O modelo monitorado dispara um alerta a cada dois meses, retreina com os rótulos daquele mês e fica a no máximo uns 0,02 da melhor acurácia possível o tempo todo.

Dois detalhes merecem atenção. O primeiro alerta (mês 2) dispara quando o modelo estático ainda vai bem (0,780, contra um melhor possível de 0,785): drift não é dano. Um alerta de PSI diz que as entradas mudaram, o que é motivo para olhar e talvez retreinar; não prova que o modelo quebrou. E o retreino precisou de rótulos. Aqui eles chegam ao fim de cada mês; com um atraso de 60 dias, o modelo monitorado estaria sempre dois meses atrás do drift.

Desalinhamento treino–serving e outras falhas de produção

O desalinhamento, medido

Etapa 5. Desalinhamento treino–serving (training–serving skew) é qualquer diferença entre como um atributo é calculado para o treino e como é calculado ao servir. Uma causa clássica é um novo aplicativo cliente que manda um campo em outra unidade. O conjunto de teste, de três jeitos:

base = loaded.predict_proba(X_test)[:, 1]
skewed = {'as in training': X_test,
          'amount in thousands': X_test.assign(credit_amount=X_test['credit_amount'] / 1000),
          'duration in years': X_test.assign(duration=X_test['duration'] / 12)}
for name, data in skewed.items():
    p = loaded.predict_proba(data)[:, 1]
    rejected = np.mean([bool(validate(r, schema)) for r in data.to_dict('records')])
    print(f"{name:20s} AUC {roc_auc_score(y_test, p):.3f} · flagged {np.mean(p > 0.5):5.1%}"
          f" · decisions changed {np.mean((p > 0.5) != (base > 0.5)):5.1%} · rejected {rejected:4.0%}")

as in training       AUC 0.783 · flagged 20.0% · decisions changed  0.0% · rejected   1%
amount in thousands  AUC 0.775 · flagged 14.8% · decisions changed  5.2% · rejected 100%
duration in years    AUC 0.751 · flagged  2.8% · decisions changed 17.2% · rejected  93%

O AUC quase não se mexe (0,783 contra 0,775 e 0,751), porque a ordenação é carregada principalmente pelos outros atributos. Um painel que acompanhasse o AUC, se houvesse rótulos, veria pouca coisa. As decisões contam outra história: com a duração em anos, 17,2% dos solicitantes recebem uma decisão diferente e a fração sinalizada como arriscada cai de 20,0% para 2,8%. A verificação do schema pega os dois bugs de uma vez, rejeitando todas as linhas em milhares e 93% das linhas em anos. Até o conjunto de teste intocado tem 1% de rejeições, porque alguns solicitantes reais caem fora da faixa de treino, um lembrete de que regras de rejeição exigem a decisão "rejeitar ou avisar" discutida acima.

Outras fontes comuns de desalinhamento:

• pré-processamento reimplementado na linguagem de serving (SQL, Java) em vez de distribuído dentro do artefato;
• atributos calculados a partir de uma tabela batch no treino, mas de um fluxo ao vivo no serving ("gasto médio nos últimos 30 dias" calculado com janelas ou fusos horários diferentes); feature stores como o Feast existem para calcular cada atributo uma única vez para os dois caminhos;
• atributos que estavam disponíveis na tabela histórica mas não no momento da predição, o vazamento de atributos que conhecem o futuro. Offline o modelo parece excelente; em produção a coluna ainda não existe. O lab mede um caso: a duração de uma ligação de vendas, conhecida só depois que a ligação termina.

Laços de realimentação

Colocando um modelo novo no ar com segurança

Um modelo retreinado só deve substituir o antigo depois de passar por um portão de avaliação: precisa vencer o modelo atual nos mesmos dados recentes e rotulados com que nenhum dos dois foi treinado, por mais do que o ruído (a comparação pareada da aula de validação). Depois ele é liberado em etapas:

Testes e documentação

Breck et al. (2017), em The ML Test Score, transformaram essas ideias numa lista de 28 testes em quatro grupos: dados, desenvolvimento do modelo, infraestrutura e monitoramento. Exemplos, parafraseados: as expectativas sobre os atributos estão registradas num schema; o treino é reprodutível; treino e serving calculam os mesmos valores de atributos; o modelo não está velho demais; a qualidade das predições nos dados servidos não regrediu. Cada teste vale meio ponto se é feito à mão e um ponto inteiro se é automatizado. A nota de um sistema é o mínimo entre os quatro grupos, já que a área mais fraca é a que causa o incidente. Quase tudo nesta aula corresponde a um desses testes: a ida e volta, o schema, as sementes, o hash dos dados, as verificações de PSI.

O próprio modelo precisa de documentação que viaje com ele. Um model card (Mitchell et al., 2019) é um documento curto com seções fixas: detalhes do modelo, uso pretendido, fatores (grupos e condições dos quais o desempenho pode depender), métricas, dados de avaliação, dados de treino, análises quantitativas, considerações éticas e, por fim, ressalvas e recomendações. O lab escreve um a partir dos seus próprios metadados e resultados de monitoramento.

Maturidade, com honestidade

A maioria das equipes não precisa de uma plataforma completa no primeiro dia. Uma escada sensata:

1. notebook → script + artefatos versionados: um pipeline salvo com seus metadados, um hash dos dados, um registro de experimentos;
2. deploy manual → CI/CD: testes (schema, ida e volta, um limiar de métrica num conjunto de validação fixo) rodam automaticamente antes de qualquer modelo ir para produção;
3. nenhuma visibilidade → logs + painéis + alertas de drift: toda predição registrada com a versão do modelo, PSI e distribuições de escore por lote;
4. retreino manual → retreino agendado ou disparado, com portões: um modelo novo só entra se vencer o atual em dados recentes, com rollback pronto.

Cada degrau elimina uma classe de incidente das 3 da manhã. Suba quando a dor justificar: um modelo que pontua uma vez por mês para um relatório interno não precisa de liberação canary.

O panorama de ferramentas

O ecossistema muda rápido; os conceitos acima, não. Um mapa, uma linha por ferramenta, para quando você precisar de uma:

Lab: colocando em produção e monitorando um modelo de marketing bancário

O cenário é real. Um banco português fez campanhas de telemarketing de um depósito a prazo entre maio de 2008 e novembro de 2010. O conjunto de dados bank-marketing da UCI registra 45.211 dessas ligações na ordem em que foram feitas: atributos do cliente, histórico de contatos e se o cliente contratou. Fazemos o papel da equipe de dados do banco: em 1º de maio de 2009, um modelo treinado com o ano anterior é implantado para ordenar os clientes antes de cada ligação; nós o acompanhamos pelos 19 meses seguintes. O drift que acontece no caminho não é simulado. Ele aconteceu, durante a crise financeira.

Os passos 1–3 constroem e validam o modelo, os passos 4–6 o empacotam e põem uma função de predição protegida na frente dele, os passos 7–9 o monitoram em produção, o passo 10 o retreina e o passo 11 escreve seu model card. Arquivos são gravados apenas numa pasta temporária. Rode os blocos em ordem, numa mesma sessão; toda saída abaixo é o que o código imprime (em poucos segundos, depois que os dados são baixados). Para trabalhar direto nele, baixe o notebook do lab e abra no Jupyter ou no Colab (Arquivo → Fazer upload de notebook).

Passo 1 — carregar as ligações, em ordem temporal

import numpy as np
import pandas as pd
from sklearn.datasets import fetch_openml

raw = fetch_openml(data_id=1461, as_frame=True).frame        # UCI bank marketing, 45.211 ligações
names = ['age', 'job', 'marital', 'education', 'default', 'balance', 'housing', 'loan',
         'contact', 'day', 'month', 'duration', 'campaign', 'pdays', 'previous', 'poutcome']
calls = raw.drop(columns='Class').set_axis(names, axis=1)
calls['y'] = (raw['Class'] == '2').astype(int)               # 1 = contratou um depósito a prazo

# as linhas estão na ordem das ligações, de maio de 2008 a novembro de 2010, mas o ano não
# é armazenado: um ano novo começa cada vez que começa uma sequência de ligações em 'jan'
month = calls['month'].astype(str)
new_run = month.ne(month.shift())
year = 2008 + (new_run & month.eq('jan')).cumsum()
calls['date'] = pd.to_datetime(year.astype(str) + '-' + month, format='%Y-%b').dt.to_period('M')

print(calls.shape, calls['date'].min(), '->', calls['date'].max())
print(calls.groupby(calls['date'].dt.year)['y'].agg(['size', 'mean']).round(3))

(45211, 18) 2008-05 -> 2010-11
       size   mean
date
2008  27729  0.051
2009  14862  0.171
2010   2620  0.516

A cópia do OpenML chama as colunas de V1 a V16, então restauramos os nomes a partir da documentação do conjunto de dados. A tabela tem o mês mas não o ano; como as linhas estão na ordem das ligações, dá para reconstruir o ano contando os janeiros. O fato mais importante do lab já aparece aqui: 5,1% das ligações de 2008 terminaram em contratação, 17,1% em 2009 e 51,6% em 2010, com muito menos ligações. O banco mudou para quem ligava e como, no meio de uma crise.

Passo 2 — a coluna que não vai existir quando o modelo for usado

from sklearn.compose import make_column_transformer
from sklearn.linear_model import LogisticRegression
from sklearn.metrics import roc_auc_score
from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import OneHotEncoder, StandardScaler

CAT = ['job', 'marital', 'education', 'default', 'housing', 'loan', 'contact', 'poutcome']
NUM = ['age', 'balance', 'campaign', 'pdays', 'previous']
for c in CAT:
    calls[c] = calls[c].astype(str)

def make_model(num=NUM):
    return make_pipeline(
        make_column_transformer((OneHotEncoder(handle_unknown='ignore'), CAT),
                                (StandardScaler(), num)),
        LogisticRegression(max_iter=2000))

history = calls[calls['date'] < '2009-05']                 # o ano anterior à implantação
fit_part = history[history['date'] < '2009-02']            # maio de 2008 - jan de 2009
val_part = history[history['date'] >= '2009-02']           # fev - abr de 2009, os mais recentes
for label, num in [('with duration', NUM + ['duration']), ('without duration', NUM)]:
    m = make_model(num).fit(fit_part[CAT + num], fit_part['y'])
    auc = roc_auc_score(val_part['y'], m.predict_proba(val_part[CAT + num])[:, 1])
    print(f"{label:17s} validation AUC {auc:.3f}")

with duration     validation AUC 0.757
without duration  validation AUC 0.619

A validação é temporal: ajuste em maio de 2008 – janeiro de 2009, validação nos três meses mais recentes, como a aula de validação recomenda para dados com ordem no tempo. duration é a duração da ligação de vendas, em segundos. Ela leva o AUC de 0,619 para 0,757 e, ainda assim, é inútil: o modelo precisa ordenar os clientes antes da ligação, quando a duração não é conhecida. A própria documentação do conjunto de dados alerta sobre isso; o lab de Gradient Boosting mediu o mesmo vazamento numa divisão aleatória. É o vazamento de um atributo que conhece o futuro. Em produção, ele não provocaria erro: a coluna simplesmente estaria ausente ou preenchida com um valor padrão. Nós a removemos, assim como day e month, que registram quando a ligação aconteceu e não quem é o cliente.

O AUC honesto, 0,619, é modesto. Prever quem vai contratar um depósito a partir da papelada de uma ligação é difícil. Não vamos tentar melhorar o modelo neste lab: o assunto é o que acontece com ele depois da implantação.

Passo 3 — treinar o modelo que o banco implanta em maio de 2009

from sklearn.metrics import brier_score_loss

FEATURES = CAT + NUM
m = make_model().fit(fit_part[FEATURES], fit_part['y'])
p_val = m.predict_proba(val_part[FEATURES])[:, 1]
validation = {'auc': round(float(roc_auc_score(val_part['y'], p_val)), 3),
              'brier': round(float(brier_score_loss(val_part['y'], p_val)), 3),
              'mean_predicted': round(float(p_val.mean()), 3),
              'actual_rate': round(float(val_part['y'].mean()), 3)}
print('validation (Feb-Apr 2009):', validation)

model = make_model().fit(history[FEATURES], history['y'])   # reajuste no ano inteiro
print(f"deployed model trained on {len(history):,} calls, {history['y'].mean():.1%} subscribed")

validation (Feb-Apr 2009): {'auc': 0.619, 'brier': 0.141, 'mean_predicted': 0.055, 'actual_rate': 0.158}
deployed model trained on 34,177 calls, 6.6% subscribed

O registro de validação tem quatro números porque eles falham de jeitos diferentes. O AUC (0,619) avalia a ordenação. O escore de Brier (0,141) e a comparação da probabilidade média prevista com a taxa observada avaliam as probabilidades. E ali o aviso já está visível antes da implantação: em fevereiro–abril de 2009 o modelo previu 5,5% de contratações em média, enquanto 15,8% contrataram. O mundo já estava se movendo dentro do ano de treino. O modelo implantado é reajustado no ano inteiro, 34.177 ligações com 6,6% de contratações.

Passo 4 — empacotar: o modelo, seus metadados, uma pasta versionada

import hashlib
import json
import os
import platform
import tempfile

import joblib
import sklearn

def data_hash(df):
    return hashlib.sha256(pd.util.hash_pandas_object(df, index=True).to_numpy().tobytes()).hexdigest()[:12]

def make_schema(X, limits):
    """Faixas e categorias do treino, mais limites rígidos que vêm do negócio."""
    schema = {}
    for col in X.columns:
        if pd.api.types.is_numeric_dtype(X[col]):
            schema[col] = {'type': 'number', 'min': float(X[col].min()), 'max': float(X[col].max())}
            if col in limits:
                schema[col]['hard_min'], schema[col]['hard_max'] = limits[col]
        else:
            schema[col] = {'type': 'category', 'values': sorted(X[col].unique())}
    return schema

store = tempfile.mkdtemp()                                   # faz o papel de um registro de modelos
version = 'bank-2009.05'
os.makedirs(os.path.join(store, version))
joblib.dump(model, os.path.join(store, version, 'model.joblib'))
metadata = {'version': version, 'trained_on': '2008-05 to 2009-04', 'rows': len(history),
            'data_hash': data_hash(history[FEATURES + ['y']]),
            'python': platform.python_version(), 'sklearn': sklearn.__version__,
            'features': FEATURES, 'schema': make_schema(history[FEATURES], {'age': (18, 100)}),
            'validation': validation}
with open(os.path.join(store, version, 'metadata.json'), 'w') as f:
    json.dump(metadata, f, indent=1)

print(sorted(os.listdir(os.path.join(store, version))))
print({k: metadata[k] for k in ['version', 'rows', 'data_hash', 'sklearn']})
print('pdays:', metadata['schema']['pdays'])
print('age:  ', metadata['schema']['age'])

['metadata.json', 'model.joblib']
{'version': 'bank-2009.05', 'rows': 34177, 'data_hash': 'b18de86a57b0', 'sklearn': '1.9.0'}
pdays: {'type': 'number', 'min': -1.0, 'max': 350.0}
age:   {'type': 'number', 'min': 19.0, 'max': 95.0, 'hard_min': 18, 'hard_max': 100}

O artefato é uma pasta com o nome da versão, contendo o pipeline e um arquivo JSON com tudo o que é preciso para confiar nele depois: o período de treino e um hash das linhas de treino, as versões das bibliotecas, a lista de atributos, o schema e as métricas de validação. O schema tem dois tipos de faixa. A faixa de treino vem dos dados (pdays, o número de dias desde o último contato com o cliente, nunca passou de 350 no ano de treino). Os limites rígidos vêm do negócio: um cliente tem entre 18 e 100 anos, digam os dados o que disserem. store faz o papel de um registro de modelos; num sistema real seria o MLflow, um bucket ou uma pasta compartilhada em que só o job de treino pode escrever.

Passo 5 — recarregar e verificar

import warnings

def load_artifact(folder):
    with open(os.path.join(folder, 'metadata.json')) as f:
        meta = json.load(f)
    if meta['sklearn'] != sklearn.__version__:
        warnings.warn(f"model built with scikit-learn {meta['sklearn']}, running {sklearn.__version__}")
    return joblib.load(os.path.join(folder, 'model.joblib')), meta

served, meta = load_artifact(os.path.join(store, version))
check = history[FEATURES].sample(2000, random_state=0)
same = np.array_equal(model.predict_proba(check), served.predict_proba(check))
print('reloaded model gives identical probabilities:', same)
print('data hash matches the metadata:', data_hash(history[FEATURES + ['y']]) == meta['data_hash'])

reloaded model gives identical probabilities: True
data hash matches the metadata: True

load_artifact é o que o processo de serving roda ao iniciar. Ele confere a versão da biblioteca com os metadados antes de confiar no pickle. Os dois prints são os testes que deveriam rodar automaticamente antes de qualquer implantação: o modelo recarregado reproduz exatamente o modelo em memória; os dados com que ele diz ter sido treinado são os que temos.

Passo 6 — uma função de predição protegida

def validate(payload, schema):
    """Erros rejeitam a requisição; avisos a deixam passar, mas são registrados."""
    errors = [f"unexpected field {k!r}" for k in payload if k not in schema]
    warns = []
    for col, rule in schema.items():
        v = payload.get(col)
        if v is None:
            errors.append(f"{col}: missing")
        elif rule['type'] == 'category':
            if v not in rule['values']:
                errors.append(f"{col}: unknown category {v!r}")
        elif isinstance(v, bool) or not isinstance(v, (int, float)):
            errors.append(f"{col}: expected a number, got {v!r}")
        elif not rule.get('hard_min', -np.inf) <= v <= rule.get('hard_max', np.inf):
            errors.append(f"{col}: {v} is impossible")
        elif not rule['min'] <= v <= rule['max']:
            warns.append(f"{col}: {v} outside training range [{rule['min']:g}, {rule['max']:g}]")
    return errors, warns

def predict(payload, model=served, meta=meta):
    errors, warns = validate(payload, meta['schema'])
    if errors:
        return {'status': 'rejected', 'errors': errors, 'model_version': meta['version']}
    p = float(model.predict_proba(pd.DataFrame([payload])[meta['features']])[0, 1])
    return {'status': 'ok', 'p_subscribe': round(p, 3), 'warnings': warns,
            'model_version': meta['version']}

client = {'job': 'technician', 'marital': 'married', 'education': 'secondary', 'default': 'no',
          'housing': 'yes', 'loan': 'no', 'contact': 'cellular', 'poutcome': 'unknown',
          'age': 41, 'balance': 1200, 'campaign': 2, 'pdays': -1, 'previous': 0}
tests = {'good client': client,
         'missing balance': {k: v for k, v in client.items() if k != 'balance'},
         'age = 212': {**client, 'age': 212},
         'unseen job': {**client, 'job': 'influencer'},
         'age as text': {**client, 'age': '41'},
         'pdays = 400': {**client, 'pdays': 400, 'previous': 1, 'poutcome': 'failure'}}
for name, payload in tests.items():
    answer = predict(payload)
    detail = answer.get('errors') or answer.get('warnings') or ''
    print(f"{name:16s} {answer['status']:8s} {answer.get('p_subscribe', ''):<6} {detail}")

good client      ok       0.063
missing balance  rejected        ['balance: missing']
age = 212        rejected        ['age: 212 is impossible']
unseen job       rejected        ["job: unknown category 'influencer'"]
age as text      rejected        ["age: expected a number, got '41'"]
pdays = 400      ok       0.024  ['pdays: 400 outside training range [-1, 350]']

Esta é a função que um serviço web chamaria (o esboço em FastAPI de Servindo predições a envolve sem mudanças). Seis requisições testam esse contrato. O cliente bom recebe uma probabilidade. Um campo ausente, uma idade impossível, uma profissão desconhecida e um número enviado como texto são rejeitados com um motivo legível. O último cliente foi contatado há 400 dias, mais tempo do que qualquer um no ano de treino, então a requisição é servida com um aviso. É a escolha certa aqui: como o próximo passo mostra, valores de pdays acima de 350 não são bug, são o que o tempo faz com este atributo.

Passo 7 — produção: monitores que não precisam de rótulos

A partir de maio de 2009 o modelo pontua toda ligação. A produção é cortada em trimestres do calendário; para cada um calculamos o que está disponível na hora: quantas requisições o validador rejeitaria ou deixaria passar com aviso; o PSI de cada atributo em relação ao ano de treino; e o PSI e a estatística KS do próprio escore do modelo.

from scipy.stats import ks_2samp

def psi(ref, cur, bins=10, eps=1e-4):
    """PSI em faixas de quantis de `ref` (números) ou nas categorias (texto)."""
    ref, cur = pd.Series(ref), pd.Series(cur)
    if pd.api.types.is_numeric_dtype(ref):
        edges = np.unique(np.quantile(ref, np.linspace(0, 1, bins + 1)))
        edges[0], edges[-1] = -np.inf, np.inf
        e = np.histogram(ref, edges)[0] / len(ref)
        a = np.histogram(cur, edges)[0] / len(cur)
    else:
        cats = sorted(set(ref) | set(cur))
        e = ref.value_counts(normalize=True).reindex(cats, fill_value=0).to_numpy()
        a = cur.value_counts(normalize=True).reindex(cats, fill_value=0).to_numpy()
    e, a = np.clip(e, eps, None), np.clip(a, eps, None)
    return float(np.sum((a - e) * np.log(a / e)))

production = calls[calls['date'] >= '2009-05'].copy()
production['batch'] = (production['date'].dt.qyear.astype(str) + 'Q'
                       + production['date'].dt.quarter.astype(str))
ref_score = served.predict_proba(history[FEATURES])[:, 1]

rows = []
for batch, g in production.groupby('batch'):
    checks = [validate(r, meta['schema']) for r in g[FEATURES].to_dict('records')]
    score = served.predict_proba(g[FEATURES])[:, 1]
    feature_psi = {c: psi(history[c], g[c]) for c in FEATURES}
    worst = max(feature_psi, key=feature_psi.get)
    rows.append({'batch': batch, 'calls': len(g),
                 'rejected': np.mean([bool(e) for e, w in checks]),
                 'warned': np.mean([bool(w) for e, w in checks]),
                 'worst feature': worst, 'its PSI': feature_psi[worst],
                 'score PSI': psi(ref_score, score),
                 'score KS D': ks_2samp(ref_score, score).statistic})
monitor = pd.DataFrame(rows).set_index('batch')
print(monitor.round(3))

        calls  rejected  warned worst feature  its PSI  score PSI  score KS D
batch
2009Q2   6217       0.0   0.147       contact    3.247      0.071       0.115
2009Q3   1261       0.0   0.053       contact    2.212      1.360       0.475
2009Q4    936       0.0   0.017       contact    2.112      1.390       0.446
2010Q1    799       0.0   0.041      poutcome    2.212      1.218       0.447
2010Q2    661       0.0   0.120      poutcome    2.145      0.891       0.391
2010Q3    865       0.0   0.148      poutcome    2.114      0.810       0.371
2010Q4    295       0.0   0.122      poutcome    2.295      1.098       0.439

Como ler, coluna por coluna:

• rejected é zero em todos os trimestres: nenhum campo ausente, nenhuma categoria desconhecida, nenhum valor impossível. Os sistemas que alimentam o modelo estão saudáveis;
• warned não é zero. Em 2009Q2, 14,7% das ligações têm um valor fora da faixa de treino, quase todos de pdays acima de 350. Os dados de treino começaram em maio de 2008, então ninguém neles poderia ter sido contatado mais de um ano antes; um ano depois, muitos clientes foram. A faixa do atributo cresce com o tempo: uma regra de rejeição rígida teria recusado uma ligação em cada sete;
• worst feature é contact em 2009, com PSI acima de 2. Em maio e junho de 2008 o canal de contato não era registrado (todo valor é unknown); a partir de julho de 2008, passou a ser. Uma mudança em como os dados são registrados é indistinguível de um drift. Um modelo treinado com a mistura de 2008 aprendeu, em parte, "canal desconhecido" como sinal. A partir de 2010 o pior atributo é poutcome, o desfecho da campanha anterior: cada vez mais clientes já tinham recebido ligações antes e muitos deles tinham dito sim;
• score PSI é pequeno em 2009Q2 (0,071) e depois salta para acima de 0,8 em todos os trimestres seguintes: a distribuição das saídas do próprio modelo se moveu, antes de chegar qualquer rótulo. É o monitor isolado mais barato e mais útil.

Passo 8 — os rótulos chegam: quão bom foi o modelo?

O desfecho de uma ligação é conhecido em poucas semanas, então ao fim de cada trimestre os rótulos já chegaram e as métricas de verdade podem ser calculadas:

perf = []
for batch, g in production.groupby('batch'):
    score = served.predict_proba(g[FEATURES])[:, 1]
    perf.append({'batch': batch, 'actual rate': g['y'].mean(), 'mean predicted': score.mean(),
                 'AUC': roc_auc_score(g['y'], score), 'Brier': brier_score_loss(g['y'], score)})
perf = pd.DataFrame(perf).set_index('batch')
print(perf.round(3))

        actual rate  mean predicted    AUC  Brier
batch
2009Q2        0.124           0.081  0.671  0.105
2009Q3        0.367           0.150  0.584  0.276
2009Q4        0.459           0.158  0.609  0.326
2010Q1        0.517           0.174  0.664  0.347
2010Q2        0.551           0.167  0.698  0.372
2010Q3        0.519           0.157  0.710  0.350
2010Q4        0.427           0.174  0.749  0.271

Os dois tipos de métrica contam histórias diferentes. O AUC fica entre 0,584 e 0,749, nunca longe do valor de validação de 0,619 (e até acima dele no fim): o modelo ainda ordena os clientes mais ou menos tão bem quanto antes. As probabilidades desabam. A fração de contratantes vai de 12,4% para acima de 50%, enquanto o modelo continua prevendo uns 16%; em 2010Q2 ele espera 16,7 contratações a cada 100 ligações onde acontecem 55,1. O escore de Brier vai de 0,105 para 0,372. Isso é principalmente prior shift, a taxa de base se movendo, somado aos atributos de histórico de contato mudando de significado. Uma equipe que monitorasse só o AUC concluiria que nada aconteceu; uma equipe que planejasse pessoal, orçamento ou volume de ligações a partir das probabilidades previstas erraria por um fator de três.

Passo 9 — regras de alerta

Monitores só servem se alguém for avisado quando disparam. As regras são escritas de antemão, com limiares; aqui, uma por monitor:

alerts = pd.DataFrame({
    'rejections > 1%': monitor['rejected'] > 0.01,
    'warnings > 5%': monitor['warned'] > 0.05,
    'feature PSI > 0.25': monitor['its PSI'] > 0.25,
    'score PSI > 0.25': monitor['score PSI'] > 0.25,
    'calibration gap > 0.05': (perf['mean predicted'] - perf['actual rate']).abs() > 0.05,
    'AUC < validation - 0.05': perf['AUC'] < validation['auc'] - 0.05,
})
print(alerts.map(lambda fired: 'ALERT' if fired else '.').T.to_string())

batch                   2009Q2 2009Q3 2009Q4 2010Q1 2010Q2 2010Q3 2010Q4
rejections > 1%              .      .      .      .      .      .      .
warnings > 5%            ALERT  ALERT      .      .  ALERT  ALERT  ALERT
feature PSI > 0.25       ALERT  ALERT  ALERT  ALERT  ALERT  ALERT  ALERT
score PSI > 0.25             .  ALERT  ALERT  ALERT  ALERT  ALERT  ALERT
calibration gap > 0.05       .  ALERT  ALERT  ALERT  ALERT  ALERT  ALERT
AUC < validation - 0.05      .      .      .      .      .      .      .

Leia como uma linha do tempo. No primeiro trimestre de produção, só os monitores de dados, que dispensam rótulos, disparam: os avisos de fora da faixa (o pdays crescente) e o PSI dos atributos (a mudança de registro em contact). A distribuição do escore e a calibração ainda são aceitáveis; quem olhasse esses alertas em meados de 2009 encontraria duas mudanças explicáveis, não um modelo quebrado. Do segundo trimestre em diante, o PSI do escore dispara; quando os rótulos chegam, o desvio de calibração o confirma. A regra do AUC nunca dispara. Escolher a métrica a monitorar é escolher o que você é capaz de perceber.

Passo 10 — retreinar com dados recentes e comparar

A correção para uma taxa de base que se move é aprendê-la de novo. Uma política simples: no início de cada trimestre, retreinar o mesmo pipeline com os últimos seis meses de ligações rotuladas e compará-lo, no trimestre, com o modelo estático implantado em maio de 2009.

from sklearn.base import clone

compare = []
for batch, g in production.groupby('batch'):
    start = g['date'].min()
    recent = calls[(calls['date'] < start) & (calls['date'] >= start - 6)]   # últimos 6 meses rotulados
    retrained = clone(model).fit(recent[FEATURES], recent['y'])
    s_old = served.predict_proba(g[FEATURES])[:, 1]
    s_new = retrained.predict_proba(g[FEATURES])[:, 1]
    compare.append({'batch': batch, 'window': len(recent),
                    'AUC static': roc_auc_score(g['y'], s_old),
                    'AUC retrained': roc_auc_score(g['y'], s_new),
                    'Brier static': brier_score_loss(g['y'], s_old),
                    'Brier retrained': brier_score_loss(g['y'], s_new),
                    'gap static': s_old.mean() - g['y'].mean(),
                    'gap retrained': s_new.mean() - g['y'].mean()})
compare = pd.DataFrame(compare).set_index('batch')
print(compare.round(3).to_string())
print(compare.drop(columns='window').mean().round(3).to_string())