Aula 9: Machine Learning para Inferência Causal

# Implementação de Causal Forest em Python usando EconML
import numpy as np
import pandas as pd
from econml.dml import CausalForestDML
from sklearn.ensemble import RandomForestRegressor, RandomForestClassifier
from sklearn.model_selection import train_test_split

def causal_forest_analysis(df, X_cols, y_col, treatment_col, 
                         heterogeneity_cols=None, n_estimators=100):
    """
    Implementa Causal Forest para estimar efeitos heterogêneos
    
    Args:
        df: DataFrame com os dados
        X_cols: Lista de colunas para covariáveis
        y_col: Nome da coluna de resultado
        treatment_col: Nome da coluna de tratamento
        heterogeneity_cols: Colunas para analisar heterogeneidade (opcional)
        n_estimators: Número de árvores na floresta
    
    Returns:
        Modelo treinado e efeitos heterogêneos
    """
    # Preparar dados
    X = df[X_cols].values
    T = df[treatment_col].values.reshape(-1, 1)
    y = df[y_col].values
    
    # Dividir dados para honestidade
    X_train, X_test, T_train, T_test, y_train, y_test = train_test_split(
        X, T, y, test_size=0.5, random_state=42)
    
    # Configurar e treinar o modelo
    cf = CausalForestDML(
        model_y=RandomForestRegressor(n_estimators=n_estimators),
        model_t=RandomForestClassifier(n_estimators=n_estimators),
        n_estimators=n_estimators,
        min_samples_leaf=5,
        max_depth=None,
        verbose=0,
        random_state=42
    )
    
    cf.fit(y=y_train, T=T_train, X=X_train)
    
    # Estimar efeitos individuais
    te_pred = cf.effect(X_test)
    
    # Análise de heterogeneidade
    heterogeneity_results = {}
    if heterogeneity_cols:
        for col in heterogeneity_cols:
            if col in df.columns:
                # Agrupar efeitos por variável de interesse
                grouped = pd.DataFrame({
                    'effect': te_pred.flatten(),
                    'group': df.loc[X_test.index, col]
                }).groupby('group')['effect'].mean()
                heterogeneity_results[col] = grouped
    
    return cf, te_pred, heterogeneity_results

# Implementação de Double Machine Learning em Python
import numpy as np
import pandas as pd
from sklearn.model_selection import KFold
from sklearn.ensemble import RandomForestRegressor, RandomForestClassifier
from sklearn.linear_model import LinearRegression

def double_ml(df, X_cols, y_col, d_col, n_splits=5, ml_model=None):
    """
    Implementa Double Machine Learning para inferência causal
    
    Args:
        df: DataFrame com os dados
        X_cols: Lista de colunas para covariáveis
        y_col: Nome da coluna de resultado
        d_col: Nome da coluna de tratamento
        n_splits: Número de folds para cross-fitting
        ml_model: Modelo ML a ser usado (default: RandomForest)
    
    Returns:
        Estimativa do parâmetro causal e estatísticas de inferência
    """
    # Preparar dados
    X = df[X_cols].values
    y = df[y_col].values
    d = df[d_col].values
    
    # Definir modelo ML padrão se não especificado
    if ml_model is None:
        ml_model_y = RandomForestRegressor(n_estimators=100, min_samples_leaf=5)
        ml_model_d = RandomForestRegressor(n_estimators=100, min_samples_leaf=5)
    else:
        ml_model_y = ml_model
        ml_model_d = ml_model
    
    # Configurar K-fold cross-fitting
    kf = KFold(n_splits=n_splits, shuffle=True, random_state=42)
    
    # Armazenar resultados por fold
    thetas = []
    
    # Loop através dos folds
    for train_idx, test_idx in kf.split(X):
        # Dividir dados
        X_train, X_test = X[train_idx], X[test_idx]
        y_train, y_test = y[train_idx], y[test_idx]
        d_train, d_test = d[train_idx], d[test_idx]
        
        # Primeira etapa: modelagem de Y
        ml_model_y.fit(X_train, y_train)
        y_pred = ml_model_y.predict(X_test)
        
        # Primeira etapa: modelagem de D
        ml_model_d.fit(X_train, d_train)
        d_pred = ml_model_d.predict(X_test)
        
        # Segunda etapa: calcular resíduos
        y_resid = y_test - y_pred
        d_resid = d_test - d_pred
        
        # Estimar theta para este fold
        theta = np.cov(y_resid, d_resid)[0, 1] / np.var(d_resid)
        thetas.append(theta)
    
    # Calcular estimativa final e estatísticas
    theta_hat = np.mean(thetas)
    theta_std = np.std(thetas) / np.sqrt(n_splits)
    theta_ci_lower = theta_hat - 1.96 * theta_std
    theta_ci_upper = theta_hat + 1.96 * theta_std
    
    results = {
        'theta_hat': theta_hat,
        'std_error': theta_std,
        'conf_int_95': [theta_ci_lower, theta_ci_upper],
        'fold_estimates': thetas
    }
    
    return results

# Análise de Efeitos Heterogêneos com Causal Forest em Python
# Rebaixamento de Rating de Crédito e Impacto no Valor das Empresas

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from econml.dml import CausalForestDML
from sklearn.ensemble import RandomForestRegressor, RandomForestClassifier
from sklearn.model_selection import train_test_split
from sklearn.inspection import partial_dependence

# 1. Carregar e preparar os dados
# Suponha que temos um DataFrame com empresas, características, tratamento (rebaixamento) e resultado (valor)
df = pd.read_csv('dados_empresas_rating.csv')

# 2. Definir variáveis de interesse
# X: covariáveis (características das empresas)
X_cols = ['tamanho', 'alavancagem', 'liquidez', 'rentabilidade', 'crescimento', 
          'idade', 'governanca', 'exportadora', 'setor_1', 'setor_2', 'setor_3']
X = df[X_cols]

# T: tratamento (1 se a empresa sofreu rebaixamento, 0 caso contrário)
T = df['rebaixamento'].values.reshape(-1, 1)

# Y: resultado (variação percentual no valor de mercado)
Y = df['var_valor_mercado'].values

# 3. Dividir em treino e teste
X_train, X_test, T_train, T_test, Y_train, Y_test = train_test_split(
    X, T, Y, test_size=0.3, random_state=42)

# 4. Configurar e treinar o modelo Causal Forest
cf_model = CausalForestDML(
    model_y=RandomForestRegressor(n_estimators=100, min_samples_leaf=10),
    model_t=RandomForestClassifier(n_estimators=100, min_samples_leaf=10),
    n_estimators=500,
    min_samples_leaf=10,
    max_depth=10,
    verbose=0,
    random_state=42
)

# Ajustar o modelo
cf_model.fit(Y_train, T_train, X=X_train)

# 5. Estimar efeitos individuais de tratamento
te_test = cf_model.effect(X_test)
df_test = X_test.copy()
df_test['efeito_estimado'] = te_test

# 6. Análise de heterogeneidade por características específicas
# a) Tamanho da empresa
tamanho_bins = pd.qcut(df_test['tamanho'], q=5, labels=False)
efeitos_por_tamanho = df_test.groupby(tamanho_bins)['efeito_estimado'].mean()

plt.figure(figsize=(10, 6))
plt.bar(range(len(efeitos_por_tamanho)), efeitos_por_tamanho.values)
plt.xlabel('Quintil de Tamanho (1=Menor, 5=Maior)')
plt.ylabel('Efeito Médio no Valor de Mercado (%)')
plt.title('Efeito Heterogêneo do Rebaixamento de Rating por Tamanho da Empresa')
plt.axhline(y=0, color='r', linestyle='-')
plt.savefig('heterogeneidade_tamanho.png')

# b) Alavancagem
alavancagem_bins = pd.qcut(df_test['alavancagem'], q=5, labels=False)
efeitos_por_alavancagem = df_test.groupby(alavancagem_bins)['efeito_estimado'].mean()

plt.figure(figsize=(10, 6))
plt.bar(range(len(efeitos_por_alavancagem)), efeitos_por_alavancagem.values)
plt.xlabel('Quintil de Alavancagem (1=Menor, 5=Maior)')
plt.ylabel('Efeito Médio no Valor de Mercado (%)')
plt.title('Efeito Heterogêneo do Rebaixamento de Rating por Alavancagem')
plt.axhline(y=0, color='r', linestyle='-')
plt.savefig('heterogeneidade_alavancagem.png')

# 7. Importância das variáveis para heterogeneidade
importancias = cf_model.feature_importances_
importancias_df = pd.DataFrame({'variavel': X_cols, 'importancia': importancias})
importancias_df = importancias_df.sort_values('importancia', ascending=False)

plt.figure(figsize=(10, 6))
plt.barh(importancias_df['variavel'], importancias_df['importancia'])
plt.xlabel('Importância Relativa')
plt.title('Importância das Variáveis para Heterogeneidade de Efeitos')
plt.tight_layout()
plt.savefig('importancia_variaveis.png')

# 8. Efeitos parciais para variáveis importantes
fig, ax = plt.subplots(figsize=(12, 8))
result = partial_dependence(cf_model.model_cate, 
                           X_test, 
                           features=['alavancagem', 'liquidez'], 
                           kind='average')
XX, YY = np.meshgrid(result["values"][0], result["values"][1])
Z = result["average"].reshape(XX.shape).T
contour = ax.contourf(XX, YY, Z, cmap='viridis')
cbar = plt.colorbar(contour)
cbar.set_label('Efeito do Rebaixamento (%)')
ax.set_xlabel('Alavancagem')
ax.set_ylabel('Liquidez')
ax.set_title('Efeito Conjunto de Alavancagem e Liquidez')
plt.savefig('efeito_conjunto.png')

# 9. Interpretação dos resultados
print("Efeito médio do rebaixamento:", np.mean(te_test))
print("Desvio padrão dos efeitos:", np.std(te_test))
print("Mínimo efeito estimado:", np.min(te_test))
print("Máximo efeito estimado:", np.max(te_test))

# 10. Função para prever o efeito para uma nova empresa
def prever_efeito_rebaixamento(empresa_caracteristicas):
    # empresa_caracteristicas deve ser um DataFrame com as mesmas colunas que X
    return cf_model.effect(empresa_caracteristicas)