Objetivo é construir modelo para predizer quando pessoa tem câncer conforme exame de sequenciamento do RNA. Os dados são reais e foram devidamente anonimizados. No dataframe, Cada linha representa amostra retirada de pessoa. Colunas são tipos de microRNA, cada registro representando intensidade com que tal microRNA está expresso. Valores de expressão variam entre [0,n]. Valores próximos a 0 indicam baixa expressão, enquanto que contrário indica alta expressão. Dados também apresentam rótulos respostas (class) sendo TP (primary solid tumor) indicando tumor e NT (normal tissue).
O notebook e dados foram referenciados do The Cancer Genome Atlas Program (TCGA) e Escola Britânica de Artes Criativas & Tecnologia (EBAC).
# Importar dependências:
import pandas as pd
import numpy as np
import seaborn as sns
import matplotlib.pyplot as plt
from sklearn.model_selection import train_test_split
from sklearn.linear_model import LogisticRegression
from sklearn.metrics import balanced_accuracy_score
from sklearn.model_selection import cross_val_score
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import StandardScaler
from sklearn.neighbors import KNeighborsClassifier
from sklearn.preprocessing import MinMaxScaler
from sklearn.feature_selection import SelectKBest
from sklearn.feature_selection import chi2
from sklearn.decomposition import PCA
# Carregar dados:
df = pd.read_csv("brca_mirnaseq.csv", sep=';', header=0, decimal=',')
df.head()
df.info()
df # 842 linhas, 898 colunas
df.shape
# Análise exploratória:
ax = sns.countplot(x="class", data=df)
df["class"].value_counts() # 755 TP (apresenta tumor), 87 NT (não apresenta tumor)
df["class"].value_counts(normalize=True) # 0.8962 TP (89.62%), 0.1038 NT (10.38%)
df.describe() # média, desvio padrão, mínimo, 25%, 50%, 75%, máximo
# Baseline (modelo simples para o problema):
X = df.drop("class", axis=1)
y = df["class"]
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.30, stratify=y, random_state=42) # 70% treino, 30% teste, random_state garante reprodutibilidade dos resultados (resultados iguais toda vez que executar o modelo)
y_train.value_counts(normalize=True) # Base treino: TP 0.896435 (89.64%), NT 0.103565 (10.36%)
y_test.value_counts(normalize=True) # Base testes: TP 0.897233 (89.72%), NT 0.102767 (10.28%)
lrc = LogisticRegression(random_state=42) # modelo de regressão logística
cv_list_lr_baseline = cross_val_score(
lrc,
X_train,
y_train,
cv=10,
scoring="balanced_accuracy"
) # 10-fold cross-validation, balanced_accuracy como métrica de avaliação
mean_cv_lr_baseline = np.mean(cv_list_lr_baseline) # média dos resultados de validação cruzada (bac) para modelo de regressão logística
std_cv_lr_baseline = np.std(cv_list_lr_baseline) # desvio padrão dos resultados de validação cruzada (bac) para modelo de regressão logística
print(f"Performance (bac): {round(mean_cv_lr_baseline, 4)} +- {round(std_cv_lr_baseline, 4)}") # Performance (bac): 0.9201 +- 0.046
# Modelagem:
knn = Pipeline(
[
('mms', MinMaxScaler()),
('skb', SelectKBest(chi2, k=10)),
('knn', KNeighborsClassifier(
n_neighbors=3,
p=2,
weights="uniform",
))
]
) # pipeline para modelo KNN, MinMaxScaler, SelectKBest e KNeighborsClassifier, k=10 (10 melhores features) e n_neighbors=3 (3 vizinhos), p=2 (distância euclidiana), weights="uniform" (pesos iguais para todos vizinhos)
cv_list_knn_euclid = cross_val_score(
knn,
X_train,
y_train,
cv=10,
scoring="balanced_accuracy"
) # 10-fold cross-validation, balanced_accuracy como métrica de avaliação para modelo KNN
mean_cv_knn_euclid = np.mean(cv_list_knn_euclid) # média dos resultados de validação cruzada (bac) para modelo KNN
std_cv_knn_euclid = np.std(cv_list_lr_baseline) # desvio padrão dos resultados de validação cruzada (bac) para modelo KNN
print(f"Performance (bac): {round(mean_cv_knn_euclid, 4)} +- {round(std_cv_knn_euclid, 4)}") # Performance (bac): 0.9703 +- 0.046
knn = Pipeline(
[
('mms', MinMaxScaler()),
('skb', SelectKBest(chi2, k=10)),
('knn', KNeighborsClassifier(
n_neighbors=3,
p=1,
weights="uniform",
))
]
) # pipeline para modelo KNN, MinMaxScaler, SelectKBest e KNeighborsClassifier, k=10 (10 melhores features) e n_neighbors=3 (3 vizinhos), p=1 (distância manhattan), weights="uniform" (pesos iguais para todos vizinhos)
# distância euclidiana utiliza linha reta para calcular distância entre 2 pontos, enquanto que distância manhattan utiliza somente linhas horizontais e verticais (mais efiiciente para dados com muitas dimensões e outliers)
cv_list_knn_man = cross_val_score(
knn,
X_train,
y_train,
cv=10,
scoring="balanced_accuracy"
) # 10-fold cross-validation, balanced_accuracy como métrica de avaliação para modelo KNN
mean_cv_knn_man = np.mean(cv_list_knn_man) # média dos resultados de validação cruzada (bac) para modelo
std_cv_knn_man = np.std(cv_list_knn_man) # desvio padrão dos resultados de validação cruzada (bac) para modelo
print(f"Performance (bac): {round(mean_cv_knn_man, 4)} +- {round(std_cv_knn_man, 4)}") # Performance (bac): 0.9638 +- 0.0407
lr = Pipeline(
[
('scaler', StandardScaler()),
('lr', LogisticRegression(
penalty="l2", # penalidade, usado para evitar overfitting
C=1, # força de regularização do modelo. Valores pequenos implicam em regularização mais forte
fit_intercept=True, # bias ou intercepto do modelo
class_weight="balanced", # peso das classes. Útil para datasets desbalanceados
random_state=42
)
)
])
cv_list_lr_l2 = cross_val_score(
lr,
X_train,
y_train,
cv=10,
scoring="balanced_accuracy"
) # 10-fold cross-validation, balanced_accuracy como métrica de avaliação para modelo de regressão logística com regularização L2
mean_cv_lr_l2 = np.mean(cv_list_lr_l2) # média dos resultados de validação cruzada (bac) para modelo de regressão logística com regularização L2
std_cv_lr_l2 = np.std(cv_list_lr_l2) # desvio padrão dos resultados de validação cruzada (bac) para modelo de regressão logística com regularização L2
print(f"Preformance (bac): {round(mean_cv_lr_l2, 4)} +- {round(std_cv_lr_l2, 4)}") # Preformance (bac): 0.9655 +- 0.0391
lr = Pipeline(
[
('scaler', StandardScaler()),
('lr', LogisticRegression(
penalty="l1", # penalidade, usado para evitar overfitting
C=1, # força de regularização do modelo. Valores pequenos implicam em regularização mais forte
fit_intercept=True, # bias ou intercepto do modelo
class_weight="balanced", # peso das classes. Útil para datasets desbalanceados
solver="liblinear",
random_state=42
)
)
])
cv_list_lr_l1 = cross_val_score(
lr,
X_train,
y_train,
cv=10,
scoring="balanced_accuracy"
) # 10-fold cross-validation, balanced_accuracy como métrica de avaliação para modelo de regressão logística com regularização L1
mean_cv_lr_l1 = np.mean(cv_list_lr_l1) # média dos resultados de validação cruzada (bac) para modelo de regressão logística com regularização L1
std_cv_lr_l1 = np.std(cv_list_lr_l1) # desvio padrão dos resultados de validação cruzada (bac) para modelo de regressão logística com regularização L1
print(f"Preformance (bac): {round(mean_cv_lr_l1, 4)} +- {round(std_cv_lr_l1, 4)}") # Preformance (bac): 0.9665 +- 0.0373
lr = Pipeline(
[
('scaler', StandardScaler()),
('pca', PCA(n_components=10)),
('lr', LogisticRegression(
penalty="l2",
C=1,
fit_intercept=True,
class_weight="balanced",
random_state=42
))
]
) # pipeline para modelo de regressão logística com PCA, StandardScaler e LogisticRegression, n_components=10 (10 componentes principais)
cv_list_lr_pca = cross_val_score(
lr,
X_train,
y_train,
cv=10,
scoring="balanced_accuracy"
) # 10-fold cross-validation, balanced_accuracy como métrica de avaliação para modelo de regressão logística com PCA e regularização L2
mean_cv_lr_pca = np.mean(cv_list_lr_pca) # média dos resultados de validação cruzada (bac) para modelo de regressão logística com PCA
std_cv_lr_pca = np.std(cv_list_lr_pca) # desvio padrão dos resultados de validação cruzada (bac) para modelo de regressão logística com PCA
print(f"Preformance (bac): {round(mean_cv_lr_pca, 4)} +- {round(std_cv_lr_pca, 4)}") # Preformance (bac): 0.9822 +- 0.0228
# Avaliação experimental:
# resultados da cross-validacao:
df_result_cv = pd.DataFrame(
[cv_list_lr_baseline, cv_list_knn_euclid, cv_list_knn_man, cv_list_lr_l2, cv_list_lr_l1, cv_list_lr_pca],
index=["baseline", "kNN-eucli", "kNN-man","LR-L2", "LR-L1", "LR-PCA"]
).T # dataframe com resultados de validação cruzada para cada modelo desenvolvido
df_result_cv # resultados de validação cruzada para cada modelo desenvolvido
df_res = df_result_cv.stack().to_frame("balanced_accuracy") # transformar dataframe em série e renomear coluna para balanced_accuracy (bac)
df_res.index.rename(["fold", "pipelines"], inplace=True) # renomear índices para fold e pipelines (modelos)
df_res = df_res.reset_index() # resetar índices do dataframe df_res para visualização dos resultados de validação cruzada para cada modelo desenvolvido (fold, pipelines, balanced_accuracy)
df_res.head(12) # visualizar primeiros 12 resultados de validação cruzada para cada modelo desenvolvido (fold, pipelines, balanced_accuracy)
plt.figure(figsize=(10,10)) # tamanho da figura
ax = sns.boxplot(x="pipelines", y="balanced_accuracy", data=df_res) # boxplot para visualização dos resultados de validação cruzada para cada modelo desenvolvido (pipelines, balanced_accuracy)
ax = sns.swarmplot(x="pipelines", y="balanced_accuracy", data=df_res, color=".40") # swarmplot para visualização dos resultados de validação cruzada para cada modelo desenvolvido (pipelines, balanced_accuracy)
plt.figure(figsize=(10,10))
sns.catplot(x="pipelines", y="balanced_accuracy", kind="violin", data=df_res) # violinplot para visualização dos resultados de validação cruzada para cada modelo desenvolvido (pipelines, balanced_accuracy)
# Avaliação de performance:
lr = Pipeline(
[
('scaler', StandardScaler()),
('pca', PCA(n_components=10)),
('lr', LogisticRegression(
penalty="l2",
C=1,
fit_intercept=True,
class_weight="balanced",
random_state=42
)
)
]) # pipeline para modelo de regressão logística com PCA, StandardScaler e LogisticRegression, n_components=10 (10 componentes principais) e regularização L2, C=1 (força de regularização), class_weight="balanced" (peso das classes), random_state=42 (garante reprodutibilidade dos resultados) e fit_intercept=True (bias ou intercepto do modelo)
lr.fit(X_train, y_train) # treinar modelo de regressão logística com PCA e regularização L2 para base de treino (X_train, y_train)
y_pred = lr.predict(X_test) # prever rótulos de resposta para base de teste (X_test), y_pred são rótulos preditos pelo modelo de regressão logística com PCA e regularização L2
lr_pca_test = balanced_accuracy_score(y_test, y_pred) # balanced_accuracy para modelo de regressão logística com PCA e regularização L2 na base de teste (X_test, y_test)
print("Performance: ", round(lr_pca_test, 4)) # Performance: 0.972 (97.2%)
from sklearn.metrics import ConfusionMatrixDisplay # Confusion matrix para avaliação de performance do modelo, matriz de confusão
ConfusionMatrixDisplay.from_estimator(lr, X_test, y_test)
plt.show() # exibir matriz de confusão para avaliação de performance do modelo
ConfusionMatrixDisplay.from_estimator(lr, X_test, y_test, normalize='true') # matriz de confusão normalizada (dados em porcentagem)
plt.show()
# Resultados confusion matrix:
true,predict
(errou) TP,NT: 4
(errou) TP,NT: 0.018 (0.018%)
(acertou) NT,NT: 25
(acertou) NT,NT: 0.96 (96%)
(acertou) TP,TP: 223
(acertou) TP,TP: 0.98 (98%)
(errou) NT,TP: 1
(errou) NT,TP: 0.038 (3.8%)
# Testar modelo com novos dados:
novo_dado_1 = np.array([[0.5, 0.8, 0.6, 0.9, 0.7, 0.8, 0.6, 0.9, 0.7, 0.8]]) # chance de apresentar tumor (TP)
predicao_1 = lr.predict(novo_dado_1)
print(f"Predição para novo dado 1: {predicao_1[0]}") # TP ou NT
novo_dado_2 = np.array([[0.1, 0.2, 0.3, 0.1, 0.2, 0.3, 0.1, 0.2, 0.3, 0.1]]) # sem chance de apresentar tumor (NT)
predicao_2 = lr.predict(novo_dado_2)
print(f"Predição para novo dado 2: {predicao_2[0]}") # TP ou NT
Elaborado por Mateus Schwede
ubsocial.github.io