<jupyter_start><jupyter_text>Extract features from T2 WSG Data<jupyter_code>import numpy as np
import matplotlib.pyplot as plt
from matplotlib.pyplot import figure
import pandas as pd
import time
from collections import OrderedDict

from sklearn.ensemble import RandomForestClassifier
from sklearn.model_selection import cross_val_score
from sklearn.model_selection import train_test_split
from sklearn import metrics
from sklearn.metrics import classification_report
from sklearn.metrics import roc_auc_score
from sklearn.metrics import roc_curve, auc
from sklearn.inspection import permutation_importance<jupyter_output><empty_output><jupyter_text>T2 BA + WSG PvL<jupyter_code>data = pd.read_csv('T2_WSG_BA_PvL.csv', index_col = 0)
data = data.T
data
y = data.iloc[:,0]
x = data.iloc[:,1:200]
X_train, X_test, y_train, y_test = train_test_split(x, y, test_size=.35, random_state=17)
ensemble_clfs = [
    ("RandomForestClassifier, max_features='sqrt'",
        RandomForestClassifier(warm_start=True, oob_score=True,
                               max_features="sqrt", 
                               random_state=42)),
    ("RandomForestClassifier, max_features='log2'",
        RandomForestClassifier(warm_start=True, max_features='log2',
                               oob_score=True, 
                               random_state=42)),
    ("RandomForestClassifier, max_features=None",
        RandomForestClassifier(warm_start=True, max_features=None,
                               oob_score=True, 
                               random_state=42))
]

# Map a classifier name to a list of (<n_estimators>, <error rate>) pairs.
error_rate = OrderedDict((label, []) for label, _ in ensemble_clfs)

# Range of `n_estimators` values to explore.
min_estimators = 100
max_estimators = 675

for label, clf in ensemble_clfs:
    for i in range(min_estimators, max_estimators + 1):
        clf.set_params(n_estimators=i)
        clf.fit(X_train, y_train)

        # Record the OOB error for each `n_estimators=i` setting.
        oob_error = 1 - clf.oob_score_
        error_rate[label].append((i, oob_error))

# Generate the "OOB error rate" vs. "n_estimators" plot.
for label, clf_err in error_rate.items():
    xs, ys = zip(*clf_err)
    plt.plot(xs, ys, label=label)

plt.xlim(min_estimators, max_estimators)
plt.xlabel("n_estimators")
plt.ylabel("OOB error rate")
plt.legend(loc="upper right")
plt.show()
clf = RandomForestClassifier(n_estimators=100, max_features=None, random_state=42)
all_accuracies = cross_val_score(estimator=clf, X=X_train, y=y_train, cv=5)
print('Mean Validation Scores: ' ,end='')
print(np.mean(all_accuracies))
clf_final = RandomForestClassifier(n_estimators=100, max_features=None, random_state=9, oob_score=True)
clf_final.fit(X_train,y_train)
y_pred = clf_final.predict(X_test)
print("Test Set Accuracy:",metrics.accuracy_score(y_test, y_pred))
# clf_final = RandomForestClassifier(n_estimators=250, max_features='sqrt', random_state=11, oob_score=True)
# clf_final.fit(X_train,y_train)
# y_pred = clf_final.predict(X_test)
# print("Test Set Accuracy:",metrics.accuracy_score(y_test, y_pred))
rf_probs = clf_final.predict_proba(X_test)[:, 1]
roc_value = roc_auc_score(y_test, rf_probs)
roc_value
# get roc/auc info
Y_score = clf_final.predict_proba(X_test)[:,1]
fpr = dict()
tpr = dict()
fpr, tpr, _ = roc_curve(y_test, Y_score)

roc_auc = dict()
roc_auc = auc(fpr, tpr)

# make the plot
plt.figure(figsize=(10,10))
plt.plot([0, 1], [0, 1], 'k--')
plt.xlim([-0.05, 1.0])
plt.ylim([0.0, 1.05])
plt.xlabel('False Positive Rate')
plt.ylabel('True Positive Rate')
plt.grid(True)
plt.plot(fpr, tpr, label='AUC = {0}'.format(roc_auc))        
plt.legend(loc="lower right", shadow=True, fancybox =True) 
plt.show()
print(clf_final.oob_score_)
feature_importance = pd.Series(clf_final.feature_importances_).sort_values(ascending=False)
feats = feature_importance[0:20]
feats.index
for i in feats.values:
    print(i)
features = []
WGS = 0
BA = 0
for i in feats.index:
    features.append(x.iloc[:,i].name)
    if i <= 26:
         WGS += 1
    else: 
        BA +=1
print(WGS)
print(BA)
features
result = permutation_importance(clf_final, X_test,y_test, n_repeats=10, random_state=42)
sorted_feature = result.importances_mean.argsort()[::-1]
top_10 = sorted_feature[0:10]
top_10
top_10
result.importances_mean[top_10]
result.importances_std[top_10]
fig,ax = plt.subplots()
ax.boxplot(result.importances[top_10].T, vert=False,labels=X_test.columns[top_10])
fig.tight_layout()
plt.show()
index = np.arange(len(top_10))
std = result.importances_std[top_10]
figure(figsize=(6,4))
plt.barh(index, result.importances_mean[top_10])
plt.yticks(index, X_test.columns[top_10], fontsize=10)
plt.title("T2 BA + WSG PvL", fontsize=20)<jupyter_output><empty_output><jupyter_text>T2 Metabolome + WSG PvL<jupyter_code>data = pd.read_csv('T2_Known_WSG_PvL.csv', index_col = 0)
data = data.T
data
y = data.iloc[:,0]
x = data.iloc[:,1:100]
X_train, X_test, y_train, y_test = train_test_split(x, y, test_size=.35, random_state=11)
ensemble_clfs = [
    ("RandomForestClassifier, max_features='sqrt'",
        RandomForestClassifier(warm_start=True, oob_score=True,
                               max_features="sqrt", 
                               random_state=42)),
    ("RandomForestClassifier, max_features='log2'",
        RandomForestClassifier(warm_start=True, max_features='log2',
                               oob_score=True, 
                               random_state=42)),
    ("RandomForestClassifier, max_features=None",
        RandomForestClassifier(warm_start=True, max_features=None,
                               oob_score=True, 
                               random_state=42))
]

# Map a classifier name to a list of (<n_estimators>, <error rate>) pairs.
error_rate = OrderedDict((label, []) for label, _ in ensemble_clfs)

# Range of `n_estimators` values to explore.
min_estimators = 100
max_estimators = 675

for label, clf in ensemble_clfs:
    for i in range(min_estimators, max_estimators + 1):
        clf.set_params(n_estimators=i)
        clf.fit(X_train, y_train)

        # Record the OOB error for each `n_estimators=i` setting.
        oob_error = 1 - clf.oob_score_
        error_rate[label].append((i, oob_error))

# Generate the "OOB error rate" vs. "n_estimators" plot.
for label, clf_err in error_rate.items():
    xs, ys = zip(*clf_err)
    plt.plot(xs, ys, label=label)

plt.xlim(min_estimators, max_estimators)
plt.xlabel("n_estimators")
plt.ylabel("OOB error rate")
plt.legend(loc="upper right")
plt.show()
clf = RandomForestClassifier(n_estimators=480, max_features='sqrt', random_state=42)
all_accuracies = cross_val_score(estimator=clf, X=X_train, y=y_train, cv=5)
print('Mean Validation Scores: ' ,end='')
print(np.mean(all_accuracies))
clf_final = RandomForestClassifier(n_estimators=480, max_features='sqrt', oob_score=True, random_state=10)
clf_final.fit(X_train,y_train)
y_pred = clf_final.predict(X_test)
print("Test Set Accuracy:",metrics.accuracy_score(y_test, y_pred))
# clf_final = RandomForestClassifier(n_estimators=150, max_features=None, oob_score=True, random_state=9)
# clf_final.fit(X_train,y_train)
# y_pred = clf_final.predict(X_test)
# print("Test Set Accuracy:",metrics.accuracy_score(y_test, y_pred))
rf_probs = clf_final.predict_proba(X_test)[:, 1]
roc_value = roc_auc_score(y_test, rf_probs)
roc_value
# get roc/auc info
Y_score = clf_final.predict_proba(X_test)[:,1]
fpr = dict()
tpr = dict()
fpr, tpr, _ = roc_curve(y_test, Y_score)

roc_auc = dict()
roc_auc = auc(fpr, tpr)

# make the plot
plt.figure(figsize=(10,10))
plt.plot([0, 1], [0, 1], 'k--')
plt.xlim([-0.05, 1.0])
plt.ylim([0.0, 1.05])
plt.xlabel('False Positive Rate')
plt.ylabel('True Positive Rate')
plt.grid(True)
plt.plot(fpr, tpr, label='AUC = {0}'.format(roc_auc))        
plt.legend(loc="lower right", shadow=True, fancybox =True) 
plt.show()
print(clf_final.oob_score_)
feature_importance = pd.Series(clf_final.feature_importances_).sort_values(ascending=False)
feats = feature_importance[0:20]
feats.index
features = []
WGS = 0
Meta = 0
for i in feats.index:
    features.append(x.iloc[:,i].name)
    if i <= 26:
         WGS += 1
    else: 
        Meta +=1
print(WGS)
print(Meta)
for i in feats.values:
    print(i)
features
result = permutation_importance(clf_final, X_test,y_test, n_repeats=10, random_state=42)
sorted_feature = result.importances_mean.argsort()[::-1]
top_10 = sorted_feature[0:10]
top_10
fig,ax = plt.subplots()
ax.boxplot(result.importances[top_10].T, vert=False,labels=X_test.columns[top_10])
fig.tight_layout()
plt.show()<jupyter_output><empty_output><jupyter_text>T2 BA + WSG LvCoL<jupyter_code>data = pd.read_csv('T2_WSG_BA_LvCoL.csv', index_col = 0)
data = data.T
data
y = data.iloc[:,0]
x = data.iloc[:,1:200]
X_train, X_test, y_train, y_test = train_test_split(x, y, test_size=.35, random_state=17)
ensemble_clfs = [
    ("RandomForestClassifier, max_features='sqrt'",
        RandomForestClassifier(warm_start=True, oob_score=True,
                               max_features="sqrt", 
                               random_state=42)),
    ("RandomForestClassifier, max_features='log2'",
        RandomForestClassifier(warm_start=True, max_features='log2',
                               oob_score=True, 
                               random_state=42)),
    ("RandomForestClassifier, max_features=None",
        RandomForestClassifier(warm_start=True, max_features=None,
                               oob_score=True, 
                               random_state=42))
]

# Map a classifier name to a list of (<n_estimators>, <error rate>) pairs.
error_rate = OrderedDict((label, []) for label, _ in ensemble_clfs)

# Range of `n_estimators` values to explore.
min_estimators = 100
max_estimators = 675

for label, clf in ensemble_clfs:
    for i in range(min_estimators, max_estimators + 1):
        clf.set_params(n_estimators=i)
        clf.fit(X_train, y_train)

        # Record the OOB error for each `n_estimators=i` setting.
        oob_error = 1 - clf.oob_score_
        error_rate[label].append((i, oob_error))

# Generate the "OOB error rate" vs. "n_estimators" plot.
for label, clf_err in error_rate.items():
    xs, ys = zip(*clf_err)
    plt.plot(xs, ys, label=label)

plt.xlim(min_estimators, max_estimators)
plt.xlabel("n_estimators")
plt.ylabel("OOB error rate")
plt.legend(loc="upper right")
plt.show()
clf = RandomForestClassifier(n_estimators=400, max_features= None, random_state=42)
all_accuracies = cross_val_score(estimator=clf, X=X_train, y=y_train, cv=5)
print('Mean Validation Scores: ' ,end='')
print(np.mean(all_accuracies))
clf_final = RandomForestClassifier(n_estimators=400, max_features=None, oob_score=True, random_state=9)
clf_final.fit(X_train,y_train)
y_pred = clf_final.predict(X_test)
print("Test Set Accuracy:",metrics.accuracy_score(y_test, y_pred))
# clf_final = RandomForestClassifier(n_estimators=100, max_features=None, oob_score=True, random_state=9)
# clf_final.fit(X_train,y_train)
# y_pred = clf_final.predict(X_test)
# print("Test Set Accuracy:",metrics.accuracy_score(y_test, y_pred))
rf_probs = clf_final.predict_proba(X_test)[:, 1]
roc_value = roc_auc_score(y_test, rf_probs)
roc_value
# get roc/auc info
Y_score = clf_final.predict_proba(X_test)[:,1]
fpr = dict()
tpr = dict()
fpr, tpr, _ = roc_curve(y_test, Y_score)

roc_auc = dict()
roc_auc = auc(fpr, tpr)

# make the plot
plt.figure(figsize=(10,10))
plt.plot([0, 1], [0, 1], 'k--')
plt.xlim([-0.05, 1.0])
plt.ylim([0.0, 1.05])
plt.xlabel('False Positive Rate')
plt.ylabel('True Positive Rate')
plt.grid(True)
plt.plot(fpr, tpr, label='AUC = {0}'.format(roc_auc))        
plt.legend(loc="lower right", shadow=True, fancybox =True) 
plt.show()
print(clf_final.oob_score_)
feature_importance = pd.Series(clf_final.feature_importances_).sort_values(ascending=False)
feats = feature_importance[0:20]
feats.index
features = []
WGS = 0
BA = 0
for i in feats.index:
    features.append(x.iloc[:,i].name)
    if i <= 26:
         WGS += 1
    else: 
        BA +=1
print(WGS)
print(BA)
for i in feats.values:
    print(i)
features
result = permutation_importance(clf_final, X_test,y_test, n_repeats=10, random_state=42)
sorted_feature = result.importances_mean.argsort()[::-1]
top_10 = sorted_feature[0:10]
top_10
result.importances_mean
index = np.arange(len(top_10))
std = result.importances_std[top_10]
figure(figsize=(6,4))
plt.barh(index, result.importances_mean[top_10])
plt.yticks(index, X_test.columns[top_10], fontsize=10)
plt.title("T2 BA + WSG LvCoL", fontsize=20)<jupyter_output><empty_output><jupyter_text>T2 Metabolome + WSG LvCoL<jupyter_code>data = pd.read_csv('T2_Known_WSG_LvCoL.csv', index_col = 0)
data = data.T
data
y = data.iloc[:,0]
x = data.iloc[:,1:100]
X_train, X_test, y_train, y_test = train_test_split(x, y, test_size=.35, random_state=11)
ensemble_clfs = [
    ("RandomForestClassifier, max_features='sqrt'",
        RandomForestClassifier(warm_start=True, oob_score=True,
                               max_features="sqrt", 
                               random_state=42)),
    ("RandomForestClassifier, max_features='log2'",
        RandomForestClassifier(warm_start=True, max_features='log2',
                               oob_score=True, 
                               random_state=42)),
    ("RandomForestClassifier, max_features=None",
        RandomForestClassifier(warm_start=True, max_features=None,
                               oob_score=True, 
                               random_state=42))
]

# Map a classifier name to a list of (<n_estimators>, <error rate>) pairs.
error_rate = OrderedDict((label, []) for label, _ in ensemble_clfs)

# Range of `n_estimators` values to explore.
min_estimators = 100
max_estimators = 675

for label, clf in ensemble_clfs:
    for i in range(min_estimators, max_estimators + 1):
        clf.set_params(n_estimators=i)
        clf.fit(X_train, y_train)

        # Record the OOB error for each `n_estimators=i` setting.
        oob_error = 1 - clf.oob_score_
        error_rate[label].append((i, oob_error))

# Generate the "OOB error rate" vs. "n_estimators" plot.
for label, clf_err in error_rate.items():
    xs, ys = zip(*clf_err)
    plt.plot(xs, ys, label=label)

plt.xlim(min_estimators, max_estimators)
plt.xlabel("n_estimators")
plt.ylabel("OOB error rate")
plt.legend(loc="upper right")
plt.show()
clf = RandomForestClassifier(n_estimators=300, max_features=None, random_state=42)
all_accuracies = cross_val_score(estimator=clf, X=X_train, y=y_train, cv=5)
print('Mean Validation Scores: ' ,end='')
print(np.mean(all_accuracies))
clf_final = RandomForestClassifier(n_estimators=300, max_features=None, oob_score=True, random_state=9)
clf_final.fit(X_train,y_train)
y_pred = clf_final.predict(X_test)
print("Test Set Accuracy:",metrics.accuracy_score(y_test, y_pred))
# clf_final = RandomForestClassifier(n_estimators=100, max_features=None, oob_score=True, random_state=9)
# clf_final.fit(X_train,y_train)
# y_pred = clf_final.predict(X_test)
# print("Test Set Accuracy:",metrics.accuracy_score(y_test, y_pred))
rf_probs = clf_final.predict_proba(X_test)[:, 1]
roc_value = roc_auc_score(y_test, rf_probs)
roc_value
# get roc/auc info
Y_score = clf_final.predict_proba(X_test)[:,1]
fpr = dict()
tpr = dict()
fpr, tpr, _ = roc_curve(y_test, Y_score)

roc_auc = dict()
roc_auc = auc(fpr, tpr)

# make the plot
plt.figure(figsize=(10,10))
plt.plot([0, 1], [0, 1], 'k--')
plt.xlim([-0.05, 1.0])
plt.ylim([0.0, 1.05])
plt.xlabel('False Positive Rate')
plt.ylabel('True Positive Rate')
plt.grid(True)
plt.plot(fpr, tpr, label='AUC = {0}'.format(roc_auc))        
plt.legend(loc="lower right", shadow=True, fancybox =True) 
plt.show()
print(clf_final.oob_score_)
feature_importance = pd.Series(clf_final.feature_importances_).sort_values(ascending=False)
feats = feature_importance[0:20]
feats.index
features = []
WGS = 0
Meta = 0
for i in feats.index:
    features.append(x.iloc[:,i].name)
    if i <= 26:
         WGS += 1
    else: 
        Meta +=1
print(WGS)
print(Meta)
for i in feats.values:
    print(i)
features<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.13081792230358694, 'fraction_numerical': 0.023023457862728063, 'mean_word_length': 4.429582210242588, 'pattern_counts': {'":': 0, '<': 33, '<?xml version=': 0, '>': 29, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '5', 'dir': 'jupyter-structured-clean-dedup', 'id': '588069', 'license': 'MIT', 'n_tokens_mistral': 6471, 'n_tokens_neox': 6136, 'n_words': 1129, 'path': 'Random_Forest/Multi-omics/T2/T2_RF_Top_Log2fold.ipynb', 'repo_name': 'bryansho/PCOS_WGS_16S_metabolome', 'size': '16114'}