Arquivo de Python - Página 52 de 94

Análise Discriminante: Encolhimento (Shrinkage)

19/12/202518/10/2025 Por antonino

Continuando nossa análise detalhada da Análise Discriminante, chegamos a um conceito crucial para aplicações práticas: o encolhimento (shrinkage). Conforme discutimos anteriormente nas formulações matemáticas do LDA e QDA, a estimação precisa das matrizes de covariância é fundamental para o desempenho destes classificadores.

O Problema da Estimação de Covariância

Analogamente ao que observamos nas implementações anteriores, quando trabalhamos com dados de alta dimensionalidade ou tamanhos de amostra limitados, a estimação das matrizes de covariância torna-se desafiadora. Primordialmente, o problema surge porque:

O número de parâmetros a estimar cresce quadraticamente com o número de features
Matrizes de covariância podem tornar-se singular ou mal-condicionada
Estimativas de máxima verossimilhança podem ter alta variância

O Conceito de Encolhimento

O encolhimento é uma técnica de regularização que combina a matriz de covariância estimada dos dados com uma estrutura mais simples e estável. Similarmente a métodos de regularização que encontramos em regressão Ridge e Lasso, o shrinkage visa melhorar a generalização do modelo.

Formulação Matemática

Para o LDA, a matriz de covariância regularizada é dada por:

\(\Sigma_{shrunk} = (1 – \gamma) \Sigma_{empírico} + \gamma \Sigma_{estruturado}\)

Onde γ é o parâmetro de encolhimento que varia entre 0 e 1.

Tipos de Estruturas de Encolhimento

Encolhimento para a Matriz Identidade

Uma abordagem comum é encolher em direção à matriz identidade:

\(\Sigma_{shrunk} = (1 – \gamma) \Sigma + \gamma \frac{\text{tr}(\Sigma)}{d} I\)

Onde d é a dimensionalidade e tr(Σ) é o traço da matriz de covariância.

Encolhimento para Covariância Diagonal

Outra abordagem utiliza uma matriz diagonal como alvo:

\(\Sigma_{shrunk} = (1 – \gamma) \Sigma + \gamma \text{diag}(\Sigma)\)

Vantagens do Encolhimento

Conforme demonstramos nas implementações anteriores, o uso de shrinkage oferece benefícios significativos:

Melhora a condição numérica das matrizes de covariância
Reduz a variância das estimativas
Permite aplicação em dados de alta dimensionalidade
Melhora a performance de generalização

Seleção do Parâmetro de Encolhimento

O parâmetro γ controla o trade-off entre viés e variância:

γ = 0: Usa apenas a covariância empírica (alto risco de overfitting)
γ = 1: Usa apenas a estrutura alvo (alto viés, baixa variância)
Valores intermediários: Balanceiam viés e variância

Implementação no scikit-learn

No scikit-learn, o encolhimento está disponível através do parâmetro shrinkage nas classes LinearDiscriminantAnalysis e QuadraticDiscriminantAnalysis. Ademais, o parâmetro shrinkage_factor pode ser estimado automaticamente usando validação cruzada.

Exemplo Prático em Python

Para ilustrar os efeitos do encolhimento, implementemos um exemplo comparativo detalhado:

import numpy as np
import matplotlib.pyplot as plt
from sklearn.discriminant_analysis import LinearDiscriminantAnalysis, QuadraticDiscriminantAnalysis
from sklearn.datasets import make_classification
from sklearn.model_selection import train_test_split, cross_val_score, GridSearchCV
from sklearn.metrics import accuracy_score, classification_report
from sklearn.preprocessing import StandardScaler
import warnings
warnings.filterwarnings('ignore')
 
'''
Demonstração dos efeitos do encolhimento (shrinkage) 
em classificadores LDA e QDA
'''
 
print("=== EFEITOS DO ENCOLHIMENTO EM LDA E QDA ===")
 
# Criando dataset com alta dimensionalidade e amostras limitadas
print("\n1. CONFIGURAÇÃO DO EXPERIMENTO")
X, y = make_classification(
    n_samples=100,  # Poucas amostras para destacar o efeito do shrinkage
    n_features=20,  # Alta dimensionalidade
    n_informative=5,
    n_redundant=10,
    n_classes=3,
    random_state=42
)
 
print(f"Dimensões do dataset: {X.shape}")
print(f"Número de classes: {len(np.unique(y))}")
 
# Padronizando os dados
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X)
 
# Dividindo em treino e teste
X_train, X_test, y_train, y_test = train_test_split(
    X_scaled, y, test_size=0.3, random_state=42
)
 
'''
Comparação de LDA com diferentes valores de shrinkage
'''
print("\n2. LDA COM DIFERENTES VALORES DE SHRINKAGE")
 
shrinkage_values = [None, 'auto', 0.1, 0.5, 0.9]
lda_results = {}
 
for shrinkage in shrinkage_values:
    if shrinkage == 'auto':
        lda = LinearDiscriminantAnalysis(solver='lsqr', shrinkage='auto')
    else:
        lda = LinearDiscriminantAnalysis(solver='lsqr', shrinkage=shrinkage)
    
    # Treinando o modelo
    lda.fit(X_train, y_train)
    
    # Fazendo previsões
    y_pred = lda.predict(X_test)
    accuracy = accuracy_score(y_test, y_pred)
    
    # Validação cruzada para avaliação robusta
    cv_scores = cross_val_score(lda, X_scaled, y, cv=5)
    
    lda_results[shrinkage] = {
        'model': lda,
        'accuracy': accuracy,
        'cv_mean': np.mean(cv_scores),
        'cv_std': np.std(cv_scores)
    }
    
    print(f"Shrinkage {shrinkage}:")
    print(f"  Acurácia teste: {accuracy:.3f}")
    print(f"  CV Score: {np.mean(cv_scores):.3f} ± {np.std(cv_scores):.3f}")
 
'''
Comparação de QDA com regularização
'''
print("\n3. QDA COM REGULARIZAÇÃO")
 
# QDA não tem parâmetro de shrinkage direto, mas usa reg_param
reg_param_values = [0.0, 0.1, 0.5, 0.9]
qda_results = {}
 
for reg_param in reg_param_values:
    qda = QuadraticDiscriminantAnalysis(reg_param=reg_param)
    
    # Treinando o modelo
    qda.fit(X_train, y_train)
    
    # Fazendo previsões
    y_pred = qda.predict(X_test)
    accuracy = accuracy_score(y_test, y_pred)
    
    # Validação cruzada
    cv_scores = cross_val_score(qda, X_scaled, y, cv=5)
    
    qda_results[reg_param] = {
        'model': qda,
        'accuracy': accuracy,
        'cv_mean': np.mean(cv_scores),
        'cv_std': np.std(cv_scores)
    }
    
    print(f"Reg_param {reg_param}:")
    print(f"  Acurácia teste: {accuracy:.3f}")
    print(f"  CV Score: {np.mean(cv_scores):.3f} ± {np.std(cv_scores):.3f}")
 
'''
Busca automática do melhor parâmetro de shrinkage
'''
print("\n4. BUSCA AUTOMÁTICA DO MELHOR SHRINKAGE")
 
# Usando GridSearchCV para encontrar o melhor shrinkage
param_grid = {
    'shrinkage': [None, 'auto', 0.0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0]
}
 
lda_gs = LinearDiscriminantAnalysis(solver='lsqr')
grid_search = GridSearchCV(lda_gs, param_grid, cv=5, scoring='accuracy')
grid_search.fit(X_scaled, y)
 
print(f"Melhor parâmetro de shrinkage: {grid_search.best_params_}")
print(f"Melhor score de validação cruzada: {grid_search.best_score_:.3f}")
 
'''
Visualização dos resultados - VERSÃO CORRIGIDA
'''
print("\n5. VISUALIZAÇÃO DOS RESULTADOS")

# Preparando dados para plotagem - CORREÇÃO: garantir tamanhos consistentes
shrinkage_labels = [str(s) for s in shrinkage_values]
lda_accuracies = [lda_results[s]['accuracy'] for s in shrinkage_values]
lda_cv_means = [lda_results[s]['cv_mean'] for s in shrinkage_values]
lda_cv_stds = [lda_results[s]['cv_std'] for s in shrinkage_values]

qda_reg_labels = [str(r) for r in reg_param_values]
qda_accuracies = [qda_results[r]['accuracy'] for r in reg_param_values]
qda_cv_means = [qda_results[r]['cv_mean'] for r in reg_param_values]
qda_cv_stds = [qda_results[r]['cv_std'] for r in reg_param_values]

fig, ((ax1, ax2), (ax3, ax4)) = plt.subplots(2, 2, figsize=(15, 12))

# Plot 1: LDA - Acurácia vs Shrinkage
ax1.errorbar(range(len(shrinkage_labels)), lda_accuracies, 
             yerr=lda_cv_stds, fmt='o-', capsize=5, label='Teste ± CV std')
ax1.set_xticks(range(len(shrinkage_labels)))
ax1.set_xticklabels(shrinkage_labels)
ax1.set_xlabel('Valor de Shrinkage')
ax1.set_ylabel('Acurácia')
ax1.set_title('LDA: Acurácia vs Shrinkage')
ax1.legend()
ax1.grid(True, alpha=0.3)

# Plot 2: QDA - Acurácia vs Regularização
ax2.errorbar(range(len(qda_reg_labels)), qda_accuracies, 
             yerr=qda_cv_stds, fmt='o-', capsize=5, label='Teste ± CV std', color='orange')
ax2.set_xticks(range(len(qda_reg_labels)))
ax2.set_xticklabels(qda_reg_labels)
ax2.set_xlabel('Valor de reg_param')
ax2.set_ylabel('Acurácia')
ax2.set_title('QDA: Acurácia vs Regularização')
ax2.legend()
ax2.grid(True, alpha=0.3)

# Plot 3: Comparação LDA vs QDA no melhor parâmetro
best_lda_shrinkage = grid_search.best_params_['shrinkage']
best_lda_score = grid_search.best_score_

best_qda_reg = max(qda_results.items(), key=lambda x: x[1]['cv_mean'])[0]
best_qda_score = qda_results[best_qda_reg]['cv_mean']

models_comparison = ['LDA (melhor)', 'QDA (melhor)']
scores_comparison = [best_lda_score, best_qda_score]

bars = ax3.bar(models_comparison, scores_comparison, color=['blue', 'orange'])
ax3.set_ylabel('Acurácia (Validação Cruzada)')
ax3.set_title('Comparação: Melhor LDA vs Melhor QDA')
ax3.grid(True, alpha=0.3)

# Adicionando valores nas barras
for bar, score in zip(bars, scores_comparison):
    ax3.text(bar.get_x() + bar.get_width()/2, bar.get_height() + 0.01,
             f'{score:.3f}', ha='center', va='bottom')

# Plot 4: Análise da estabilidade dos modelos - CORREÇÃO: garantir tamanhos iguais
shrinkage_stability = [lda_results[s]['cv_std'] for s in shrinkage_values]

# CORREÇÃO: Usar apenas os valores de reg_param que temos para QDA
reg_stability = [qda_results[r]['cv_std'] for r in reg_param_values]

# CORREÇÃO: Criar arrays com o mesmo tamanho para comparação
min_length = min(len(shrinkage_labels), len(qda_reg_labels))
x_pos = np.arange(min_length)
width = 0.35

ax4.bar(x_pos - width/2, shrinkage_stability[:min_length], width, label='LDA', alpha=0.7)
ax4.bar(x_pos + width/2, reg_stability[:min_length], width, label='QDA', alpha=0.7)
ax4.set_xticks(x_pos)
ax4.set_xticklabels(shrinkage_labels[:min_length])
ax4.set_xlabel('Parâmetro de Regularização')
ax4.set_ylabel('Desvio Padrão (CV)')
ax4.set_title('Estabilidade: Desvio Padrão na Validação Cruzada')
ax4.legend()
ax4.grid(True, alpha=0.3)

plt.tight_layout()
plt.show()

'''
Análise de casos extremos - VERSÃO CORRIGIDA
'''
print("\n6. ANÁLISE DE CASOS EXTREMOS")

# Caso com dimensionalidade ainda mais alta
print("Caso com alta dimensionalidade:")
X_high, y_high = make_classification(
    n_samples=50,   # Muito poucas amostras
    n_features=30,  # Dimensionalidade muito alta
    n_informative=5,
    n_classes=2,
    random_state=42
)

X_high_scaled = StandardScaler().fit_transform(X_high)

# Testando LDA sem e com shrinkage - CORREÇÃO: usar solver apropriado
try:
    lda_no_shrinkage = LinearDiscriminantAnalysis(solver='lsqr', shrinkage=None)
    lda_no_shrinkage.fit(X_high_scaled, y_high)
    score_no_shrink = lda_no_shrinkage.score(X_high_scaled, y_high)
    print(f"  LDA sem shrinkage: OK (score: {score_no_shrink:.3f})")
except Exception as e:
    print(f"  LDA sem shrinkage falhou: {e}")

lda_with_shrinkage = LinearDiscriminantAnalysis(solver='lsqr', shrinkage='auto')
lda_with_shrinkage.fit(X_high_scaled, y_high)
score_with_shrink = lda_with_shrinkage.score(X_high_scaled, y_high)
print(f"  LDA com shrinkage automático: OK (score: {score_with_shrink:.3f})")

'''
Interpretação prática dos resultados - VERSÃO ATUALIZADA
'''
print("\n7. INTERPRETAÇÃO PRÁTICA")
print("Os resultados demonstram que:")
print(f"  - Melhor shrinkage para LDA: {best_lda_shrinkage}")
print(f"  - Melhor reg_param para QDA: {best_qda_reg}")
print(f"  - Melhor score LDA: {best_lda_score:.3f}")
print(f"  - Melhor score QDA: {best_qda_score:.3f}")
print("  - O shrinkage melhora a estabilidade numérica do LDA")
print("  - Valores moderados de shrinkage geralmente oferecem o melhor trade-off")

print("\nRecomendações práticas:")
print("  - Use shrinkage quando tiver muitas features ou poucas amostras")
print("  - Experimente com 'auto' primeiro para uma configuração rápida")
print("  - Use validação cruzada para tuning fino do parâmetro")
print("  - Considere QDA com reg_param quando as covariâncias forem diferentes")

100

101

102

103

104

105

106

107

108

109

110

111

112

113

114

115

116

117

118

119

120

121

122

123

124

125

126

127

128

129

130

131

132

133

134

135

136

137

138

139

140

141

142

143

144

145

146

147

148

149

150

151

152

153

154

155

156

157

158

159

160

161

162

163

164

165

166

167

168

169

170

171

172

173

174

175

176

177

178

179

180

181

182

183

184

185

186

187

188

189

190

191

192

193

194

195

196

197

198

199

200

201

202

203

204

205

206

207

208

209

210

211

212

213

214

215

216

217

218

219

220

221

222

223

224

225

226

227

228

229

230

231

232

233

234

235

236

237

238

239

240

241

242

243

244

245

246

247

248

249

250

251

252

253

254

255

256

257

import numpy as np

import matplotlib.pyplot as plt

from sklearn.discriminant_analysis import LinearDiscriminantAnalysis, QuadraticDiscriminantAnalysis

from sklearn.datasets import make_classification

from sklearn.model_selection import train_test_split, cross_val_score, GridSearchCV

from sklearn.metrics import accuracy_score, classification_report

from sklearn.preprocessing import StandardScaler

import warnings

warnings.filterwarnings('ignore')

'''

Demonstração dos efeitos do encolhimento (shrinkage)

em classificadores LDA e QDA

'''

print("=== EFEITOS DO ENCOLHIMENTO EM LDA E QDA ===")

# Criando dataset com alta dimensionalidade e amostras limitadas

print("\n1. CONFIGURAÇÃO DO EXPERIMENTO")

X, y = make_classification(

n_samples=100, # Poucas amostras para destacar o efeito do shrinkage

n_features=20, # Alta dimensionalidade

n_informative=5,

n_redundant=10,

n_classes=3,

random_state=42

)

print(f"Dimensões do dataset: {X.shape}")

print(f"Número de classes: {len(np.unique(y))}")

# Padronizando os dados

scaler = StandardScaler()

X_scaled = scaler.fit_transform(X)

# Dividindo em treino e teste

X_train, X_test, y_train, y_test = train_test_split(

X_scaled, y, test_size=0.3, random_state=42

)

'''

Comparação de LDA com diferentes valores de shrinkage

'''

print("\n2. LDA COM DIFERENTES VALORES DE SHRINKAGE")

shrinkage_values = [None, 'auto', 0.1, 0.5, 0.9]

lda_results = {}

for shrinkage in shrinkage_values:

if shrinkage == 'auto':

lda = LinearDiscriminantAnalysis(solver='lsqr', shrinkage='auto')

else:

lda = LinearDiscriminantAnalysis(solver='lsqr', shrinkage=shrinkage)

# Treinando o modelo

lda.fit(X_train, y_train)

# Fazendo previsões

y_pred = lda.predict(X_test)

accuracy = accuracy_score(y_test, y_pred)

# Validação cruzada para avaliação robusta

cv_scores = cross_val_score(lda, X_scaled, y, cv=5)

lda_results[shrinkage] = {

'model': lda,

'accuracy': accuracy,

'cv_mean': np.mean(cv_scores),

'cv_std': np.std(cv_scores)

}

print(f"Shrinkage {shrinkage}:")

print(f" Acurácia teste: {accuracy:.3f}")

print(f" CV Score: {np.mean(cv_scores):.3f} ± {np.std(cv_scores):.3f}")

'''

Comparação de QDA com regularização

'''

print("\n3. QDA COM REGULARIZAÇÃO")

# QDA não tem parâmetro de shrinkage direto, mas usa reg_param

reg_param_values = [0.0, 0.1, 0.5, 0.9]

qda_results = {}

for reg_param in reg_param_values:

qda = QuadraticDiscriminantAnalysis(reg_param=reg_param)

# Treinando o modelo

qda.fit(X_train, y_train)

# Fazendo previsões

y_pred = qda.predict(X_test)

accuracy = accuracy_score(y_test, y_pred)

# Validação cruzada

cv_scores = cross_val_score(qda, X_scaled, y, cv=5)

qda_results[reg_param] = {

'model': qda,

'accuracy': accuracy,

'cv_mean': np.mean(cv_scores),

'cv_std': np.std(cv_scores)

}

print(f"Reg_param {reg_param}:")

print(f" Acurácia teste: {accuracy:.3f}")

print(f" CV Score: {np.mean(cv_scores):.3f} ± {np.std(cv_scores):.3f}")

'''

Busca automática do melhor parâmetro de shrinkage

'''

print("\n4. BUSCA AUTOMÁTICA DO MELHOR SHRINKAGE")

# Usando GridSearchCV para encontrar o melhor shrinkage

param_grid = {

'shrinkage': [None, 'auto', 0.0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0]

}

lda_gs = LinearDiscriminantAnalysis(solver='lsqr')

grid_search = GridSearchCV(lda_gs, param_grid, cv=5, scoring='accuracy')

grid_search.fit(X_scaled, y)

print(f"Melhor parâmetro de shrinkage: {grid_search.best_params_}")

print(f"Melhor score de validação cruzada: {grid_search.best_score_:.3f}")

'''

Visualização dos resultados - VERSÃO CORRIGIDA

'''

print("\n5. VISUALIZAÇÃO DOS RESULTADOS")

# Preparando dados para plotagem - CORREÇÃO: garantir tamanhos consistentes

shrinkage_labels = [str(s) for s in shrinkage_values]

lda_accuracies = [lda_results[s]['accuracy'] for s in shrinkage_values]

lda_cv_means = [lda_results[s]['cv_mean'] for s in shrinkage_values]

lda_cv_stds = [lda_results[s]['cv_std'] for s in shrinkage_values]

qda_reg_labels = [str(r) for r in reg_param_values]

qda_accuracies = [qda_results[r]['accuracy'] for r in reg_param_values]

qda_cv_means = [qda_results[r]['cv_mean'] for r in reg_param_values]

qda_cv_stds = [qda_results[r]['cv_std'] for r in reg_param_values]

fig, ((ax1, ax2), (ax3, ax4)) = plt.subplots(2, 2, figsize=(15, 12))

# Plot 1: LDA - Acurácia vs Shrinkage

ax1.errorbar(range(len(shrinkage_labels)), lda_accuracies,

yerr=lda_cv_stds, fmt='o-', capsize=5, label='Teste ± CV std')

ax1.set_xticks(range(len(shrinkage_labels)))

ax1.set_xticklabels(shrinkage_labels)

ax1.set_xlabel('Valor de Shrinkage')

ax1.set_ylabel('Acurácia')

ax1.set_title('LDA: Acurácia vs Shrinkage')

ax1.legend()

ax1.grid(True, alpha=0.3)

# Plot 2: QDA - Acurácia vs Regularização

ax2.errorbar(range(len(qda_reg_labels)), qda_accuracies,

yerr=qda_cv_stds, fmt='o-', capsize=5, label='Teste ± CV std', color='orange')

ax2.set_xticks(range(len(qda_reg_labels)))

ax2.set_xticklabels(qda_reg_labels)

ax2.set_xlabel('Valor de reg_param')

ax2.set_ylabel('Acurácia')

ax2.set_title('QDA: Acurácia vs Regularização')

ax2.legend()

ax2.grid(True, alpha=0.3)

# Plot 3: Comparação LDA vs QDA no melhor parâmetro

best_lda_shrinkage = grid_search.best_params_['shrinkage']

best_lda_score = grid_search.best_score_

best_qda_reg = max(qda_results.items(), key=lambda x: x[1]['cv_mean'])[0]

best_qda_score = qda_results[best_qda_reg]['cv_mean']

models_comparison = ['LDA (melhor)', 'QDA (melhor)']

scores_comparison = [best_lda_score, best_qda_score]

bars = ax3.bar(models_comparison, scores_comparison, color=['blue', 'orange'])

ax3.set_ylabel('Acurácia (Validação Cruzada)')

ax3.set_title('Comparação: Melhor LDA vs Melhor QDA')

ax3.grid(True, alpha=0.3)

# Adicionando valores nas barras

for bar, score in zip(bars, scores_comparison):

ax3.text(bar.get_x() + bar.get_width()/2, bar.get_height() + 0.01,

f'{score:.3f}', ha='center', va='bottom')

# Plot 4: Análise da estabilidade dos modelos - CORREÇÃO: garantir tamanhos iguais

shrinkage_stability = [lda_results[s]['cv_std'] for s in shrinkage_values]

# CORREÇÃO: Usar apenas os valores de reg_param que temos para QDA

reg_stability = [qda_results[r]['cv_std'] for r in reg_param_values]

# CORREÇÃO: Criar arrays com o mesmo tamanho para comparação

min_length = min(len(shrinkage_labels), len(qda_reg_labels))

x_pos = np.arange(min_length)

width = 0.35

ax4.bar(x_pos - width/2, shrinkage_stability[:min_length], width, label='LDA', alpha=0.7)

ax4.bar(x_pos + width/2, reg_stability[:min_length], width, label='QDA', alpha=0.7)

ax4.set_xticks(x_pos)

ax4.set_xticklabels(shrinkage_labels[:min_length])

ax4.set_xlabel('Parâmetro de Regularização')

ax4.set_ylabel('Desvio Padrão (CV)')

ax4.set_title('Estabilidade: Desvio Padrão na Validação Cruzada')

ax4.legend()

ax4.grid(True, alpha=0.3)

plt.tight_layout()

plt.show()

'''

Análise de casos extremos - VERSÃO CORRIGIDA

'''

print("\n6. ANÁLISE DE CASOS EXTREMOS")

# Caso com dimensionalidade ainda mais alta

print("Caso com alta dimensionalidade:")

X_high, y_high = make_classification(

n_samples=50, # Muito poucas amostras

n_features=30, # Dimensionalidade muito alta

n_informative=5,

n_classes=2,

random_state=42

)

X_high_scaled = StandardScaler().fit_transform(X_high)

# Testando LDA sem e com shrinkage - CORREÇÃO: usar solver apropriado

try:

lda_no_shrinkage = LinearDiscriminantAnalysis(solver='lsqr', shrinkage=None)

lda_no_shrinkage.fit(X_high_scaled, y_high)

score_no_shrink = lda_no_shrinkage.score(X_high_scaled, y_high)

print(f" LDA sem shrinkage: OK (score: {score_no_shrink:.3f})")

except Exception as e:

print(f" LDA sem shrinkage falhou: {e}")

lda_with_shrinkage = LinearDiscriminantAnalysis(solver='lsqr', shrinkage='auto')

lda_with_shrinkage.fit(X_high_scaled, y_high)

score_with_shrink = lda_with_shrinkage.score(X_high_scaled, y_high)

print(f" LDA com shrinkage automático: OK (score: {score_with_shrink:.3f})")

'''

Interpretação prática dos resultados - VERSÃO ATUALIZADA

'''

print("\n7. INTERPRETAÇÃO PRÁTICA")

print("Os resultados demonstram que:")

print(f" - Melhor shrinkage para LDA: {best_lda_shrinkage}")

print(f" - Melhor reg_param para QDA: {best_qda_reg}")

print(f" - Melhor score LDA: {best_lda_score:.3f}")

print(f" - Melhor score QDA: {best_qda_score:.3f}")

print(" - O shrinkage melhora a estabilidade numérica do LDA")

print(" - Valores moderados de shrinkage geralmente oferecem o melhor trade-off")

print("\nRecomendações práticas:")

print(" - Use shrinkage quando tiver muitas features ou poucas amostras")

print(" - Experimente com 'auto' primeiro para uma configuração rápida")

print(" - Use validação cruzada para tuning fino do parâmetro")

print(" - Considere QDA com reg_param quando as covariâncias forem diferentes")

Interpretação dos Resultados

Analisando os experimentos, podemos observar padrões importantes:

O encolhimento melhora significativamente a estabilidade dos classificadores
Valores intermediários de shrinkage geralmente produzem os melhores resultados
O parâmetro ‘auto’ oferece uma solução prática e eficiente
Em casos de alta dimensionalidade, o shrinkage torna-se essencial

Conexões com Tópicos Anteriores

Similarmente aos conceitos de regularização que discutimos em regressão linear, o encolhimento no LDA e QDA:

Balanceia o trade-off entre viés e variância
Melhora a generalização em dados limitados
Lida com problemas de multicolinearidade
Oferece robustez contra overfitting

Considerações Práticas

Quando Usar Encolhimento

Número de features maior que número de amostras
Matrizes de covariância singular ou mal-condicionada
Dados com alta correlação entre features
Aplicações onde robustez é mais importante que precisão máxima

Escolha do Parâmetro

Inegavelmente, a seleção do parâmetro de encolhimento deve considerar:

Usar validação cruzada para tuning
Começar com valores padrão ou ‘auto’
Considerar o trade-off entre performance e complexidade computacional
Avaliar a sensibilidade do modelo a diferentes valores

Conclusão

O encolhimento representa uma ferramenta essencial para aplicar LDA e QDA em cenários práticos, especialmente quando lidamos com as limitações de dados do mundo real. Embora adicione complexidade ao modelo, os benefícios em termos de robustez e generalização frequentemente justificam seu uso.

Portanto, ao implementar classificadores discriminantes em aplicações reais, o uso criterioso de técnicas de encolhimento pode significar a diferença entre um modelo que funciona bem na prática e um que falha devido a problemas numéricos ou overfitting.

Referência

Este post explora o item 1.2.4. Encolhimento da documentação do scikit-learn:

https://scikit-learn.org/0.21/modules/lda_qda.html

Análise Discriminante: Formulação Matemática da Redução de Dimensionalidade LDA

19/12/202518/10/2025 Por antonino

Após explorarmos as formulações dos classificadores LDA e QDA, aprofundemo-nos agora na fundamentação matemática da redução de dimensionalidade via LDA. Primordialmente, esta abordagem representa uma aplicação elegante dos conceitos discriminativos para transformação de dados em espaços de menor dimensão.

O Problema de Otimização

Conforme discutimos anteriormente, o LDA para redução de dimensionalidade busca encontrar projeções que maximizem a separação entre classes. Analogamente, podemos formular este problema como uma otimização que maximiza o critério de Fisher:

\(J(W) = \frac{W^T S_B W}{W^T S_W W}\)

Onde W é a matriz de transformação que projeta os dados em um espaço de dimensão reduzida.

Matrizes de Dispersão

Matriz de Dispersão Entre Classes (Between-Class Scatter)

Esta matriz captura a dispersão das médias das classes em relação à média global:

\(S_B = \sum_{k=1}^K N_k (\mu_k – \mu)(\mu_k – \mu)^T\)

Onde:

\(K\) é o número de classes
\(N_k\) é o número de amostras na classe k
\(μ_k\) é a média da classe k
\(μ\) é a média global dos dados

Matriz de Dispersão Dentro das Classes (Within-Class Scatter)

Esta matriz quantifica a dispersão dos dados dentro de cada classe:

\(S_W = \sum_{k=1}^K \sum_{x_i \in C_k} (x_i – \mu_k)(x_i – \mu_k)^T\)

Alternativamente, podemos expressar S_W em termos das matrizes de covariância de cada classe:

\(S_W = \sum_{k=1}^K N_k \Sigma_k\)

Solução do Problema de Autovalor Generalizado

Maximizar J(W) equivale a resolver o problema de autovalor generalizado:

\(S_B w = \lambda S_W w\)

Os autovetores correspondentes aos maiores autovalores formam as direções de projeção ótimas. Ademais, o número máximo de componentes discriminantes é limitado por:

\(L \leq \min(K-1, d)\)

Onde K é o número de classes e d é a dimensionalidade original.

Interpretação Geométrica

Projeção Ótima

Os autovetores w_i que resolvem o problema de autovalor definem um novo sistema de coordenadas onde:

\(y = W^T x\)

Neste espaço projetado, a razão entre a dispersão entre classes e dentro das classes é maximizada.

Variância Explicada

Os autovalores λ_i associados a cada autovetor indicam o poder discriminativo de cada componente:

\(\text{Variância explicada} = \frac{\lambda_i}{\sum_{j=1}^L \lambda_j}\)

Relação com Análise de Variância (ANOVA)

Similarmente à ANOVA, o LDA decompõe a variabilidade total em componentes:

\(S_T = S_W + S_B\)

Onde \(S_T\) é a matriz de dispersão total. Inegavelmente, esta decomposição revela a estrutura fundamental dos dados.

Regularização e Estabilidade Numérica

Quando \(S_W\) é singular ou mal-condicionada, utilizamos regularização:

\(S_W^{reg} = S_W + \gamma I\)

Esta abordagem melhora a estabilidade numérica sem comprometer significativamente o poder discriminativo.

Exemplo Prático em Python

Para ilustrar estas formulações matemáticas, implementemos um exemplo detalhado que demonstra os cálculos fundamentais:

import numpy as np
import matplotlib.pyplot as plt
from sklearn.discriminant_analysis import LinearDiscriminantAnalysis
from sklearn.datasets import load_iris, make_classification
from sklearn.preprocessing import StandardScaler
import scipy.linalg as la

'''
Implementação manual da redução de dimensionalidade LDA
para demonstrar as formulações matemáticas
'''

class ManualLDAReducer:
    '''Implementação manual do LDA para redução de dimensionalidade'''
    
    def fit(self, X, y, n_components=None):
        """
        Implementa a formulação matemática completa do LDA
        para redução de dimensionalidade
        """
        self.classes_ = np.unique(y)
        self.n_classes = len(self.classes_)
        self.n_features = X.shape[1]
        
        if n_components is None:
            n_components = min(self.n_classes - 1, self.n_features)
        self.n_components = n_components
        
        # Calculando estatísticas básicas
        self.means_ = self._calculate_class_means(X, y)
        self.global_mean_ = np.mean(X, axis=0)
        self.priors_ = self._calculate_priors(y)
        
        # Calculando matrizes de dispersão
        self.S_B = self._calculate_between_class_scatter()
        self.S_W = self._calculate_within_class_scatter(X, y)
        
        # Resolvendo problema de autovalor generalizado
        self.eigenvalues_, self.eigenvectors_ = self._solve_eigenproblem()
        
        # Selecionando componentes principais
        self.components_ = self.eigenvectors_[:, :self.n_components]
        
        return self
    
    def _calculate_class_means(self, X, y):
        """Calcula as médias de cada classe"""
        means = []
        for k in self.classes_:
            class_mask = (y == k)
            means.append(np.mean(X[class_mask], axis=0))
        return np.array(means)
    
    def _calculate_priors(self, y):
        """Calcula as probabilidades a priori de cada classe"""
        priors = []
        for k in self.classes_:
            priors.append(np.mean(y == k))
        return np.array(priors)
    
    def _calculate_between_class_scatter(self):
        """Calcula a matriz de dispersão entre classes S_B"""
        S_B = np.zeros((self.n_features, self.n_features))
        for k, mean_k in enumerate(self.means_):
            n_k = np.sum(self.priors_[k] * len(self.classes_))  # Aproximação
            diff = (mean_k - self.global_mean_).reshape(-1, 1)
            S_B += n_k * (diff @ diff.T)
        return S_B
    
    def _calculate_within_class_scatter(self, X, y):
        """Calcula a matriz de dispersão dentro das classes S_W"""
        S_W = np.zeros((self.n_features, self.n_features))
        for k, mean_k in enumerate(self.means_):
            class_mask = (y == self.classes_[k])
            X_class = X[class_mask]
            X_centered = X_class - mean_k
            S_W += X_centered.T @ X_centered
        return S_W
    
    def _solve_eigenproblem(self):
        """
        Resolve o problema de autovalor generalizado S_B w = λ S_W w
        usando decomposição simultânea das matrizes
        """
        # Decomposição de Cholesky de S_W para estabilidade numérica
        try:
            L = la.cholesky(self.S_W, lower=True)
            Linv = la.inv(L)
            
            # Transformando o problema em problema de autovalor regular
            S_B_transformed = Linv @ self.S_B @ Linv.T
            
            # Resolvendo problema de autovalor simétrico
            eigvals, eigvecs_transformed = la.eigh(S_B_transformed)
            
            # Transformando autovetores de volta para o espaço original
            eigenvectors = Linv.T @ eigvecs_transformed
            
            # Ordenando autovalores e autovetores em ordem decrescente
            idx = np.argsort(eigvals)[::-1]
            eigenvalues = eigvals[idx]
            eigenvectors = eigenvectors[:, idx]
            
        except la.LinAlgError:
            # Fallback: usar pseudo-inversa se Cholesky falhar
            print("Cholesky failed, using pseudo-inverse approach")
            S_W_pinv = la.pinv(self.S_W)
            S_combined = S_W_pinv @ self.S_B
            eigenvalues, eigenvectors = la.eig(S_combined)
            eigenvalues = eigenvalues.real
            eigenvectors = eigenvectors.real
            
            # Ordenando
            idx = np.argsort(eigenvalues)[::-1]
            eigenvalues = eigenvalues[idx]
            eigenvectors = eigenvectors[:, idx]
        
        return eigenvalues, eigenvectors
    
    def transform(self, X):
        """Projeta os dados no espaço LDA"""
        return X @ self.components_
    
    def explained_variance_ratio(self):
        """Calcula a razão de variância explicada por cada componente"""
        total_variance = np.sum(self.eigenvalues_)
        return self.eigenvalues_[:self.n_components] / total_variance

'''
Demonstração completa da formulação matemática
'''
print("=== FORMULAÇÃO MATEMÁTICA DA REDUÇÃO LDA ===")

# Carregando dataset para demonstração
iris = load_iris()
X, y = iris.data, iris.target

# Padronizando os dados
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X)

print(f"Dados originais: {X_scaled.shape}")
print(f"Número de classes: {len(np.unique(y))}")

# Aplicando implementação manual
print("\n1. IMPLEMENTAÇÃO MANUAL")
manual_lda = ManualLDAReducer()
manual_lda.fit(X_scaled, y, n_components=2)
X_manual = manual_lda.transform(X_scaled)

print(f"Matriz S_B (Between-class scatter): {manual_lda.S_B.shape}")
print(f"Matriz S_W (Within-class scatter): {manual_lda.S_W.shape}")
print(f"Autovalores: {manual_lda.eigenvalues_}")
print(f"Razão de variância explicada: {manual_lda.explained_variance_ratio()}")
print(f"Variância total explicada: {np.sum(manual_lda.explained_variance_ratio()):.3f}")

# Comparando com scikit-learn
print("\n2. COMPARAÇÃO COM SCIKIT-LEARN")
sklearn_lda = LinearDiscriminantAnalysis(n_components=2)
X_sklearn = sklearn_lda.fit_transform(X_scaled, y)

print(f"Autovalores (scikit-learn): {sklearn_lda.explained_variance_ratio_}")
print(f"Variância total explicada: {np.sum(sklearn_lda.explained_variance_ratio_):.3f}")

'''
Análise das matrizes de dispersão
'''
print("\n3. ANÁLISE DAS MATRIZES DE DISPERSÃO")

# Verificando propriedades matemáticas
S_T = manual_lda.S_W + manual_lda.S_B  # Decomposição total
print(f"Decomposição S_T = S_W + S_B:")
print(f"  Norma de S_T - (S_W + S_B): {np.linalg.norm(S_T - (manual_lda.S_W + manual_lda.S_B)):.6f}")

# Análise de rank
print(f"Rank das matrizes:")
print(f"  Rank(S_B): {np.linalg.matrix_rank(manual_lda.S_B)}")
print(f"  Rank(S_W): {np.linalg.matrix_rank(manual_lda.S_W)}")
print(f"  Rank(S_T): {np.linalg.matrix_rank(S_T)}")

'''
Visualização comparativa
'''
print("\n4. VISUALIZAÇÃO DOS RESULTADOS")

fig, ((ax1, ax2), (ax3, ax4)) = plt.subplots(2, 2, figsize=(15, 12))

# Plot 1: Dados originais (duas primeiras features)
for i, class_label in enumerate(np.unique(y)):
    ax1.scatter(X_scaled[y == class_label, 0], X_scaled[y == class_label, 1], 
                alpha=0.7, label=f'Classe {class_label}')
ax1.set_title('Dados Originais (2 Features)')
ax1.set_xlabel('Feature 1')
ax1.set_ylabel('Feature 2')
ax1.legend()
ax1.grid(True, alpha=0.3)

# Plot 2: Projeção manual LDA
for i, class_label in enumerate(np.unique(y)):
    ax2.scatter(X_manual[y == class_label, 0], X_manual[y == class_label, 1], 
                alpha=0.7, label=f'Classe {class_label}')
ax2.set_title('Projeção LDA Manual')
ax2.set_xlabel('Componente LDA 1')
ax2.set_ylabel('Componente LDA 2')
ax2.legend()
ax2.grid(True, alpha=0.3)

# Plot 3: Projeção scikit-learn LDA
for i, class_label in enumerate(np.unique(y)):
    ax3.scatter(X_sklearn[y == class_label, 0], X_sklearn[y == class_label, 1], 
                alpha=0.7, label=f'Classe {class_label}')
ax3.set_title('Projeção LDA Scikit-learn')
ax3.set_xlabel('Componente LDA 1')
ax3.set_ylabel('Componente LDA 2')
ax3.legend()
ax3.grid(True, alpha=0.3)

# Plot 4: Comparação das projeções
ax4.scatter(X_manual[:, 0], X_manual[:, 1], alpha=0.5, label='Manual', marker='o')
ax4.scatter(X_sklearn[:, 0], X_sklearn[:, 1], alpha=0.5, label='Scikit-learn', marker='x')
ax4.set_title('Comparação das Projeções')
ax4.set_xlabel('Componente LDA 1')
ax4.set_ylabel('Componente LDA 2')
ax4.legend()
ax4.grid(True, alpha=0.3)

plt.tight_layout()
plt.show()

'''
Análise dos autovetores e direções discriminantes
'''
print("\n5. ANÁLISE DAS DIREÇÕES DISCRIMINANTES")

print("Autovetores (componentes LDA) - Manual:")
print(manual_lda.components_)

if hasattr(sklearn_lda, 'scalings_'):
    print("\nAutovetores (componentes LDA) - Scikit-learn:")
    print(sklearn_lda.scalings_)

# Calculando ângulo entre os componentes
if hasattr(sklearn_lda, 'scalings_') and manual_lda.components_.shape == sklearn_lda.scalings_.shape:
    for i in range(manual_lda.n_components):
        v1 = manual_lda.components_[:, i]
        v2 = sklearn_lda.scalings_[:, i]
        cosine_sim = np.abs(np.dot(v1, v2) / (np.linalg.norm(v1) * np.linalg.norm(v2)))
        print(f"Similaridade do cosseno componente {i+1}: {cosine_sim:.6f}")

'''
Exemplo com dados de alta dimensionalidade
'''
print("\n6. APLICAÇÃO EM ALTA DIMENSIONALIDADE")

# Criando dataset com mais features
X_high, y_high = make_classification(
    n_samples=200, n_features=10, n_informative=5, 
    n_redundant=3, n_classes=4, random_state=42
)

X_high_scaled = StandardScaler().fit_transform(X_high)

manual_lda_high = ManualLDAReducer()
manual_lda_high.fit(X_high_scaled, y_high, n_components=3)
X_high_manual = manual_lda_high.transform(X_high_scaled)

print(f"Dimensionalidade original: {X_high_scaled.shape}")
print(f"Dimensionalidade após LDA: {X_high_manual.shape}")
print(f"Autovalores: {manual_lda_high.eigenvalues_[:3]}")
print(f"Variância explicada: {manual_lda_high.explained_variance_ratio()}")

'''
Interpretação matemática dos resultados
'''
print("\n7. INTERPRETAÇÃO MATEMÁTICA")
print("A formulação matemática demonstra que:")
print("  - O LDA maximiza a razão entre dispersão entre classes e dentro delas")
print("  - Os autovetores definem as direções de máxima separabilidade")
print("  - Os autovalores quantificam o poder discriminativo de cada componente")
print("  - A decomposição S_T = S_W + S_B revela a estrutura de variabilidade")
print("  - A implementação manual valida as fundamentações teóricas")

100

101

102

103

104

105

106

107

108

109

110

111

112

113

114

115

116

117

118

119

120

121

122

123

124

125

126

127

128

129

130

131

132

133

134

135

136

137

138

139

140

141

142

143

144

145

146

147

148

149

150

151

152

153

154

155

156

157

158

159

160

161

162

163

164

165

166

167

168

169

170

171

172

173

174

175

176

177

178

179

180

181

182

183

184

185

186

187

188

189

190

191

192

193

194

195

196

197

198

199

200

201

202

203

204

205

206

207

208

209

210

211

212

213

214

215

216

217

218

219

220

221

222

223

224

225

226

227

228

229

230

231

232

233

234

235

236

237

238

239

240

241

242

243

244

245

246

247

248

249

250

251

252

253

254

255

256

257

258

259

260

261

262

263

264

265

266

267

268

269

270

271

272

273

274

275

276

277

278

279

280

281

import numpy as np

import matplotlib.pyplot as plt

from sklearn.discriminant_analysis import LinearDiscriminantAnalysis

from sklearn.datasets import load_iris, make_classification

from sklearn.preprocessing import StandardScaler

import scipy.linalg as la

'''

Implementação manual da redução de dimensionalidade LDA

para demonstrar as formulações matemáticas

'''

class ManualLDAReducer:

'''Implementação manual do LDA para redução de dimensionalidade'''

def fit(self, X, y, n_components=None):

"""

Implementa a formulação matemática completa do LDA

para redução de dimensionalidade

"""

self.classes_ = np.unique(y)

self.n_classes = len(self.classes_)

self.n_features = X.shape[1]

if n_components is None:

n_components = min(self.n_classes - 1, self.n_features)

self.n_components = n_components

# Calculando estatísticas básicas

self.means_ = self._calculate_class_means(X, y)

self.global_mean_ = np.mean(X, axis=0)

self.priors_ = self._calculate_priors(y)

# Calculando matrizes de dispersão

self.S_B = self._calculate_between_class_scatter()

self.S_W = self._calculate_within_class_scatter(X, y)

# Resolvendo problema de autovalor generalizado

self.eigenvalues_, self.eigenvectors_ = self._solve_eigenproblem()

# Selecionando componentes principais

self.components_ = self.eigenvectors_[:, :self.n_components]

return self

def _calculate_class_means(self, X, y):

"""Calcula as médias de cada classe"""

means = []

for k in self.classes_:

class_mask = (y == k)

means.append(np.mean(X[class_mask], axis=0))

return np.array(means)

def _calculate_priors(self, y):

"""Calcula as probabilidades a priori de cada classe"""

priors = []

for k in self.classes_:

priors.append(np.mean(y == k))

return np.array(priors)

def _calculate_between_class_scatter(self):

"""Calcula a matriz de dispersão entre classes S_B"""

S_B = np.zeros((self.n_features, self.n_features))

for k, mean_k in enumerate(self.means_):

n_k = np.sum(self.priors_[k] * len(self.classes_)) # Aproximação

diff = (mean_k - self.global_mean_).reshape(-1, 1)

S_B += n_k * (diff @ diff.T)

return S_B

def _calculate_within_class_scatter(self, X, y):

"""Calcula a matriz de dispersão dentro das classes S_W"""

S_W = np.zeros((self.n_features, self.n_features))

for k, mean_k in enumerate(self.means_):

class_mask = (y == self.classes_[k])

X_class = X[class_mask]

X_centered = X_class - mean_k

S_W += X_centered.T @ X_centered

return S_W

def _solve_eigenproblem(self):

"""

Resolve o problema de autovalor generalizado S_B w = λ S_W w

usando decomposição simultânea das matrizes

"""

# Decomposição de Cholesky de S_W para estabilidade numérica

try:

L = la.cholesky(self.S_W, lower=True)

Linv = la.inv(L)

# Transformando o problema em problema de autovalor regular

S_B_transformed = Linv @ self.S_B @ Linv.T

# Resolvendo problema de autovalor simétrico

eigvals, eigvecs_transformed = la.eigh(S_B_transformed)

# Transformando autovetores de volta para o espaço original

eigenvectors = Linv.T @ eigvecs_transformed

# Ordenando autovalores e autovetores em ordem decrescente

idx = np.argsort(eigvals)[::-1]

eigenvalues = eigvals[idx]

eigenvectors = eigenvectors[:, idx]

except la.LinAlgError:

# Fallback: usar pseudo-inversa se Cholesky falhar

print("Cholesky failed, using pseudo-inverse approach")

S_W_pinv = la.pinv(self.S_W)

S_combined = S_W_pinv @ self.S_B

eigenvalues, eigenvectors = la.eig(S_combined)

eigenvalues = eigenvalues.real

eigenvectors = eigenvectors.real

# Ordenando

idx = np.argsort(eigenvalues)[::-1]

eigenvalues = eigenvalues[idx]

eigenvectors = eigenvectors[:, idx]

return eigenvalues, eigenvectors

def transform(self, X):

"""Projeta os dados no espaço LDA"""

return X @ self.components_

def explained_variance_ratio(self):

"""Calcula a razão de variância explicada por cada componente"""

total_variance = np.sum(self.eigenvalues_)

return self.eigenvalues_[:self.n_components] / total_variance

'''

Demonstração completa da formulação matemática

'''

print("=== FORMULAÇÃO MATEMÁTICA DA REDUÇÃO LDA ===")

# Carregando dataset para demonstração

iris = load_iris()

X, y = iris.data, iris.target

# Padronizando os dados

scaler = StandardScaler()

X_scaled = scaler.fit_transform(X)

print(f"Dados originais: {X_scaled.shape}")

print(f"Número de classes: {len(np.unique(y))}")

# Aplicando implementação manual

print("\n1. IMPLEMENTAÇÃO MANUAL")

manual_lda = ManualLDAReducer()

manual_lda.fit(X_scaled, y, n_components=2)

X_manual = manual_lda.transform(X_scaled)

print(f"Matriz S_B (Between-class scatter): {manual_lda.S_B.shape}")

print(f"Matriz S_W (Within-class scatter): {manual_lda.S_W.shape}")

print(f"Autovalores: {manual_lda.eigenvalues_}")

print(f"Razão de variância explicada: {manual_lda.explained_variance_ratio()}")

print(f"Variância total explicada: {np.sum(manual_lda.explained_variance_ratio()):.3f}")

# Comparando com scikit-learn

print("\n2. COMPARAÇÃO COM SCIKIT-LEARN")

sklearn_lda = LinearDiscriminantAnalysis(n_components=2)

X_sklearn = sklearn_lda.fit_transform(X_scaled, y)

print(f"Autovalores (scikit-learn): {sklearn_lda.explained_variance_ratio_}")

print(f"Variância total explicada: {np.sum(sklearn_lda.explained_variance_ratio_):.3f}")

'''

Análise das matrizes de dispersão

'''

print("\n3. ANÁLISE DAS MATRIZES DE DISPERSÃO")

# Verificando propriedades matemáticas

S_T = manual_lda.S_W + manual_lda.S_B # Decomposição total

print(f"Decomposição S_T = S_W + S_B:")

print(f" Norma de S_T - (S_W + S_B): {np.linalg.norm(S_T - (manual_lda.S_W + manual_lda.S_B)):.6f}")

# Análise de rank

print(f"Rank das matrizes:")

print(f" Rank(S_B): {np.linalg.matrix_rank(manual_lda.S_B)}")

print(f" Rank(S_W): {np.linalg.matrix_rank(manual_lda.S_W)}")

print(f" Rank(S_T): {np.linalg.matrix_rank(S_T)}")

'''

Visualização comparativa

'''

print("\n4. VISUALIZAÇÃO DOS RESULTADOS")

fig, ((ax1, ax2), (ax3, ax4)) = plt.subplots(2, 2, figsize=(15, 12))

# Plot 1: Dados originais (duas primeiras features)

for i, class_label in enumerate(np.unique(y)):

ax1.scatter(X_scaled[y == class_label, 0], X_scaled[y == class_label, 1],

alpha=0.7, label=f'Classe {class_label}')

ax1.set_title('Dados Originais (2 Features)')

ax1.set_xlabel('Feature 1')

ax1.set_ylabel('Feature 2')

ax1.legend()

ax1.grid(True, alpha=0.3)

# Plot 2: Projeção manual LDA

for i, class_label in enumerate(np.unique(y)):

ax2.scatter(X_manual[y == class_label, 0], X_manual[y == class_label, 1],

alpha=0.7, label=f'Classe {class_label}')

ax2.set_title('Projeção LDA Manual')

ax2.set_xlabel('Componente LDA 1')

ax2.set_ylabel('Componente LDA 2')

ax2.legend()

ax2.grid(True, alpha=0.3)

# Plot 3: Projeção scikit-learn LDA

for i, class_label in enumerate(np.unique(y)):

ax3.scatter(X_sklearn[y == class_label, 0], X_sklearn[y == class_label, 1],

alpha=0.7, label=f'Classe {class_label}')

ax3.set_title('Projeção LDA Scikit-learn')

ax3.set_xlabel('Componente LDA 1')

ax3.set_ylabel('Componente LDA 2')

ax3.legend()

ax3.grid(True, alpha=0.3)

# Plot 4: Comparação das projeções

ax4.scatter(X_manual[:, 0], X_manual[:, 1], alpha=0.5, label='Manual', marker='o')

ax4.scatter(X_sklearn[:, 0], X_sklearn[:, 1], alpha=0.5, label='Scikit-learn', marker='x')

ax4.set_title('Comparação das Projeções')

ax4.set_xlabel('Componente LDA 1')

ax4.set_ylabel('Componente LDA 2')

ax4.legend()

ax4.grid(True, alpha=0.3)

plt.tight_layout()

plt.show()

'''

Análise dos autovetores e direções discriminantes

'''

print("\n5. ANÁLISE DAS DIREÇÕES DISCRIMINANTES")

print("Autovetores (componentes LDA) - Manual:")

print(manual_lda.components_)

if hasattr(sklearn_lda, 'scalings_'):

print("\nAutovetores (componentes LDA) - Scikit-learn:")

print(sklearn_lda.scalings_)

# Calculando ângulo entre os componentes

if hasattr(sklearn_lda, 'scalings_') and manual_lda.components_.shape == sklearn_lda.scalings_.shape:

for i in range(manual_lda.n_components):

v1 = manual_lda.components_[:, i]

v2 = sklearn_lda.scalings_[:, i]

cosine_sim = np.abs(np.dot(v1, v2) / (np.linalg.norm(v1) * np.linalg.norm(v2)))

print(f"Similaridade do cosseno componente {i+1}: {cosine_sim:.6f}")

'''

Exemplo com dados de alta dimensionalidade

'''

print("\n6. APLICAÇÃO EM ALTA DIMENSIONALIDADE")

# Criando dataset com mais features

X_high, y_high = make_classification(

n_samples=200, n_features=10, n_informative=5,

n_redundant=3, n_classes=4, random_state=42

)

X_high_scaled = StandardScaler().fit_transform(X_high)

manual_lda_high = ManualLDAReducer()

manual_lda_high.fit(X_high_scaled, y_high, n_components=3)

X_high_manual = manual_lda_high.transform(X_high_scaled)

print(f"Dimensionalidade original: {X_high_scaled.shape}")

print(f"Dimensionalidade após LDA: {X_high_manual.shape}")

print(f"Autovalores: {manual_lda_high.eigenvalues_[:3]}")

print(f"Variância explicada: {manual_lda_high.explained_variance_ratio()}")

'''

Interpretação matemática dos resultados

'''

print("\n7. INTERPRETAÇÃO MATEMÁTICA")

print("A formulação matemática demonstra que:")

print(" - O LDA maximiza a razão entre dispersão entre classes e dentro delas")

print(" - Os autovetores definem as direções de máxima separabilidade")

print(" - Os autovalores quantificam o poder discriminativo de cada componente")

print(" - A decomposição S_T = S_W + S_B revela a estrutura de variabilidade")

print(" - A implementação manual valida as fundamentações teóricas")

Interpretação dos Resultados Matemáticos

Analisando a implementação e os resultados, podemos observar que:

As matrizes S_B e S_W capturam adequadamente a estrutura de dispersão dos dados
Os autovalores refletem o poder discriminativo de cada componente LDA
A decomposição da variabilidade total é matematicamente consistente
As projeções preservam a separabilidade entre classes mesmo com redução dimensional

Propriedades Matemáticas Importantes

Ortogonalidade dos Componentes

No espaço transformado, os componentes LDA são ortogonais em relação à métrica definida por S_W:

\(w_i^T S_W w_j = 0 \quad \text{para } i \neq j\)

Maximização Sequencial

Os componentes são encontrados sequencialmente, onde cada novo componente maximiza a separação sujeito à ortogonalidade com os anteriores.

Considerações Numéricas

Na prática, várias considerações numéricas são importantes:

O problema de autovalor generalizado requer cuidado com matrizes singulares
Técnicas de regularização melhoram a estabilidade numérica
Decomposições como Cholesky podem ser preferíveis para eficiência
A seleção do número de componentes deve considerar tanto aspectos matemáticos quanto práticos

Relação com Outras Técnicas

Similarmente a outras técnicas de redução de dimensionalidade, o LDA compartilha conceitos fundamentais:

Como o PCA, busca direções de máxima variância, mas com foco na variância entre classes
Como a ANOVA, decompõe a variabilidade total em componentes estruturados
Como métodos de projeção linear, transforma dados preservando propriedades desejadas

Portanto, a formulação matemática do LDA para redução de dimensionalidade não apenas fornece uma ferramenta prática, mas também revela insights profundos sobre a estrutura dos dados e as relações entre observações de diferentes classes.