CCLE_Transfer.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "VsAnMTpNAV-h"
   },
   "source": [
    "## Transfer learning from CCLE to MMRF (V6)\n",
    "\n",
    "Last Modified: 05/19/2021\n",
    "\n",
    "Written By: Sumi Thakur and Qingyang Xu\n",
    "\n",
    "- Pretrain on CCLE cell line data to predict IC50\n",
    "\n",
    "- Transfer to patient RNA-seq data to predict PFS\n",
    "\n",
    "References\n",
    "\n",
    "- Download patient genomic data (e.g. `MMRF_CoMMpass_IA15a_CNA_Exome_PerGene_LargestSegment.txt`)\n",
    "\n",
    "https://research.themmrf.org/\n",
    "\n",
    "- Download DevMap cell line data (e.g. `CCLE_expression.csv`)\n",
    "\n",
    "https://depmap.org/portal/download/"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/"
    },
    "id": "Dl5gEjZlpx9L",
    "outputId": "2e23cc79-87bb-4beb-94b4-b643f4a629ac"
   },
   "outputs": [],
   "source": [
    "#from google.colab import drive\n",
    "import os\n",
    "import pandas as pd\n",
    "from sklearn.model_selection import train_test_split\n",
    "import tensorflow as tf\n",
    "from kerastuner import HyperModel\n",
    "from keras import models, layers,regularizers\n",
    "from tensorflow.keras.optimizers import Adam\n",
    "from kerastuner.tuners.bayesian import BayesianOptimization\n",
    "from tqdm import tqdm\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from sklearn.metrics import r2_score\n",
    "from scipy.stats import pearsonr\n",
    "import pickle"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "id": "1JXuEqQZAV-r"
   },
   "outputs": [],
   "source": [
    "from sklearn.linear_model import LinearRegression, LogisticRegression\n",
    "from sklearn.ensemble import RandomForestClassifier, GradientBoostingClassifier\n",
    "from sklearn import neural_network\n",
    "from sklearn.model_selection import GridSearchCV, RandomizedSearchCV\n",
    "\n",
    "from sklearn import metrics as skmetrics\n",
    "from sklearn.decomposition import PCA"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 52,
   "metadata": {
    "colab": {
     "base_uri": "https://localhost:8080/",
     "height": 310
    },
    "id": "faHIU_mAAV-r",
    "outputId": "36584a3a-5532-470a-d78a-c8ce66b472cb"
   },
   "outputs": [],
   "source": [
    "from utils import MMRF_utils_V2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import sklearn \n",
    "print(tf.__version__, sklearn.__version__)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "BR2-PSkBqB7K",
    "outputId": "97e69f46-c455-4807-f3de-90e336c6d94b"
   },
   "source": [
    "## 1. Pretrain on CCLE IC50 data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "id": "gK362ObPAV-s"
   },
   "outputs": [],
   "source": [
    "def normalize_features(train, test):\n",
    "    \n",
    "    print('Normalizing input features...')\n",
    "    nsample, nfeature = train.shape\n",
    "    \n",
    "    assert nfeature == test.shape[1]\n",
    "\n",
    "    for i in range(nfeature):\n",
    "        mu = np.mean(train[:,i])\n",
    "        sigma = np.std(train[:,i])\n",
    "        train[:,i] -= mu\n",
    "        test[:,i] -= mu\n",
    "        \n",
    "        if sigma > 0: # some genes may have zero variance\n",
    "            train[:,i] /= sigma\n",
    "            test[:,i] /= sigma\n",
    "    \n",
    "    return train, test"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "id": "N29jLNJ4qMi6"
   },
   "outputs": [],
   "source": [
    "def preprocess_CCLE(ccle_exp, sample_info, sanger_dose_response, mapping):\n",
    "    \n",
    "    print('Preprocessing CCLE data...')\n",
    "    \n",
    "    ccle_exp = ccle_exp.rename(columns={'Unnamed: 0':'DepMap_ID'})\n",
    "    cols_to_keep = ['DepMap_ID'] + list(pd.unique(mapping.HGNC_ID))\n",
    "    ccle_chosen = ccle_exp[cols_to_keep]\n",
    "\n",
    "    sanger_dose_response_filt = sanger_dose_response[(sanger_dose_response.DRUG_NAME.str.contains(\"BORTEZOMIB\"))\n",
    "                                                | (sanger_dose_response.DRUG_NAME.str.contains(\"LENALIDOMIDE\"))\n",
    "                                                ]\n",
    "\n",
    "    sanger_dose_response_filt = sanger_dose_response_filt.rename(columns={'ARXSPAN_ID':'DepMap_ID'})\n",
    "    \n",
    "    merged_df = sanger_dose_response_filt.merge(ccle_chosen, on= 'DepMap_ID')[list(ccle_chosen.columns)+['DRUG_NAME','IC50_PUBLISHED']]\n",
    "    merged_df['log(IC_50)'] = merged_df.IC50_PUBLISHED.apply(np.log10)\n",
    "    \n",
    "    merged_df_bort = merged_df[merged_df.DRUG_NAME == 'BORTEZOMIB'].drop_duplicates()\n",
    "    merged_df_lenal = merged_df[merged_df.DRUG_NAME == 'LENALIDOMIDE'].drop_duplicates()\n",
    "    bort_labels = merged_df_bort['log(IC_50)']\n",
    "    lenal_labels = merged_df_lenal['log(IC_50)']\n",
    "    bort_data = merged_df_bort.drop(columns = ['DepMap_ID','IC50_PUBLISHED','DRUG_NAME','log(IC_50)'])\n",
    "    lenal_data = merged_df_lenal.drop(columns = ['DepMap_ID','IC50_PUBLISHED','DRUG_NAME','log(IC_50)'])\n",
    "    \n",
    "    return bort_data, bort_labels, lenal_data, lenal_labels"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "id": "v5P_TPRVqj1g"
   },
   "outputs": [],
   "source": [
    "def split_train_test_CCLE(X, y, normalize=True, nPCA=0):\n",
    "    \n",
    "    X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=0, train_size=0.7)\n",
    "    \n",
    "    X_train = np.array(X_train)\n",
    "    y_train = np.array(y_train)\n",
    "    \n",
    "    X_test = np.array(X_test)\n",
    "    y_test = np.array(y_test)\n",
    "    \n",
    "    X_train, X_test = normalize_features(X_train, X_test)\n",
    "    \n",
    "    pca = None\n",
    "    # try running PCA on raw RNAseq data\n",
    "    if nPCA > 0:\n",
    "        print('Running PCA-%d on gene expressions...'%nPCA)\n",
    "        pca = PCA(n_components=nPCA, random_state=0)\n",
    "        X_train = pca.fit_transform(X_train)\n",
    "        X_test = pca.transform(X_test)\n",
    "    \n",
    "    assert X_train.shape[0] == y_train.shape[0]\n",
    "    \n",
    "    assert X_test.shape[0] == y_test.shape[0]\n",
    "    \n",
    "    print('Training set shape:')\n",
    "    print(X_train.shape)\n",
    "    \n",
    "    print('Test set shape:')\n",
    "    print(X_test.shape)\n",
    "    \n",
    "    return X_train, y_train, X_test, y_test, pca"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "id": "QBp-hX0IAV-u"
   },
   "outputs": [],
   "source": [
    "def linear_regression_IC50(X_train, y_train, X_test, y_test):\n",
    "    print('Running linear regression to predict IC50...')\n",
    "    reg = LinearRegression().fit(X_train, y_train)\n",
    "    y_pred = reg.predict(X_test)\n",
    "\n",
    "    R2 = r2_score(y_test, y_pred)\n",
    "    print('R2: %0.2f'%R2)\n",
    "\n",
    "    rho, p = pearsonr(y_test, y_pred)\n",
    "    print('Correlation: %0.2f'%rho)\n",
    "    print('p-value: %0.5f'%p)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "class RegressionHyperModel(HyperModel):\n",
    "    def __init__(self, input_shape, enc_layer_size=5):\n",
    "        self.input_shape = input_shape\n",
    "        self.enc_layer_size = enc_layer_size\n",
    "    def build(self, hp):\n",
    "        model = tf.keras.Sequential()\n",
    "        model.add(\n",
    "            tf.keras.layers.Dense(\n",
    "                units=hp.Int('units', 40, 100, step=5, default=80),\n",
    "                activation=hp.Choice(\n",
    "                    'dense_activation',\n",
    "                    values=['relu', 'tanh', 'sigmoid'],\n",
    "                    default='relu'),\n",
    "                name = 'dense1',\n",
    "                input_shape = self.input_shape\n",
    "            )\n",
    "        )\n",
    "        \n",
    "        model.add(\n",
    "            tf.keras.layers.Dropout(\n",
    "                hp.Float(\n",
    "                    'dropout',\n",
    "                    min_value=0.0,\n",
    "                    max_value=0.6,\n",
    "                    default=0.3,\n",
    "                    step=0.1),\n",
    "                name='dropout1'\n",
    "            )\n",
    "        )\n",
    "        \n",
    "        model.add(\n",
    "            tf.keras.layers.Dense(\n",
    "                units=hp.Int('units', 10, 50, step=5, default=40),\n",
    "                activation=hp.Choice(\n",
    "                    'dense_activation',\n",
    "                    values=['relu', 'tanh', 'sigmoid'],\n",
    "                    default='relu'),\n",
    "                name='dense2'\n",
    "            )\n",
    "        )\n",
    "        \n",
    "        model.add(\n",
    "            tf.keras.layers.Dropout(\n",
    "                hp.Float(\n",
    "                    'dropout',\n",
    "                    min_value=0.0,\n",
    "                    max_value=0.6,\n",
    "                    default=0.3,\n",
    "                    step=0.1),\n",
    "                name='dropout2'\n",
    "            )\n",
    "        )\n",
    "        \n",
    "        model.add(tf.keras.layers.Dense(self.enc_layer_size, name='encoder'))\n",
    "        \n",
    "        model.add(tf.keras.layers.Dense(1))\n",
    "        \n",
    "        model.compile(\n",
    "            optimizer='Adam',loss='mse',metrics=['mse']\n",
    "        )\n",
    "        \n",
    "        return model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "id": "4l4DavAnwphv"
   },
   "outputs": [],
   "source": [
    "def pretrain_MLP_model(X_train, y_train, X_test, y_test, params):\n",
    "    \n",
    "    print('Training MLP model...')\n",
    "\n",
    "    enc_layer_size = params['enc_layer']\n",
    "    \n",
    "    input_shape = (X_train.shape[1],)\n",
    "    print(input_shape)\n",
    "    hypermodel = RegressionHyperModel(input_shape, enc_layer_size)    \n",
    "    \n",
    "    ### using bayesian hyperparamter optimization\n",
    "    tuner_bo = BayesianOptimization(\n",
    "                hypermodel,\n",
    "                objective='mse',\n",
    "                max_trials=10,\n",
    "                seed=42,\n",
    "                executions_per_trial=2\n",
    "            )\n",
    "    \n",
    "    ### set epochs\n",
    "    tuner_bo.search(X_train, y_train, epochs=10, validation_split=0.2, verbose=1)\n",
    "    best_model = tuner_bo.get_best_models(num_models=1)[0]\n",
    "    best_model.get_config()\n",
    "  \n",
    "    return best_model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "id": "6VpQBcyEAV-w"
   },
   "outputs": [],
   "source": [
    "def plot_R2_correlation(R2_his, rho_his):\n",
    "\n",
    "    fig, ax1 = plt.subplots()\n",
    "\n",
    "    color = 'tab:red'\n",
    "    ax1.set_xlabel('Epoch')\n",
    "    ax1.set_ylabel(r'$R^2$', color=color)\n",
    "    ax1.plot(R2_his, color=color)\n",
    "    ax1.tick_params(axis='y', labelcolor=color)\n",
    "\n",
    "    ax2 = ax1.twinx()  # instantiate a second axes that shares the same x-axis\n",
    "\n",
    "    color = 'tab:blue'\n",
    "    ax2.set_ylabel(r'$\\rho$', color=color)  # we already handled the x-label with ax1\n",
    "    ax2.plot(rho_his, color=color)\n",
    "    ax2.tick_params(axis='y', labelcolor=color)\n",
    "\n",
    "    fig.tight_layout()  # otherwise the right y-label is slightly clipped\n",
    "    plt.title(r'$R^2$ and $\\rho$ during training')\n",
    "    plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "id": "jSBPPvqXAV-x"
   },
   "outputs": [],
   "source": [
    "def plot_losses(train_loss,val_loss):\n",
    "\n",
    "    plt.yscale('log')\n",
    "    plt.plot(train_loss[1:], label='Train Loss')\n",
    "    plt.plot(val_loss[1:], label='Valid Loss')\n",
    "    plt.legend()\n",
    "    plt.title('Train vs. Valid Losses during Training')\n",
    "    plt.xlabel('Epoch')\n",
    "    plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "ofFt5StIAV-x"
   },
   "source": [
    "## 2. Transfer to MMRF patient data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "id": "VOhK6SyXAV-x"
   },
   "outputs": [],
   "source": [
    "def preprocess_MMRF_genomic(genomic_fn, CCLE_cols, mapping):\n",
    "\n",
    "    print('\\nPreprocess patient genomic data in MMRF...')\n",
    "    genomic_df = MMRF_utils_V2.preprocess_genomic_data(genomic_fn, nPCA=0)\n",
    "\n",
    "    mapping_dict = {}\n",
    "\n",
    "    # create mapping dict from Ensembl to HGNC\n",
    "    for i in range(mapping.shape[0]):\n",
    "        eid = mapping['Ensembl_ID'].iloc[i]\n",
    "        hid = mapping['HGNC_ID'].iloc[i]\n",
    "        mapping_dict[eid] = hid\n",
    "        \n",
    "    # average over columns with the same Ensembl ID\n",
    "    gene_df = genomic_df[mapping['Ensembl_ID']]\n",
    "    pids = genomic_df['pids']\n",
    "\n",
    "    HGNC_cols = [mapping_dict[eid] for eid in gene_df.columns]\n",
    "    gene_df.columns = HGNC_cols\n",
    "    # average over columns with the same HGNC name\n",
    "    gene_df = gene_df.groupby(by=gene_df.columns, axis=1).mean()\n",
    "    \n",
    "    gene_df = gene_df[CCLE_cols]\n",
    "    \n",
    "    # normalize each gene\n",
    "    for col in gene_df.columns:\n",
    "        mu = gene_df[col].mean()\n",
    "        sigma = gene_df[col].std()\n",
    "        gene_df[col] -= mu\n",
    "        if sigma>0: gene_df[col] /= sigma\n",
    "    \n",
    "    return gene_df, pids"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "p6YkwuwUAV-y"
   },
   "source": [
    "## 3. Comparing AUC of three approaches"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "id": "E79nyLAsAV-y"
   },
   "outputs": [],
   "source": [
    "# no grid search\n",
    "def predict_PFS(clf, model_name, train_X, train_y, test_X, test_y, training=True):\n",
    "    \n",
    "    if training:\n",
    "        print('Fitting model %s...'%model_name)\n",
    "        clf.fit(train_X, train_y)\n",
    "\n",
    "    # predict PoS in test set\n",
    "    print('Predicting on test set...')\n",
    "    y_pred = clf.predict_proba(test_X)[:,1]\n",
    "\n",
    "    auc = skmetrics.roc_auc_score(test_y, y_pred)\n",
    "    f1 = skmetrics.f1_score(test_y, y_pred.round())\n",
    "    pres = skmetrics.precision_score(test_y, y_pred.round())\n",
    "    rec = skmetrics.recall_score(test_y, y_pred.round())\n",
    "    acc = skmetrics.accuracy_score(test_y, y_pred.round())\n",
    "\n",
    "    print('Accuracy: %0.2f'%acc)\n",
    "    print('Precision: %0.2f'%pres)\n",
    "    print('Recall: %0.2f'%rec)\n",
    "    print('F1 score: %0.2f'%f1)\n",
    "    print('AUC: %0.2f'%auc)\n",
    "    print('')\n",
    "    \n",
    "    return clf"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "id": "yTYigO1PAV-z"
   },
   "outputs": [],
   "source": [
    "def predict_PFS_grid_search(clf, model_name, train_X, train_y, test_X, test_y):\n",
    "    \n",
    "    print('Commencing grid search for %s'%model_name)\n",
    "\n",
    "    grid = {\"max_depth\": [10,20,50],\n",
    "            \"n_estimators\":[1000,5000]}\n",
    "        \n",
    "    clf_search = RandomizedSearchCV(clf, param_distributions=grid, cv=5, random_state=50, scoring='roc_auc')\n",
    "    clf_search.fit(train_X, train_y)\n",
    "\n",
    "    print('Random search complete, fitting model %s...'%model_name)\n",
    "    best_model =  clf_search.best_estimator_\n",
    "    best_model.fit(train_X, train_y)\n",
    "\n",
    "    # predict PoS in test set\n",
    "    print('Predicting on test set...')\n",
    "    y_pred = best_model.predict_proba(test_X)[:,1]\n",
    "\n",
    "    auc = skmetrics.roc_auc_score(test_y, y_pred)\n",
    "    f1 = skmetrics.f1_score(test_y, y_pred.round())\n",
    "    pres = skmetrics.precision_score(test_y, y_pred.round())\n",
    "    rec = skmetrics.recall_score(test_y, y_pred.round())\n",
    "    acc = skmetrics.accuracy_score(test_y, y_pred.round())\n",
    "\n",
    "    print('Accuracy: %0.2f'%acc)\n",
    "    print('Precision: %0.2f'%pres)\n",
    "    print('Recall: %0.2f'%rec)\n",
    "    print('F1 score: %0.2f'%f1)\n",
    "    print('AUC: %0.2f'%auc)\n",
    "    print('')\n",
    "    \n",
    "    return best_model"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "K9xa7NVIAV-y"
   },
   "source": [
    "### 3.1 Only use patient data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "def stratify_patients(data, drug='All', option='', mmrf_asct=None):\n",
    "    \n",
    "    if drug != 'All':\n",
    "        data = data[data[drug+'1']==1]\n",
    "        \n",
    "    if option == 'ASCT':\n",
    "        data = data.merge(mmrf_asct[mmrf_asct['ASCT']==1], left_on='pids', right_on='pids')\n",
    "        data.pop('ASCT')\n",
    "    \n",
    "    if option == 'NON-ASCT':\n",
    "        data = data.merge(mmrf_asct[mmrf_asct['ASCT']==0], left_on='pids', right_on='pids')\n",
    "        data.pop('ASCT')\n",
    "\n",
    "    if 'age' in option:\n",
    "        option = option.split('_')[1]\n",
    "        if '+' in option:\n",
    "            age = int(option[:-1])\n",
    "            print('Age above %d'%age)            \n",
    "            data = data[data['age']>= age]\n",
    "            \n",
    "        if '-' in option:\n",
    "            age = int(option[:-1])\n",
    "            print('Age below %d'%age)            \n",
    "            data = data[data['age']< age]\n",
    "        \n",
    "        \n",
    "    return data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "id": "2LiVWQ6eAV-z"
   },
   "outputs": [],
   "source": [
    "def predict_pfs_patient_only(train_df,valid_df,test_df, clf, model_name='',\n",
    "                             train_drug='All', test_drug='All', option='', mmrf_asct=None, grid_search=False, training=True):\n",
    "\n",
    "    print('Predicting PFS using patient clinical data only...')\n",
    "    \n",
    "    assert train_drug in ['Bor','Len','All'] and test_drug in ['Bor','Len','All']\n",
    "    \n",
    "    print('Predicting patients with (%s, %s) as first line treatment...'%(train_drug,test_drug))\n",
    "    \n",
    "    # no genomic data\n",
    "    train = pd.concat([train_df,valid_df])\n",
    "    test = test_df\n",
    "    \n",
    "    train = stratify_patients(train, drug=train_drug, option=option, mmrf_asct=mmrf_asct)\n",
    "    test = stratify_patients(test, drug=test_drug, option=option, mmrf_asct=mmrf_asct)\n",
    "        \n",
    "    print(train.shape)    \n",
    "    print(test.shape)\n",
    "    \n",
    "    pfs_thresh = 12\n",
    "    train_y = np.array(train.pop('pfs'))\n",
    "    train_y = np.array([int(pfs>pfs_thresh) for pfs in train_y])\n",
    "    p1 = train.pop('pids')\n",
    "    train_X = np.array(train)\n",
    "\n",
    "    test_y = np.array(test.pop('pfs'))\n",
    "    test_y = np.array([int(pfs>pfs_thresh) for pfs in test_y])\n",
    "    p0 = test.pop('pids')\n",
    "    test_X = np.array(test)\n",
    "    \n",
    "    if grid_search:\n",
    "        clf_fit = predict_PFS_grid_search(clf, model_name, train_X, train_y, test_X, test_y)\n",
    "    else:\n",
    "        clf_fit = predict_PFS(clf, model_name, train_X, train_y, test_X, test_y, training)\n",
    "        \n",
    "    return clf_fit"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "3ijqqr15AV-0"
   },
   "source": [
    "### 3.2 PCA on raw RNA"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "id": "Dc8dyvDGAV-0"
   },
   "outputs": [],
   "source": [
    "def predict_pfs_patient_PCA(genomic_fn, nPCA, train_df,valid_df, test_df, clf, model_name='',\n",
    "                            train_drug='All', test_drug='Bor', option='', mmrf_asct=None, grid_search=False, training=True):\n",
    "    \n",
    "    print('Predicting PFS using patient clinical data and PCA-%d on RNAseq...'%nPCA)\n",
    "    \n",
    "    assert nPCA > 0    \n",
    "    assert train_drug in ['Bor','Len','All'] and test_drug in ['Bor','Len','All']\n",
    "    \n",
    "    print('Predicting patients with (%s, %s) as first line treatment...'%(train_drug,test_drug))\n",
    "    \n",
    "    print('\\nPreprocess patient genomic data in MMRF...')\n",
    "    \n",
    "    genomic_df = MMRF_utils_V2.preprocess_genomic_data(genomic_fn, nPCA=nPCA)\n",
    "    print(genomic_df.shape)\n",
    "    \n",
    "    print('\\nMerging patient data with genomic data...')\n",
    "    train_patient_genomic = train_df.merge(genomic_df, left_on='pids', right_on='pids')\n",
    "    print(train_patient_genomic.shape)\n",
    "\n",
    "    valid_patient_genomic = valid_df.merge(genomic_df, left_on='pids', right_on='pids')\n",
    "    print(valid_patient_genomic.shape)\n",
    "\n",
    "    test_patient_genomic = test_df.merge(genomic_df, left_on='pids', right_on='pids')\n",
    "    print(test_patient_genomic.shape)\n",
    "    \n",
    "    train = pd.concat([train_patient_genomic,valid_patient_genomic])\n",
    "    train = train\n",
    "    test = test_patient_genomic\n",
    "    \n",
    "    train = stratify_patients(train, drug=train_drug, option=option, mmrf_asct=mmrf_asct)\n",
    "    test = stratify_patients(test, drug=test_drug, option=option, mmrf_asct=mmrf_asct)\n",
    "\n",
    "    print(train.shape)\n",
    "    print(test.shape)\n",
    "\n",
    "    pfs_thresh = 12\n",
    "    train_y = np.array(train.pop('pfs'))\n",
    "    train_y = np.array([int(pfs>pfs_thresh) for pfs in train_y])\n",
    "    p1 = train.pop('pids')\n",
    "    train_X = np.array(train)\n",
    "\n",
    "    test_y = np.array(test.pop('pfs'))\n",
    "    test_y = np.array([int(pfs>pfs_thresh) for pfs in test_y])\n",
    "    p0 = test.pop('pids')\n",
    "    test_X = np.array(test)\n",
    "    \n",
    "    if grid_search:\n",
    "        clf_fit = predict_PFS_grid_search(clf, model_name, train_X, train_y, test_X, test_y)\n",
    "    else:\n",
    "        clf_fit = predict_PFS(clf, model_name, train_X, train_y, test_X, test_y, training)\n",
    "\n",
    "    return clf_fit"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "ERhw6vN_AV-1"
   },
   "source": [
    "### 3.3 Transfer from CCLE"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "id": "YiO9FF0OAV-1"
   },
   "outputs": [],
   "source": [
    "def predict_pfs_patient_CCLE(enc_gene_df, train_df,valid_df, test_df, clf, model_name='',\n",
    "                            train_drug='All', test_drug='Bor', option='', mmrf_asct=None, grid_search=False, training=True):\n",
    "    \n",
    "    print('\\nPredicing PFS with transfer model from CCLE...')\n",
    "    \n",
    "    print('\\nMerging patient data with genomic data...')\n",
    "    train_patient_genomic = train_df.merge(enc_gene_df, left_on='pids', right_on='pids')\n",
    "    print(train_patient_genomic.shape)\n",
    "\n",
    "    valid_patient_genomic = valid_df.merge(enc_gene_df, left_on='pids', right_on='pids')\n",
    "    print(valid_patient_genomic.shape)\n",
    "\n",
    "    test_patient_genomic = test_df.merge(enc_gene_df, left_on='pids', right_on='pids')\n",
    "    print(test_patient_genomic.shape)\n",
    "    \n",
    "    train = pd.concat([train_patient_genomic,valid_patient_genomic])\n",
    "    test = test_patient_genomic\n",
    "       \n",
    "    train = stratify_patients(train, drug=train_drug, option=option, mmrf_asct=mmrf_asct)\n",
    "    test = stratify_patients(test, drug=test_drug, option=option, mmrf_asct=mmrf_asct)\n",
    "    \n",
    "    print(train.shape)\n",
    "    print(test.shape)\n",
    "    \n",
    "    pfs_thresh = 12\n",
    "    train_y = np.array(train.pop('pfs'))\n",
    "    train_y = np.array([int(pfs>pfs_thresh) for pfs in train_y])\n",
    "    p1 = train.pop('pids')\n",
    "    train_X = np.array(train)\n",
    "\n",
    "    test_y = np.array(test.pop('pfs'))\n",
    "    test_y = np.array([int(pfs>pfs_thresh) for pfs in test_y])\n",
    "    p0 = test.pop('pids')\n",
    "    test_X = np.array(test)\n",
    "    \n",
    "    if grid_search:\n",
    "        clf_fit = predict_PFS_grid_search(clf, model_name, train_X, train_y, test_X, test_y)\n",
    "    else:\n",
    "        clf_fit = predict_PFS(clf, model_name, train_X, train_y, test_X, test_y, training)\n",
    "\n",
    "    return clf_fit"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 3.4 Auto-encoder for RNA-seq"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "def predict_pfs_patient_AE(ae_dict, ae_type, train_df, valid_df, test_df, clf, model_name='',\n",
    "                             train_drug='All', test_drug='Bor', option='', mmrf_asct=None, grid_search=False, training=True):\n",
    "        \n",
    "    print('\\nPredicing PFS with auto-encoded RNAseq...')\n",
    "\n",
    "    print('Using auto-encoder %s'%ae_type)\n",
    "    \n",
    "    train_ae = ae_dict['/mmrf_train_%s'%ae_type]\n",
    "    valid_ae = ae_dict['/mmrf_valid_%s'%ae_type]\n",
    "    test_ae = ae_dict['/mmrf_test_%s'%ae_type]\n",
    "\n",
    "    train_ae_df = train_df.merge(train_ae, left_on='pids', right_on='pids')\n",
    "    print(train_ae_df.shape)\n",
    "\n",
    "    valid_ae_df = valid_df.merge(valid_ae, left_on='pids', right_on='pids')\n",
    "    print(valid_ae_df.shape)\n",
    "\n",
    "    test_ae_df = test_df.merge(test_ae, left_on='pids', right_on='pids')\n",
    "    print(test_ae_df.shape)\n",
    "    \n",
    "    train = pd.concat([train_ae_df,valid_ae_df])\n",
    "    test = test_ae_df\n",
    "\n",
    "    train = stratify_patients(train, drug=train_drug, option=option, mmrf_asct=mmrf_asct)\n",
    "    test = stratify_patients(test, drug=test_drug, option=option, mmrf_asct=mmrf_asct)\n",
    "    \n",
    "    print(train.shape)\n",
    "    print(test.shape)\n",
    "    \n",
    "    pfs_thresh = 12\n",
    "    train_y = np.array(train.pop('pfs'))\n",
    "    train_y = np.array([int(pfs>pfs_thresh) for pfs in train_y])\n",
    "    p1 = train.pop('pids')\n",
    "    train_X = np.array(train)\n",
    "\n",
    "    test_y = np.array(test.pop('pfs'))\n",
    "    test_y = np.array([int(pfs>pfs_thresh) for pfs in test_y])\n",
    "    p0 = test.pop('pids')\n",
    "    test_X = np.array(test)\n",
    "    \n",
    "    if grid_search:\n",
    "        clf_fit = predict_PFS_grid_search(clf, model_name, train_X, train_y, test_X, test_y)\n",
    "    else:\n",
    "        clf_fit = predict_PFS(clf, model_name, train_X, train_y, test_X, test_y, training)\n",
    "\n",
    "    return clf_fit"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4. Run experiments"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Testing transfer learning pipeline from CCLE to MMRF...\n",
      "Reading data files...\n"
     ]
    }
   ],
   "source": [
    "print('Testing transfer learning pipeline from CCLE to MMRF...')\n",
    "    \n",
    "print('Reading data files...')\n",
    "    \n",
    "data_dir = './data/' # Enter the path to your folder here\n",
    "    \n",
    "# have the following files in your data folder\n",
    "ind = 1\n",
    "mmrf_filename = data_dir+'cleaned_mm%d_2mos_pfs_ind.pkl'%(ind)\n",
    "    \n",
    "ccle_exp = pd.read_csv(data_dir+\"CCLE_expression.csv\")\n",
    "sample_info = pd.read_csv(data_dir+\"sample_info.csv\")\n",
    "sanger_dose_response = pd.read_csv(data_dir+\"sanger-dose-response.csv\")\n",
    "mapping = pd.read_csv(data_dir+'Ensembl_HGNC_map_042421.csv')\n",
    "   \n",
    "genomic_fn = data_dir+'MMRF_CoMMpass_IA15a_E74GTF_Salmon_Gene_TPM.txt'\n",
    "\n",
    "mmrf_asct = pd.read_csv('./data/MMRF_ASCT.csv')\n",
    "mmrf_asct = mmrf_asct[['pids','ASCT']]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Part I. Pretraining on CCLE IC50 data\n",
      "Preprocessing CCLE data...\n",
      "Normalizing input features...\n",
      "Running PCA-100 on gene expressions...\n",
      "Training set shape:\n",
      "(560, 100)\n",
      "Test set shape:\n",
      "(241, 100)\n"
     ]
    }
   ],
   "source": [
    "print('\\nPart I. Pretraining on CCLE IC50 data')\n",
    "    \n",
    "bort_data, bort_labels, lenal_data, lenal_labels = preprocess_CCLE(ccle_exp, sample_info, sanger_dose_response, mapping)\n",
    "\n",
    "X_train, y_train, X_test, y_test, pca = split_train_test_CCLE(bort_data, bort_labels, normalize=True, nPCA=100)\n",
    "\n",
    "#linear_regression_IC50(X_train, y_train, X_test, y_test) # linear regression baseline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# pretrain a model with encoding layer size 5\n",
    "#MLP_params = {'niter':30, 'nepoch':10,'batch_size':32,'hidden_layers':[50,20,10],'lr':1e-3, 'enc_layer':5}\n",
    "#pretrain_model = pretrain_MLP_model(X_train, y_train, X_test, y_test, MLP_params)\n",
    "#pretrain_model.save(\"./pretrain_model_051921\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [],
   "source": [
    "# load pretrained NN model for exact reproducibility \n",
    "# model with encoding layer size 5\n",
    "pretrain_model = tf.keras.models.load_model(\"./CCLE_Transfer_Pretrain_Model_latent5\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'name': 'sequential',\n",
       " 'layers': [{'class_name': 'InputLayer',\n",
       "   'config': {'batch_input_shape': (None, 100),\n",
       "    'dtype': 'float32',\n",
       "    'sparse': False,\n",
       "    'ragged': False,\n",
       "    'name': 'dense1_input'}},\n",
       "  {'class_name': 'Dense',\n",
       "   'config': {'name': 'dense1',\n",
       "    'trainable': True,\n",
       "    'batch_input_shape': (None, 100),\n",
       "    'dtype': 'float32',\n",
       "    'units': 45,\n",
       "    'activation': 'sigmoid',\n",
       "    'use_bias': True,\n",
       "    'kernel_initializer': {'class_name': 'GlorotUniform',\n",
       "     'config': {'seed': None}},\n",
       "    'bias_initializer': {'class_name': 'Zeros', 'config': {}},\n",
       "    'kernel_regularizer': None,\n",
       "    'bias_regularizer': None,\n",
       "    'activity_regularizer': None,\n",
       "    'kernel_constraint': None,\n",
       "    'bias_constraint': None}},\n",
       "  {'class_name': 'Dropout',\n",
       "   'config': {'name': 'dropout1',\n",
       "    'trainable': True,\n",
       "    'dtype': 'float32',\n",
       "    'rate': 0.4,\n",
       "    'noise_shape': None,\n",
       "    'seed': None}},\n",
       "  {'class_name': 'Dense',\n",
       "   'config': {'name': 'dense2',\n",
       "    'trainable': True,\n",
       "    'dtype': 'float32',\n",
       "    'units': 45,\n",
       "    'activation': 'sigmoid',\n",
       "    'use_bias': True,\n",
       "    'kernel_initializer': {'class_name': 'GlorotUniform',\n",
       "     'config': {'seed': None}},\n",
       "    'bias_initializer': {'class_name': 'Zeros', 'config': {}},\n",
       "    'kernel_regularizer': None,\n",
       "    'bias_regularizer': None,\n",
       "    'activity_regularizer': None,\n",
       "    'kernel_constraint': None,\n",
       "    'bias_constraint': None}},\n",
       "  {'class_name': 'Dropout',\n",
       "   'config': {'name': 'dropout2',\n",
       "    'trainable': True,\n",
       "    'dtype': 'float32',\n",
       "    'rate': 0.4,\n",
       "    'noise_shape': None,\n",
       "    'seed': None}},\n",
       "  {'class_name': 'Dense',\n",
       "   'config': {'name': 'encoder',\n",
       "    'trainable': True,\n",
       "    'dtype': 'float32',\n",
       "    'units': 5,\n",
       "    'activation': 'linear',\n",
       "    'use_bias': True,\n",
       "    'kernel_initializer': {'class_name': 'GlorotUniform',\n",
       "     'config': {'seed': None}},\n",
       "    'bias_initializer': {'class_name': 'Zeros', 'config': {}},\n",
       "    'kernel_regularizer': None,\n",
       "    'bias_regularizer': None,\n",
       "    'activity_regularizer': None,\n",
       "    'kernel_constraint': None,\n",
       "    'bias_constraint': None}},\n",
       "  {'class_name': 'Dense',\n",
       "   'config': {'name': 'dense',\n",
       "    'trainable': True,\n",
       "    'dtype': 'float32',\n",
       "    'units': 1,\n",
       "    'activation': 'linear',\n",
       "    'use_bias': True,\n",
       "    'kernel_initializer': {'class_name': 'GlorotUniform',\n",
       "     'config': {'seed': None}},\n",
       "    'bias_initializer': {'class_name': 'Zeros', 'config': {}},\n",
       "    'kernel_regularizer': None,\n",
       "    'bias_regularizer': None,\n",
       "    'activity_regularizer': None,\n",
       "    'kernel_constraint': None,\n",
       "    'bias_constraint': None}}]}"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pretrain_model.get_config()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# pretrain a second model with encoding layer size 10\n",
    "# IMPORTANT: remove the training records in \"untitled project\" before running this cell\n",
    "#MLP_params = {'niter':30, 'nepoch':10,'batch_size':32,'hidden_layers':[50,20,10],'lr':1e-3, 'enc_layer':10}\n",
    "#pretrain_model2 = pretrain_MLP_model(X_train, y_train, X_test, y_test, MLP_params)\n",
    "#pretrain_model2.save(\"./pretrain_model2_051921\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 60,
   "metadata": {},
   "outputs": [],
   "source": [
    "# load pretrained NN model for exact reproducibility\n",
    "# model with encoding layer size 5\n",
    "pretrain_model2 = tf.keras.models.load_model(\"./CCLE_Transfer_Pretrain_Model_latent10\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'name': 'sequential',\n",
       " 'layers': [{'class_name': 'InputLayer',\n",
       "   'config': {'batch_input_shape': (None, 100),\n",
       "    'dtype': 'float32',\n",
       "    'sparse': False,\n",
       "    'ragged': False,\n",
       "    'name': 'dense1_input'}},\n",
       "  {'class_name': 'Dense',\n",
       "   'config': {'name': 'dense1',\n",
       "    'trainable': True,\n",
       "    'batch_input_shape': (None, 100),\n",
       "    'dtype': 'float32',\n",
       "    'units': 100,\n",
       "    'activation': 'sigmoid',\n",
       "    'use_bias': True,\n",
       "    'kernel_initializer': {'class_name': 'GlorotUniform',\n",
       "     'config': {'seed': None}},\n",
       "    'bias_initializer': {'class_name': 'Zeros', 'config': {}},\n",
       "    'kernel_regularizer': None,\n",
       "    'bias_regularizer': None,\n",
       "    'activity_regularizer': None,\n",
       "    'kernel_constraint': None,\n",
       "    'bias_constraint': None}},\n",
       "  {'class_name': 'Dropout',\n",
       "   'config': {'name': 'dropout1',\n",
       "    'trainable': True,\n",
       "    'dtype': 'float32',\n",
       "    'rate': 0.0,\n",
       "    'noise_shape': None,\n",
       "    'seed': None}},\n",
       "  {'class_name': 'Dense',\n",
       "   'config': {'name': 'dense2',\n",
       "    'trainable': True,\n",
       "    'dtype': 'float32',\n",
       "    'units': 100,\n",
       "    'activation': 'sigmoid',\n",
       "    'use_bias': True,\n",
       "    'kernel_initializer': {'class_name': 'GlorotUniform',\n",
       "     'config': {'seed': None}},\n",
       "    'bias_initializer': {'class_name': 'Zeros', 'config': {}},\n",
       "    'kernel_regularizer': None,\n",
       "    'bias_regularizer': None,\n",
       "    'activity_regularizer': None,\n",
       "    'kernel_constraint': None,\n",
       "    'bias_constraint': None}},\n",
       "  {'class_name': 'Dropout',\n",
       "   'config': {'name': 'dropout2',\n",
       "    'trainable': True,\n",
       "    'dtype': 'float32',\n",
       "    'rate': 0.0,\n",
       "    'noise_shape': None,\n",
       "    'seed': None}},\n",
       "  {'class_name': 'Dense',\n",
       "   'config': {'name': 'encoder',\n",
       "    'trainable': True,\n",
       "    'dtype': 'float32',\n",
       "    'units': 10,\n",
       "    'activation': 'linear',\n",
       "    'use_bias': True,\n",
       "    'kernel_initializer': {'class_name': 'GlorotUniform',\n",
       "     'config': {'seed': None}},\n",
       "    'bias_initializer': {'class_name': 'Zeros', 'config': {}},\n",
       "    'kernel_regularizer': None,\n",
       "    'bias_regularizer': None,\n",
       "    'activity_regularizer': None,\n",
       "    'kernel_constraint': None,\n",
       "    'bias_constraint': None}},\n",
       "  {'class_name': 'Dense',\n",
       "   'config': {'name': 'dense',\n",
       "    'trainable': True,\n",
       "    'dtype': 'float32',\n",
       "    'units': 1,\n",
       "    'activation': 'linear',\n",
       "    'use_bias': True,\n",
       "    'kernel_initializer': {'class_name': 'GlorotUniform',\n",
       "     'config': {'seed': None}},\n",
       "    'bias_initializer': {'class_name': 'Zeros', 'config': {}},\n",
       "    'kernel_regularizer': None,\n",
       "    'bias_regularizer': None,\n",
       "    'activity_regularizer': None,\n",
       "    'kernel_constraint': None,\n",
       "    'bias_constraint': None}}]}"
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pretrain_model2.get_config()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Part II. Preprocess MMRF genomic data\n",
      "\n",
      "Preprocess patient genomic data in MMRF...\n",
      "Reading ./data/MMRF_CoMMpass_IA15a_E74GTF_Salmon_Gene_TPM.txt\n",
      "Raw data shape:\n",
      "(57997, 922)\n",
      "0 out of 57997 genes have missing data \n",
      "Using raw genomic data...\n",
      "Created new dataframe...\n"
     ]
    }
   ],
   "source": [
    "print('\\nPart II. Preprocess MMRF genomic data')\n",
    "\n",
    "# run PCA on MMRF genomic data and use IC50 model in Part I to encode\n",
    "gene_df, pids = preprocess_MMRF_genomic(genomic_fn, bort_data.columns, mapping)\n",
    "pca_gene = pca.transform(gene_df)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(774, 5)"
      ]
     },
     "execution_count": 27,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# encoded genes for the first CCLE IC50 model\n",
    "layer_name = 'encoder'\n",
    "encoder = tf.keras.Model(inputs=pretrain_model.input,\n",
    "                         outputs=pretrain_model.get_layer(layer_name).output)\n",
    "enc_gene = encoder.predict(pca_gene)\n",
    "    \n",
    "enc_gene_df = pd.DataFrame(data=enc_gene, columns=['Enc%d'%(i+1) for i in range(enc_gene.shape[1])])\n",
    "enc_gene_df['pids'] = pids\n",
    "\n",
    "enc_gene.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 55,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[-0.0362463 , -1.2298235 , -0.9147611 , -1.3867134 ,  0.3141091 ],\n",
       "       [ 0.07769116, -1.2051928 , -0.75010157, -1.0937194 ,  0.28480598],\n",
       "       [ 0.03127058, -1.1820287 , -0.7730207 , -0.9916511 ,  0.44245133],\n",
       "       [ 0.20903537, -1.188896  , -0.7983977 , -1.1671392 ,  0.2296299 ],\n",
       "       [-0.02217743, -1.2201385 , -0.79128355, -0.99069536,  0.39180973],\n",
       "       [-0.06209853, -1.2417276 , -0.77394193, -1.0590153 ,  0.3051305 ],\n",
       "       [ 0.03451161, -1.0500154 , -0.79868567, -1.0617337 ,  0.5624187 ],\n",
       "       [ 0.09678172, -1.1868719 , -0.7550749 , -0.9623866 ,  0.31102258],\n",
       "       [ 0.03811919, -1.2442173 , -0.8186367 , -0.8882976 ,  0.38186038],\n",
       "       [ 0.06550027, -1.1912589 , -0.95859486, -1.1243361 ,  0.30956605]],\n",
       "      dtype=float32)"
      ]
     },
     "execution_count": 55,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "enc_gene[:10,:]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(774, 10)"
      ]
     },
     "execution_count": 28,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# encoded genes for the second CCLE IC50 model\n",
    "layer_name = 'encoder'\n",
    "encoder2 = tf.keras.Model(inputs=pretrain_model2.input,\n",
    "                         outputs=pretrain_model2.get_layer(layer_name).output)\n",
    "enc_gene2 = encoder2.predict(pca_gene)\n",
    "\n",
    "enc_gene_df2 = pd.DataFrame(data=enc_gene2, columns=['Enc%d'%(i+1) for i in range(enc_gene2.shape[1])])\n",
    "enc_gene_df2['pids'] = pids\n",
    "\n",
    "enc_gene2.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Part III. Predict PFS on MMRF patient data\n",
      "loading from: ./data/cleaned_mm1_2mos_pfs_ind.pkl\n",
      "pids\n",
      "(494,)\n",
      "x\n",
      "(494, 33, 16)\n",
      "m\n",
      "(494, 33, 16)\n",
      "feature_names_x\n",
      "(16,)\n",
      "['cbc_abs_neut' 'chem_albumin' 'chem_bun' 'chem_calcium' 'chem_creatinine'\n",
      " 'chem_glucose' 'cbc_hemoglobin' 'serum_kappa' 'serum_m_protein'\n",
      " 'cbc_platelet' 'chem_totprot' 'cbc_wbc' 'serum_iga' 'serum_igg'\n",
      " 'serum_igm' 'serum_lambda']\n",
      "ys_seq\n",
      "(494, 1)\n",
      "ce\n",
      "(494, 1)\n",
      "feature_names_y\n",
      "(1,)\n",
      "['progression free survival (all)']\n",
      "b\n",
      "(494, 16)\n",
      "feature_names\n",
      "(16,)\n",
      "Index(['iss', 'age', 'gender', 'ecog', 'serum_beta2_microglobulin', 'PC1',\n",
      "       'PC2', 'PC3', 'PC4', 'PC5', 'heavy_chain', 'igg_type', 'iga_type',\n",
      "       'igm_type', 'kappa_type', 'lambda_type'],\n",
      "      dtype='object')\n",
      "a\n",
      "(494, 33, 9)\n",
      "m_a\n",
      "(494, 33, 6)\n",
      "feature_names_a\n",
      "(9,)\n",
      "['local_clock' 'Bor' 'Car' 'Cyc' 'Dex' 'Len' 'line1' 'line2' 'line3plus']\n",
      "\n",
      "Preprocess patient data in MMRF\n",
      "Keep first 1 clinical visits.\n",
      "(461, 27)\n",
      "Keep first 1 clinical visits.\n",
      "(161, 27)\n",
      "Keep first 1 clinical visits.\n",
      "(214, 27)\n",
      "Index(['pids', 'pfs', 'iss', 'age', 'gender', 'ecog',\n",
      "       'serum_beta2_microglobulin', 'PC1', 'PC2', 'PC3', 'PC4', 'PC5',\n",
      "       'heavy_chain', 'igg_type', 'iga_type', 'igm_type', 'kappa_type',\n",
      "       'lambda_type', 'local_clock1', 'Bor1', 'Car1', 'Cyc1', 'Dex1', 'Len1',\n",
      "       'line11', 'line21', 'line3plus1'],\n",
      "      dtype='object')\n"
     ]
    }
   ],
   "source": [
    "#data_dir = './data/'\n",
    "#ind = 1\n",
    "#mmrf_filename = data_dir+'cleaned_mm%d_2mos_pfs_ind.pkl'%(ind)\n",
    "\n",
    "print('\\nPart III. Predict PFS on MMRF patient data')\n",
    "train, test, valid = MMRF_utils_V2.get_train_test_valid(mmrf_filename, ind, show_features=True)\n",
    "\n",
    "print('\\nPreprocess patient data in MMRF')\n",
    "# include the first n clinical visit data for each patient (default 0)\n",
    "train_df = MMRF_utils_V2.preprocess_patient_data(train, num_clin_visits=1)  \n",
    "print(train_df.shape)\n",
    "\n",
    "valid_df = MMRF_utils_V2.preprocess_patient_data(valid, num_clin_visits=1)  \n",
    "print(valid_df.shape)\n",
    "\n",
    "test_df = MMRF_utils_V2.preprocess_patient_data(test, num_clin_visits=1)  \n",
    "print(test_df.shape)\n",
    "\n",
    "print(train_df.columns)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 4.1 Only use patient data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Comparing prediction accuracy of different models... \n",
      "\n",
      "\n",
      " Part I \n",
      "\n",
      "Predicting PFS using patient clinical data only...\n",
      "Predicting patients with (Bor, Bor) as first line treatment...\n",
      "(526, 27)\n",
      "(170, 27)\n",
      "Commencing grid search for RF\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\qingy\\Anaconda3\\envs\\data_v2\\lib\\site-packages\\sklearn\\model_selection\\_search.py:289: UserWarning: The total space of parameters 6 is smaller than n_iter=10. Running 6 iterations. For exhaustive searches, use GridSearchCV.\n",
      "  % (grid_size, self.n_iter, grid_size), UserWarning)\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Random search complete, fitting model RF...\n",
      "Predicting on test set...\n",
      "Accuracy: 0.59\n",
      "Precision: 0.60\n",
      "Recall: 0.86\n",
      "F1 score: 0.71\n",
      "AUC: 0.58\n",
      "\n"
     ]
    }
   ],
   "source": [
    "print('\\nComparing prediction accuracy of different models... \\n')\n",
    "    \n",
    "# Model 1: patient data only\n",
    "print('\\n Part I \\n')\n",
    "clf1 = RandomForestClassifier(n_estimators=1000, max_depth=10, random_state=0)    \n",
    "clf1 = predict_pfs_patient_only(train_df,valid_df,test_df, clf1, model_name='RF',\n",
    "                             train_drug='Bor', test_drug='Bor', option='', mmrf_asct=None, grid_search=True, training=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Predicting patients with ASCT\n",
      "Predicting PFS using patient clinical data only...\n",
      "Predicting patients with (Bor, Bor) as first line treatment...\n",
      "(328, 27)\n",
      "(96, 27)\n",
      "Predicting on test set...\n",
      "Accuracy: 0.67\n",
      "Precision: 0.67\n",
      "Recall: 0.93\n",
      "F1 score: 0.78\n",
      "AUC: 0.54\n",
      "\n",
      "Predicting patients without ASCT\n",
      "Predicting PFS using patient clinical data only...\n",
      "Predicting patients with (Bor, Bor) as first line treatment...\n",
      "(198, 27)\n",
      "(74, 27)\n",
      "Predicting on test set...\n",
      "Accuracy: 0.50\n",
      "Precision: 0.49\n",
      "Recall: 0.75\n",
      "F1 score: 0.59\n",
      "AUC: 0.57\n",
      "\n"
     ]
    }
   ],
   "source": [
    "print('Predicting patients with ASCT')\n",
    "_ = predict_pfs_patient_only(train_df,valid_df,test_df, clf1, model_name='RF',\n",
    "                             train_drug='Bor', test_drug='Bor', option='ASCT', mmrf_asct=mmrf_asct, training=False)\n",
    "\n",
    "print('Predicting patients without ASCT')\n",
    "_ = predict_pfs_patient_only(train_df,valid_df,test_df, clf1, model_name='RF',\n",
    "                             train_drug='Bor', test_drug='Bor', option='NON-ASCT', mmrf_asct=mmrf_asct, training=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Best model for Part I:\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "{'bootstrap': True,\n",
       " 'ccp_alpha': 0.0,\n",
       " 'class_weight': None,\n",
       " 'criterion': 'gini',\n",
       " 'max_depth': 10,\n",
       " 'max_features': 'auto',\n",
       " 'max_leaf_nodes': None,\n",
       " 'max_samples': None,\n",
       " 'min_impurity_decrease': 0.0,\n",
       " 'min_impurity_split': None,\n",
       " 'min_samples_leaf': 1,\n",
       " 'min_samples_split': 2,\n",
       " 'min_weight_fraction_leaf': 0.0,\n",
       " 'n_estimators': 5000,\n",
       " 'n_jobs': None,\n",
       " 'oob_score': False,\n",
       " 'random_state': 0,\n",
       " 'verbose': 0,\n",
       " 'warm_start': False}"
      ]
     },
     "execution_count": 32,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "print('Best model for Part I:')\n",
    "clf1.get_params()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 4.2 Use patient data and PCA on raw RNA-seq"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      " Part II \n",
      "\n",
      "Predicting PFS using patient clinical data and PCA-50 on RNAseq...\n",
      "Predicting patients with (Bor, Bor) as first line treatment...\n",
      "\n",
      "Preprocess patient genomic data in MMRF...\n",
      "Reading ./data/MMRF_CoMMpass_IA15a_E74GTF_Salmon_Gene_TPM.txt\n",
      "Raw data shape:\n",
      "(57997, 922)\n",
      "0 out of 57997 genes have missing data \n",
      "Running PCA with 50 principal components...\n",
      "Finished PCA!\n",
      "Created new dataframe...\n",
      "(774, 51)\n",
      "\n",
      "Merging patient data with genomic data...\n",
      "(313, 77)\n",
      "(97, 77)\n",
      "(142, 77)\n",
      "(345, 77)\n",
      "(112, 77)\n",
      "Commencing grid search for RF\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\qingy\\Anaconda3\\envs\\data_v2\\lib\\site-packages\\sklearn\\model_selection\\_search.py:289: UserWarning: The total space of parameters 6 is smaller than n_iter=10. Running 6 iterations. For exhaustive searches, use GridSearchCV.\n",
      "  % (grid_size, self.n_iter, grid_size), UserWarning)\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Random search complete, fitting model RF...\n",
      "Predicting on test set...\n",
      "Accuracy: 0.59\n",
      "Precision: 0.58\n",
      "Recall: 0.95\n",
      "F1 score: 0.72\n",
      "AUC: 0.65\n",
      "\n"
     ]
    }
   ],
   "source": [
    "# Model 2: patient data with PCA on MMRF RNAseq\n",
    "print('\\n Part II \\n')\n",
    "\n",
    "nPCA = 50\n",
    "genomic_fn = data_dir+'MMRF_CoMMpass_IA15a_E74GTF_Salmon_Gene_TPM.txt'\n",
    "\n",
    "clf2 = RandomForestClassifier(n_estimators=1000, max_depth=10, random_state=0)\n",
    "clf2 = predict_pfs_patient_PCA(genomic_fn, nPCA, train_df, valid_df, test_df, clf2, model_name='RF',\n",
    "                            train_drug='Bor', test_drug='Bor', option='', mmrf_asct=None, grid_search=True, training=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Predicting PFS using patient clinical data and PCA-50 on RNAseq...\n",
      "Predicting patients with (Bor, Bor) as first line treatment...\n",
      "\n",
      "Preprocess patient genomic data in MMRF...\n",
      "Reading ./data/MMRF_CoMMpass_IA15a_E74GTF_Salmon_Gene_TPM.txt\n",
      "Raw data shape:\n",
      "(57997, 922)\n",
      "0 out of 57997 genes have missing data \n",
      "Running PCA with 50 principal components...\n",
      "Finished PCA!\n",
      "Created new dataframe...\n",
      "(774, 51)\n",
      "\n",
      "Merging patient data with genomic data...\n",
      "(313, 77)\n",
      "(97, 77)\n",
      "(142, 77)\n",
      "(218, 77)\n",
      "(58, 77)\n",
      "Predicting on test set...\n",
      "Accuracy: 0.67\n",
      "Precision: 0.68\n",
      "Recall: 0.95\n",
      "F1 score: 0.79\n",
      "AUC: 0.65\n",
      "\n"
     ]
    }
   ],
   "source": [
    "_ = predict_pfs_patient_PCA(genomic_fn, nPCA, train_df, valid_df, test_df, clf2, model_name='RF',\n",
    "                            train_drug='Bor', test_drug='Bor', option='ASCT', mmrf_asct=mmrf_asct, grid_search=False, training=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Predicting PFS using patient clinical data and PCA-50 on RNAseq...\n",
      "Predicting patients with (Bor, Bor) as first line treatment...\n",
      "\n",
      "Preprocess patient genomic data in MMRF...\n",
      "Reading ./data/MMRF_CoMMpass_IA15a_E74GTF_Salmon_Gene_TPM.txt\n",
      "Raw data shape:\n",
      "(57997, 922)\n",
      "0 out of 57997 genes have missing data \n",
      "Running PCA with 50 principal components...\n",
      "Finished PCA!\n",
      "Created new dataframe...\n",
      "(774, 51)\n",
      "\n",
      "Merging patient data with genomic data...\n",
      "(313, 77)\n",
      "(97, 77)\n",
      "(142, 77)\n",
      "(127, 77)\n",
      "(54, 77)\n",
      "Predicting on test set...\n",
      "Accuracy: 0.50\n",
      "Precision: 0.47\n",
      "Recall: 0.96\n",
      "F1 score: 0.63\n",
      "AUC: 0.66\n",
      "\n"
     ]
    }
   ],
   "source": [
    "_ = predict_pfs_patient_PCA(genomic_fn, nPCA, train_df, valid_df, test_df, clf2, model_name='RF',\n",
    "                            train_drug='Bor', test_drug='Bor', option='NON-ASCT', mmrf_asct=mmrf_asct, grid_search=False, training=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Best model for Part II:\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "{'bootstrap': True,\n",
       " 'ccp_alpha': 0.0,\n",
       " 'class_weight': None,\n",
       " 'criterion': 'gini',\n",
       " 'max_depth': 20,\n",
       " 'max_features': 'auto',\n",
       " 'max_leaf_nodes': None,\n",
       " 'max_samples': None,\n",
       " 'min_impurity_decrease': 0.0,\n",
       " 'min_impurity_split': None,\n",
       " 'min_samples_leaf': 1,\n",
       " 'min_samples_split': 2,\n",
       " 'min_weight_fraction_leaf': 0.0,\n",
       " 'n_estimators': 1000,\n",
       " 'n_jobs': None,\n",
       " 'oob_score': False,\n",
       " 'random_state': 0,\n",
       " 'verbose': 0,\n",
       " 'warm_start': False}"
      ]
     },
     "execution_count": 36,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "print('Best model for Part II:')\n",
    "\n",
    "clf2.get_params()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 4.3 Transfer from CCLE IC50 Predictions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      " Part III \n",
      "\n",
      "First IC50 model\n",
      "\n",
      "\n",
      "Predicing PFS with transfer model from CCLE...\n",
      "\n",
      "Merging patient data with genomic data...\n",
      "(313, 32)\n",
      "(97, 32)\n",
      "(142, 32)\n",
      "(345, 32)\n",
      "(112, 32)\n",
      "Commencing grid search for RF\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\qingy\\Anaconda3\\envs\\data_v2\\lib\\site-packages\\sklearn\\model_selection\\_search.py:289: UserWarning: The total space of parameters 6 is smaller than n_iter=10. Running 6 iterations. For exhaustive searches, use GridSearchCV.\n",
      "  % (grid_size, self.n_iter, grid_size), UserWarning)\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Random search complete, fitting model RF...\n",
      "Predicting on test set...\n",
      "Accuracy: 0.63\n",
      "Precision: 0.61\n",
      "Recall: 0.92\n",
      "F1 score: 0.74\n",
      "AUC: 0.68\n",
      "\n"
     ]
    }
   ],
   "source": [
    "print('\\n Part III \\n')\n",
    "print('First IC50 model\\n')\n",
    "    \n",
    "clf3 = RandomForestClassifier(n_estimators=1000, max_depth=10, random_state=0)\n",
    "clf3 = predict_pfs_patient_CCLE(enc_gene_df, train_df, valid_df, test_df, clf3, model_name='RF',\n",
    "                            train_drug='Bor', test_drug='Bor', option='', mmrf_asct=None, grid_search=True, training=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Predicing PFS with transfer model from CCLE...\n",
      "\n",
      "Merging patient data with genomic data...\n",
      "(313, 32)\n",
      "(97, 32)\n",
      "(142, 32)\n",
      "(218, 32)\n",
      "(58, 32)\n",
      "Predicting on test set...\n",
      "Accuracy: 0.72\n",
      "Precision: 0.72\n",
      "Recall: 0.95\n",
      "F1 score: 0.82\n",
      "AUC: 0.68\n",
      "\n",
      "\n",
      "Predicing PFS with transfer model from CCLE...\n",
      "\n",
      "Merging patient data with genomic data...\n",
      "(313, 32)\n",
      "(97, 32)\n",
      "(142, 32)\n",
      "(127, 32)\n",
      "(54, 32)\n",
      "Predicting on test set...\n",
      "Accuracy: 0.54\n",
      "Precision: 0.49\n",
      "Recall: 0.88\n",
      "F1 score: 0.63\n",
      "AUC: 0.66\n",
      "\n"
     ]
    }
   ],
   "source": [
    "_ = predict_pfs_patient_CCLE(enc_gene_df, train_df, valid_df, test_df, clf3, model_name='RF',\n",
    "                            train_drug='Bor', test_drug='Bor', option='ASCT', mmrf_asct=mmrf_asct, grid_search=False, training=False)\n",
    "\n",
    "_ = predict_pfs_patient_CCLE(enc_gene_df, train_df, valid_df, test_df, clf3, model_name='RF',\n",
    "                            train_drug='Bor', test_drug='Bor', option='NON-ASCT', mmrf_asct=mmrf_asct, grid_search=False, training=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'bootstrap': True,\n",
       " 'ccp_alpha': 0.0,\n",
       " 'class_weight': None,\n",
       " 'criterion': 'gini',\n",
       " 'max_depth': 20,\n",
       " 'max_features': 'auto',\n",
       " 'max_leaf_nodes': None,\n",
       " 'max_samples': None,\n",
       " 'min_impurity_decrease': 0.0,\n",
       " 'min_impurity_split': None,\n",
       " 'min_samples_leaf': 1,\n",
       " 'min_samples_split': 2,\n",
       " 'min_weight_fraction_leaf': 0.0,\n",
       " 'n_estimators': 1000,\n",
       " 'n_jobs': None,\n",
       " 'oob_score': False,\n",
       " 'random_state': 0,\n",
       " 'verbose': 0,\n",
       " 'warm_start': False}"
      ]
     },
     "execution_count": 39,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "clf3.get_params()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      " Part III \n",
      "\n",
      "Second IC50 model\n",
      "\n",
      "\n",
      "Predicing PFS with transfer model from CCLE...\n",
      "\n",
      "Merging patient data with genomic data...\n",
      "(313, 37)\n",
      "(97, 37)\n",
      "(142, 37)\n",
      "(345, 37)\n",
      "(112, 37)\n",
      "Commencing grid search for RF\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\qingy\\Anaconda3\\envs\\data_v2\\lib\\site-packages\\sklearn\\model_selection\\_search.py:289: UserWarning: The total space of parameters 6 is smaller than n_iter=10. Running 6 iterations. For exhaustive searches, use GridSearchCV.\n",
      "  % (grid_size, self.n_iter, grid_size), UserWarning)\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Random search complete, fitting model RF...\n",
      "Predicting on test set...\n",
      "Accuracy: 0.62\n",
      "Precision: 0.61\n",
      "Recall: 0.85\n",
      "F1 score: 0.71\n",
      "AUC: 0.65\n",
      "\n"
     ]
    }
   ],
   "source": [
    "print('\\n Part III \\n')\n",
    "print('Second IC50 model\\n')\n",
    "clf3_2 = RandomForestClassifier(n_estimators=1000, max_depth=10, random_state=0)\n",
    "clf3_2 = predict_pfs_patient_CCLE(enc_gene_df2, train_df, valid_df, test_df, clf3_2, model_name='RF',\n",
    "                            train_drug='Bor', test_drug='Bor', option='', mmrf_asct=None, grid_search=True, training=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Predicing PFS with transfer model from CCLE...\n",
      "\n",
      "Merging patient data with genomic data...\n",
      "(313, 37)\n",
      "(97, 37)\n",
      "(142, 37)\n",
      "(218, 37)\n",
      "(58, 37)\n",
      "Predicting on test set...\n",
      "Accuracy: 0.71\n",
      "Precision: 0.72\n",
      "Recall: 0.89\n",
      "F1 score: 0.80\n",
      "AUC: 0.64\n",
      "\n",
      "\n",
      "Predicing PFS with transfer model from CCLE...\n",
      "\n",
      "Merging patient data with genomic data...\n",
      "(313, 37)\n",
      "(97, 37)\n",
      "(142, 37)\n",
      "(127, 37)\n",
      "(54, 37)\n",
      "Predicting on test set...\n",
      "Accuracy: 0.52\n",
      "Precision: 0.47\n",
      "Recall: 0.79\n",
      "F1 score: 0.59\n",
      "AUC: 0.63\n",
      "\n"
     ]
    }
   ],
   "source": [
    "_ = predict_pfs_patient_CCLE(enc_gene_df2, train_df, valid_df, test_df, clf3_2, model_name='RF',\n",
    "                            train_drug='Bor', test_drug='Bor', option='ASCT', mmrf_asct=mmrf_asct, grid_search=False, training=False)\n",
    "\n",
    "_ = predict_pfs_patient_CCLE(enc_gene_df2, train_df, valid_df, test_df, clf3_2, model_name='RF',\n",
    "                            train_drug='Bor', test_drug='Bor', option='NON-ASCT', mmrf_asct=mmrf_asct, grid_search=False, training=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'bootstrap': True,\n",
       " 'ccp_alpha': 0.0,\n",
       " 'class_weight': None,\n",
       " 'criterion': 'gini',\n",
       " 'max_depth': 20,\n",
       " 'max_features': 'auto',\n",
       " 'max_leaf_nodes': None,\n",
       " 'max_samples': None,\n",
       " 'min_impurity_decrease': 0.0,\n",
       " 'min_impurity_split': None,\n",
       " 'min_samples_leaf': 1,\n",
       " 'min_samples_split': 2,\n",
       " 'min_weight_fraction_leaf': 0.0,\n",
       " 'n_estimators': 5000,\n",
       " 'n_jobs': None,\n",
       " 'oob_score': False,\n",
       " 'random_state': 0,\n",
       " 'verbose': 0,\n",
       " 'warm_start': False}"
      ]
     },
     "execution_count": 42,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "clf3_2.get_params()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 4.4 Transfer from auto-encoded RNA-seq"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Part IV. Predict PFS with auto-encoded RNA-seq\n",
      "\n",
      "Predicing PFS with auto-encoded RNAseq...\n",
      "Using auto-encoder combinedae\n",
      "(313, 91)\n",
      "(97, 91)\n",
      "(142, 91)\n",
      "(345, 91)\n",
      "(112, 91)\n",
      "Commencing grid search for RF\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\qingy\\Anaconda3\\envs\\data_v2\\lib\\site-packages\\sklearn\\model_selection\\_search.py:289: UserWarning: The total space of parameters 6 is smaller than n_iter=10. Running 6 iterations. For exhaustive searches, use GridSearchCV.\n",
      "  % (grid_size, self.n_iter, grid_size), UserWarning)\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Random search complete, fitting model RF...\n",
      "Predicting on test set...\n",
      "Accuracy: 0.57\n",
      "Precision: 0.57\n",
      "Recall: 0.87\n",
      "F1 score: 0.69\n",
      "AUC: 0.62\n",
      "\n"
     ]
    }
   ],
   "source": [
    "# Model 4: Auto-encoder for CCLE/MMRF\n",
    "print('\\nPart IV. Predict PFS with auto-encoded RNA-seq')\n",
    "\n",
    "ae_fn = data_dir + 'normalized_embeddings_64.pkl' # auto-encoder of RNA\n",
    "ae_dict = pickle.load(open(ae_fn, \"rb\" ))\n",
    "ae_type = 'combinedae' # other types: ['mmrfae', 'ccleae', 'combinedae']\n",
    "\n",
    "clf4 = RandomForestClassifier(n_estimators=1000, max_depth=50, random_state=0)\n",
    "\n",
    "clf4 = predict_pfs_patient_AE(ae_dict, ae_type, train_df, valid_df, test_df, clf4, model_name='RF',\n",
    "                        train_drug='Bor', test_drug='Bor', option='', mmrf_asct=None, grid_search=True, training=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 44,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Predicing PFS with auto-encoded RNAseq...\n",
      "Using auto-encoder combinedae\n",
      "(313, 91)\n",
      "(97, 91)\n",
      "(142, 91)\n",
      "(218, 91)\n",
      "(58, 91)\n",
      "Predicting on test set...\n",
      "Accuracy: 0.66\n",
      "Precision: 0.69\n",
      "Recall: 0.87\n",
      "F1 score: 0.77\n",
      "AUC: 0.63\n",
      "\n"
     ]
    }
   ],
   "source": [
    "_ = predict_pfs_patient_AE(ae_dict, ae_type, train_df, valid_df, test_df, clf4, model_name='RF',\n",
    "                        train_drug='Bor', test_drug='Bor', option='ASCT', mmrf_asct=mmrf_asct, grid_search=False, training=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Predicing PFS with auto-encoded RNAseq...\n",
      "Using auto-encoder combinedae\n",
      "(313, 91)\n",
      "(97, 91)\n",
      "(142, 91)\n",
      "(127, 91)\n",
      "(54, 91)\n",
      "Predicting on test set...\n",
      "Accuracy: 0.48\n",
      "Precision: 0.46\n",
      "Recall: 0.88\n",
      "F1 score: 0.60\n",
      "AUC: 0.60\n",
      "\n"
     ]
    }
   ],
   "source": [
    "_ = predict_pfs_patient_AE(ae_dict, ae_type, train_df, valid_df, test_df, clf4, model_name='RF',\n",
    "                        train_drug='Bor', test_drug='Bor', option='NON-ASCT', mmrf_asct=mmrf_asct, grid_search=False, training=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Best model for Part IV (Latent 64):\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "{'bootstrap': True,\n",
       " 'ccp_alpha': 0.0,\n",
       " 'class_weight': None,\n",
       " 'criterion': 'gini',\n",
       " 'max_depth': 20,\n",
       " 'max_features': 'auto',\n",
       " 'max_leaf_nodes': None,\n",
       " 'max_samples': None,\n",
       " 'min_impurity_decrease': 0.0,\n",
       " 'min_impurity_split': None,\n",
       " 'min_samples_leaf': 1,\n",
       " 'min_samples_split': 2,\n",
       " 'min_weight_fraction_leaf': 0.0,\n",
       " 'n_estimators': 5000,\n",
       " 'n_jobs': None,\n",
       " 'oob_score': False,\n",
       " 'random_state': 0,\n",
       " 'verbose': 0,\n",
       " 'warm_start': False}"
      ]
     },
     "execution_count": 46,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "print('Best model for Part IV (Latent 64):')\n",
    "\n",
    "clf4.get_params()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      " Part IV \n",
      "\n",
      "\n",
      "Predicing PFS with auto-encoded RNAseq...\n",
      "Using auto-encoder combinedae\n",
      "(313, 155)\n",
      "(97, 155)\n",
      "(142, 155)\n",
      "(345, 155)\n",
      "(112, 155)\n",
      "Commencing grid search for RF\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\qingy\\Anaconda3\\envs\\data_v2\\lib\\site-packages\\sklearn\\model_selection\\_search.py:289: UserWarning: The total space of parameters 6 is smaller than n_iter=10. Running 6 iterations. For exhaustive searches, use GridSearchCV.\n",
      "  % (grid_size, self.n_iter, grid_size), UserWarning)\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Random search complete, fitting model RF...\n",
      "Predicting on test set...\n",
      "Accuracy: 0.59\n",
      "Precision: 0.59\n",
      "Recall: 0.87\n",
      "F1 score: 0.70\n",
      "AUC: 0.65\n",
      "\n"
     ]
    }
   ],
   "source": [
    "# Model 4: Auto-encoder for CCLE/MMRF\n",
    "# Try the AE with latent dimension 128\n",
    "print('\\n Part IV \\n')\n",
    "\n",
    "ae_fn = data_dir + 'normalized_embeddings_128.pkl' # auto-encoder of RNA\n",
    "ae_dict = pickle.load(open(ae_fn, \"rb\" ))\n",
    "ae_type = 'combinedae' # other types: ['mmrfae', 'ccleae', 'combinedae']\n",
    "\n",
    "clf4_2 = RandomForestClassifier(n_estimators=1000, max_depth=50, random_state=0)\n",
    "\n",
    "clf4_2 = predict_pfs_patient_AE(ae_dict, ae_type, train_df, valid_df, test_df, clf4_2, model_name='RF',\n",
    "                        train_drug='Bor', test_drug='Bor', option='', mmrf_asct=None, grid_search=True, training=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 48,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Predicing PFS with auto-encoded RNAseq...\n",
      "Using auto-encoder combinedae\n",
      "(313, 155)\n",
      "(97, 155)\n",
      "(142, 155)\n",
      "(218, 155)\n",
      "(58, 155)\n",
      "Predicting on test set...\n",
      "Accuracy: 0.64\n",
      "Precision: 0.68\n",
      "Recall: 0.84\n",
      "F1 score: 0.75\n",
      "AUC: 0.66\n",
      "\n"
     ]
    }
   ],
   "source": [
    "_ = predict_pfs_patient_AE(ae_dict, ae_type, train_df, valid_df, test_df, clf4_2, model_name='RF',\n",
    "                        train_drug='Bor', test_drug='Bor', option='ASCT', mmrf_asct=mmrf_asct, grid_search=False, training=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 49,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Predicing PFS with auto-encoded RNAseq...\n",
      "Using auto-encoder combinedae\n",
      "(313, 155)\n",
      "(97, 155)\n",
      "(142, 155)\n",
      "(127, 155)\n",
      "(54, 155)\n",
      "Predicting on test set...\n",
      "Accuracy: 0.54\n",
      "Precision: 0.49\n",
      "Recall: 0.92\n",
      "F1 score: 0.64\n",
      "AUC: 0.63\n",
      "\n"
     ]
    }
   ],
   "source": [
    "_ = predict_pfs_patient_AE(ae_dict, ae_type, train_df, valid_df, test_df, clf4_2, model_name='RF',\n",
    "                        train_drug='Bor', test_drug='Bor', option='NON-ASCT', mmrf_asct=mmrf_asct, grid_search=False, training=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 50,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Best model for Part IV (Latent Dim 128):\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "{'bootstrap': True,\n",
       " 'ccp_alpha': 0.0,\n",
       " 'class_weight': None,\n",
       " 'criterion': 'gini',\n",
       " 'max_depth': 10,\n",
       " 'max_features': 'auto',\n",
       " 'max_leaf_nodes': None,\n",
       " 'max_samples': None,\n",
       " 'min_impurity_decrease': 0.0,\n",
       " 'min_impurity_split': None,\n",
       " 'min_samples_leaf': 1,\n",
       " 'min_samples_split': 2,\n",
       " 'min_weight_fraction_leaf': 0.0,\n",
       " 'n_estimators': 5000,\n",
       " 'n_jobs': None,\n",
       " 'oob_score': False,\n",
       " 'random_state': 0,\n",
       " 'verbose': 0,\n",
       " 'warm_start': False}"
      ]
     },
     "execution_count": 50,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "print('Best model for Part IV (Latent Dim 128):')\n",
    "\n",
    "clf4_2.get_params()"
   ]
  }
 ],
 "metadata": {
  "colab": {
   "name": "Task1_CCLE_IC50_V_hypermodel.ipynb",
   "provenance": []
  },
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}