Predicción de ventas, segmentación de clientes y detección de anomalías en retail automotriz

pip install pandas numpy scikit-learn xgboost statsmodels scipy matplotlib seaborn holidays openpyxl

proyecto-retail/
├── data/
│   ├── raw/                        # datos originales del cliente (confidencial)
│   └── data_productos_sint.csv     # datos sintéticos generados
├── notebooks/
│   ├── 01_generacion_sintetica.ipynb
│   ├── 02_prediccion_ventas.ipynb
│   ├── 03_anomalias.ipynb
│   └── 04_segmentacion_rfm.ipynb
├── src/
│   ├── data_generator.py
│   ├── forecasting.py
│   ├── anomaly_detection.py
│   └── segmentation.py
└── requirements.txt

import pandas as pd
import holidays

fechas = pd.date_range(start="2020-01-01", end="2024-12-31", freq="D")
festivos_peru = holidays.country_holidays("PE", years=range(2020, 2025))

fechas_validas = [
    f for f in fechas
    if f.weekday() != 6
    and f not in festivos_peru
]

print(f"Días hábiles 2020-2024: {len(fechas_validas)}")

import numpy as np
import pandas as pd

synthetic_data["Fecha_Emision"] = pd.to_datetime(synthetic_data["Fecha_Emision"])

inicio_pandemia = pd.to_datetime("2020-03-01")
fin_pandemia    = pd.to_datetime("2021-03-31")

mascara = (
    (synthetic_data["Fecha_Emision"] >= inicio_pandemia) &
    (synthetic_data["Fecha_Emision"] <= fin_pandemia)
)

factor_pandemia = np.random.uniform(0.3, 0.6, mascara.sum())
synthetic_data.loc[mascara, "Cantidad"] = (
    synthetic_data.loc[mascara, "Cantidad"] * factor_pandemia
).astype(int)

import pandas as pd

data = pd.read_csv("data/data_productos_sint.csv")
data["Fecha_Emision"] = pd.to_datetime(data["Fecha_Emision"])

data_productos = data[["Fecha_Emision", "Codigo_Producto", "Cantidad"]].copy()

productos_agrupados = (
    data_productos
    .groupby("Codigo_Producto")["Cantidad"]
    .sum()
    .reset_index()
    .sort_values(by="Cantidad", ascending=False)
)

top10 = productos_agrupados.head(10)
codigos_top10 = top10["Codigo_Producto"].tolist()
print(f"Top 10 productos: {codigos_top10}")

from statsmodels.graphics.tsaplots import plot_acf
import matplotlib.pyplot as plt

data_productos["Mes"] = data_productos["Fecha_Emision"].dt.to_period("M")
diccionario_dataframes = {}

for codigo in codigos_top10:
    df_prod = (
        data_productos[data_productos["Codigo_Producto"] == codigo]
        .groupby("Mes")["Cantidad"]
        .sum()
        .reset_index()
    )
    df_prod.columns = ["Fecha", "Cantidad"]
    diccionario_dataframes[codigo] = df_prod

fig, axes = plt.subplots(5, 2, figsize=(14, 18))
for ax, (codigo, df) in zip(axes.flatten(), diccionario_dataframes.items()):
    plot_acf(df["Cantidad"], lags=12, ax=ax, title=f"ACF — {codigo}")
plt.tight_layout()
plt.show()

NUM_LAGS = 5

for codigo, df in diccionario_dataframes.items():
    for lag in range(1, NUM_LAGS + 1):
        df[f"lag_{lag}"] = df["Cantidad"].shift(lag)
    df.fillna(0, inplace=True)
    diccionario_dataframes[codigo] = df

print(diccionario_dataframes[codigos_top10[0]].head(7))

from datetime import timedelta
from scipy import stats
import pandas as pd

df_clientes = data[["Fecha_Emision", "Cliente", "Valor_Venta", "Documento"]].copy()

snapshot_date = df_clientes["Fecha_Emision"].max() + timedelta(days=1)

customers = df_clientes.groupby("Cliente").agg(
    Ultimo    = ("Fecha_Emision", lambda x: (snapshot_date - x.max()).days),
    Frecuencia= ("Documento",     "nunique"),
    Monto     = ("Valor_Venta",   "sum")
).reset_index()

customers_fix = pd.DataFrame()
customers_fix["Ultimo"]     = stats.boxcox(customers["Ultimo"])[0]
customers_fix["Frecuencia"] = stats.boxcox(customers["Frecuencia"])[0]
customers_fix["Monto"]      = stats.boxcox(customers["Monto"])[0]

print(customers_fix.describe())

from sklearn.tree import DecisionTreeRegressor
from sklearn.ensemble import GradientBoostingRegressor, AdaBoostRegressor, RandomForestRegressor
from sklearn.model_selection import train_test_split
from sklearn.metrics import mean_absolute_percentage_error
from xgboost import XGBRegressor
import numpy as np

modelos = {
    "Decision Tree":      DecisionTreeRegressor(random_state=42),
    "Random Forest":      RandomForestRegressor(random_state=42),
    "Gradient Boosting":  GradientBoostingRegressor(random_state=42),
    "XGBoost":            XGBRegressor(random_state=42, verbosity=0),
    "AdaBoost":           AdaBoostRegressor(random_state=42),
}

errores_mape = {}

for codigo, df in diccionario_dataframes.items():
    features = [f"lag_{i}" for i in range(1, NUM_LAGS + 1)]
    X = df[features]
    y = df["Cantidad"]

    X_train, X_test, y_train, y_test = train_test_split(
        X, y, test_size=0.2, shuffle=False
    )

    errores_producto = {}
    for nombre, modelo in modelos.items():
        modelo.fit(X_train, y_train)
        y_pred = modelo.predict(X_test)
        mape = mean_absolute_percentage_error(y_test, y_pred) * 100
        errores_producto[nombre] = round(mape, 2)

    errores_mape[codigo] = errores_producto

import pandas as pd
df_errores = pd.DataFrame(errores_mape).T
print("MAPE promedio por modelo:")
print(df_errores.mean().sort_values())

diccionario_anomalias_iqr = {}

for codigo, df in diccionario_dataframes.items():
    Q1  = df["Cantidad"].quantile(0.25)
    Q3  = df["Cantidad"].quantile(0.75)
    IQR = Q3 - Q1

    limite_inferior = Q1 - 1.5 * IQR
    limite_superior = Q3 + 1.5 * IQR

    anomalias = df[
        (df["Cantidad"] < limite_inferior) |
        (df["Cantidad"] > limite_superior)
    ]
    diccionario_anomalias_iqr[codigo] = anomalias
    print(f"{codigo}: {len(anomalias)} anomalías detectadas con IQR")

import numpy as np
from sklearn.cluster import DBSCAN
from sklearn.neighbors import NearestNeighbors
from sklearn.preprocessing import StandardScaler

def calcular_epsilon(X_scaled, k=5):
    nbrs = NearestNeighbors(n_neighbors=k).fit(X_scaled)
    distancias, _ = nbrs.kneighbors(X_scaled)
    dist_k = np.sort(distancias[:, -1])
    diff = np.diff(dist_k)
    epsilon = dist_k[np.argmax(diff)]
    return epsilon

diccionario_anomalias_dbscan = {}

for codigo, df in diccionario_dataframes.items():
    X = df[["Cantidad"]].values
    scaler = StandardScaler()
    X_scaled = scaler.fit_transform(X)

    epsilon = calcular_epsilon(X_scaled)
    db = DBSCAN(eps=epsilon, min_samples=3).fit(X_scaled)

    df_result = df.copy()
    df_result["cluster"] = db.labels_
    anomalias = df_result[df_result["cluster"] == -1]
    diccionario_anomalias_dbscan[codigo] = anomalias
    print(f"{codigo}: {len(anomalias)} anomalías con DBSCAN (ε={epsilon:.4f})")

from sklearn.cluster import KMeans
from sklearn.preprocessing import StandardScaler
from sklearn.metrics import silhouette_score
import matplotlib.pyplot as plt

scaler = StandardScaler()
customers_normalized = scaler.fit_transform(customers_fix)

sse = {}
for k in range(1, 11):
    kmeans = KMeans(n_clusters=k, random_state=42, n_init=10)
    kmeans.fit(customers_normalized)
    sse[k] = kmeans.inertia_

plt.figure(figsize=(8, 4))
plt.plot(list(sse.keys()), list(sse.values()), marker="o")
plt.xlabel("Número de clusters (k)")
plt.ylabel("SSE")
plt.title("Método del codo — Segmentación RFM")
plt.grid(True, alpha=0.3)
plt.show()

for k in range(2, 6):
    km = KMeans(n_clusters=k, random_state=42, n_init=10).fit(customers_normalized)
    score = silhouette_score(customers_normalized, km.labels_)
    print(f"k={k}: Silhouette = {score:.4f}")

model = KMeans(n_clusters=3, random_state=42, n_init=10)
model.fit(customers_normalized)
customers["Cluster"] = model.labels_

perfil_clusters = customers.groupby("Cluster").agg(
    Ultimo_promedio    = ("Ultimo",     "mean"),
    Frecuencia_promedio= ("Frecuencia", "mean"),
    Monto_promedio     = ("Monto",      "mean"),
    N_clientes         = ("Cliente",    "count")
).round(1)

print(perfil_clusters)
# Cluster 0: Reciente, alta frecuencia, alto monto  → Clientes VIP
# Cluster 1: Alejados, baja frecuencia, bajo monto  → En riesgo de churn
# Cluster 2: Frecuencia media, monto medio          → Clientes regulares