scbert_baselines_LR.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/opt/conda/rpeyser/envs/rna-pretrain/lib/python3.9/site-packages/tqdm/auto.py:22: TqdmWarning: IProgress not found. Please update jupyter and ipywidgets. See https://ipywidgets.readthedocs.io/en/stable/user_install.html\n",
      "  from .autonotebook import tqdm as notebook_tqdm\n"
     ]
    }
   ],
   "source": [
    "# -*- coding: utf-8 -*-\n",
    "import os\n",
    "import gc\n",
    "import argparse\n",
    "import json\n",
    "import random\n",
    "import math\n",
    "import random\n",
    "from functools import reduce\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "from scipy import sparse\n",
    "from sklearn.model_selection import train_test_split, ShuffleSplit, StratifiedShuffleSplit, StratifiedKFold\n",
    "from sklearn.metrics import accuracy_score, f1_score, confusion_matrix, precision_recall_fscore_support, classification_report\n",
    "import torch\n",
    "from torch import nn\n",
    "from torch.optim import Adam, SGD, AdamW\n",
    "from torch.nn import functional as F\n",
    "from torch.optim.lr_scheduler import StepLR, CosineAnnealingWarmRestarts, CyclicLR\n",
    "from torch.utils.data import DataLoader, Dataset\n",
    "from torch.utils.data.distributed import DistributedSampler\n",
    "from torch.nn.parallel import DistributedDataParallel as DDP\n",
    "import torch.distributed as dist\n",
    "from tqdm import tqdm\n",
    "\n",
    "from performer_pytorch import PerformerLM\n",
    "import scanpy as sc\n",
    "import anndata as ad\n",
    "from utils import *\n",
    "from datetime import datetime\n",
    "from time import time\n",
    "import torch.multiprocessing as mp\n",
    "from torch.nn.modules.utils import consume_prefix_in_state_dict_if_present\n",
    "from torch.utils.tensorboard import SummaryWriter\n",
    "\n",
    "from sklearn.model_selection import cross_validate\n",
    "from sklearn.linear_model import LogisticRegression\n",
    "from sklearn.metrics import f1_score"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [],
   "source": [
    "SEED = 2021\n",
    "NREPS = 3\n",
    "SAMPLING_FRACS = [1.0, 0.75, 0.5]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Control sources of randomness\n",
    "torch.manual_seed(SEED)\n",
    "random.seed(SEED)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "AnnData object with n_obs × n_vars = 68450 × 16906\n",
       "    obs: 'TSNE.1', 'TSNE.2', 'celltype', 'n_genes'\n",
       "    uns: 'log1p'"
      ]
     },
     "execution_count": 36,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#read Zheng data\n",
    "zheng_data = sc.read_h5ad(\"/data/rna_rep_learning/scBERT/Zheng68K.h5ad\")\n",
    "zheng_data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {},
   "outputs": [],
   "source": [
    "data = zheng_data.X\n",
    "label = zheng_data.obs.celltype"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "frac 1.0, rep 0\n",
      "Loaded data...\n",
      "c=0.001\n",
      "c=0.01\n",
      "c=0.1\n",
      "c=1\n",
      "{0.001: 0.6885135135135135, 0.01: 0.7837837837837838, 0.1: 0.8014061358655953, 1: 0.7621621621621621}\n",
      "best c=0.001\n",
      "train set accuracy: 0.6993\n"
     ]
    },
    {
     "ename": "NameError",
     "evalue": "name 'y_test' is not defined",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
      "Input \u001b[0;32mIn [45]\u001b[0m, in \u001b[0;36m<cell line: 1>\u001b[0;34m()\u001b[0m\n\u001b[1;32m     30\u001b[0m lr\u001b[38;5;241m.\u001b[39mfit(X_train, y_train)\n\u001b[1;32m     31\u001b[0m \u001b[38;5;28mprint\u001b[39m(\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mtrain set accuracy: \u001b[39m\u001b[38;5;124m\"\u001b[39m \u001b[38;5;241m+\u001b[39m \u001b[38;5;28mstr\u001b[39m(np\u001b[38;5;241m.\u001b[39maround(lr\u001b[38;5;241m.\u001b[39mscore(X_train, y_train), \u001b[38;5;241m4\u001b[39m)))\n\u001b[0;32m---> 32\u001b[0m \u001b[38;5;28mprint\u001b[39m(\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mtest set accuracy: \u001b[39m\u001b[38;5;124m\"\u001b[39m \u001b[38;5;241m+\u001b[39m \u001b[38;5;28mstr\u001b[39m(np\u001b[38;5;241m.\u001b[39maround(lr\u001b[38;5;241m.\u001b[39mscore(X_test, \u001b[43my_test\u001b[49m), \u001b[38;5;241m4\u001b[39m)))\n\u001b[1;32m     33\u001b[0m val_macro_f1 \u001b[38;5;241m=\u001b[39m sklearn\u001b[38;5;241m.\u001b[39mmetrics\u001b[38;5;241m.\u001b[39mf1_score(y_test, lr\u001b[38;5;241m.\u001b[39mpredict(X_test), average\u001b[38;5;241m=\u001b[39m\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mmacro\u001b[39m\u001b[38;5;124m\"\u001b[39m)\n\u001b[1;32m     34\u001b[0m \u001b[38;5;28mprint\u001b[39m(\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mtest set macro F1: \u001b[39m\u001b[38;5;124m\"\u001b[39m \u001b[38;5;241m+\u001b[39m \u001b[38;5;28mstr\u001b[39m(np\u001b[38;5;241m.\u001b[39maround(val_macro_f1, \u001b[38;5;241m4\u001b[39m)))\n",
      "\u001b[0;31mNameError\u001b[0m: name 'y_test' is not defined"
     ]
    }
   ],
   "source": [
    "for k in np.arange(NREPS):\n",
    "    for frac in SAMPLING_FRACS:\n",
    "        print(\"frac {}, rep {}\".format(frac, k))\n",
    "        #downsample training set\n",
    "        sss = StratifiedShuffleSplit(n_splits=1, test_size=0.2, random_state=k) #same train/val set split for each frac in k\n",
    "        for index_train, index_val in sss.split(data, label):\n",
    "            np.random.seed(k)\n",
    "            index_train_small = np.random.choice(index_train, round(index_train.shape[0]*frac), replace=False)\n",
    "            X_train, y_train = data[index_train_small], label[index_train_small]\n",
    "            X_test, y_test = data[index_val], label[index_val]\n",
    "\n",
    "        print(\"Loaded data...\")\n",
    "\n",
    "        #train on train_dataset\n",
    "        \n",
    "        #hyperparameter tune using k-fold val on training data\n",
    "        cv_results = {}\n",
    "        for c in [1e-3, 1e-2, 1e-1, 1]:\n",
    "            print(\"c={}\".format(c))\n",
    "            lr = LogisticRegression(random_state=0, penalty=\"l1\", C=c, solver=\"liblinear\")\n",
    "            res = cross_validate(lr, X_train, y_train, scoring=['accuracy'])\n",
    "            cv_results[c] = np.mean(res['test_accuracy'])\n",
    "        print(cv_results)\n",
    "\n",
    "        #choose best c and calc performance on val_dataset\n",
    "        best_ind = np.argmax(list(cv_results.values()))\n",
    "        c = list(cv_results.keys())[best_ind]\n",
    "        print(\"best c={}\".format(c))\n",
    "        lr = LogisticRegression(random_state=0, penalty=\"l1\", C=c, solver=\"liblinear\")\n",
    "        lr.fit(X_train, y_train)\n",
    "        print(\"train set accuracy: \" + str(np.around(lr.score(X_train, y_train), 4)))\n",
    "        print(\"test set accuracy: \" + str(np.around(lr.score(X_test, y_test), 4)))\n",
    "        val_macro_f1 = sklearn.metrics.f1_score(y_test, lr.predict(X_test), average=\"macro\")\n",
    "        print(\"test set macro F1: \" + str(np.around(val_macro_f1, 4)))\n",
    "        \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 59,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "frac 0.75, rep 0\n",
      "Loaded data...\n",
      "c=0.001\n",
      "c=0.01\n",
      "c=0.1\n",
      "c=1\n",
      "{0.001: 0.6780618456294131, 0.01: 0.7754808862916971, 0.1: 0.7914049184319454, 1: 0.7393474555636718}\n",
      "best c=0.1\n",
      "train set accuracy: 0.9522\n",
      "test set accuracy: 0.7981\n"
     ]
    },
    {
     "ename": "NameError",
     "evalue": "name 'sklearn' is not defined",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
      "Input \u001b[0;32mIn [59]\u001b[0m, in \u001b[0;36m<cell line: 1>\u001b[0;34m()\u001b[0m\n\u001b[1;32m     33\u001b[0m \u001b[38;5;28mprint\u001b[39m(\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mtrain set accuracy: \u001b[39m\u001b[38;5;124m\"\u001b[39m \u001b[38;5;241m+\u001b[39m \u001b[38;5;28mstr\u001b[39m(np\u001b[38;5;241m.\u001b[39maround(lr\u001b[38;5;241m.\u001b[39mscore(X_train, y_train), \u001b[38;5;241m4\u001b[39m)))\n\u001b[1;32m     34\u001b[0m \u001b[38;5;28mprint\u001b[39m(\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mtest set accuracy: \u001b[39m\u001b[38;5;124m\"\u001b[39m \u001b[38;5;241m+\u001b[39m \u001b[38;5;28mstr\u001b[39m(np\u001b[38;5;241m.\u001b[39maround(lr\u001b[38;5;241m.\u001b[39mscore(X_test, y_test), \u001b[38;5;241m4\u001b[39m)))\n\u001b[0;32m---> 35\u001b[0m val_macro_f1 \u001b[38;5;241m=\u001b[39m \u001b[43msklearn\u001b[49m\u001b[38;5;241m.\u001b[39mmetrics\u001b[38;5;241m.\u001b[39mf1_score(y_test, lr\u001b[38;5;241m.\u001b[39mpredict(X_test), average\u001b[38;5;241m=\u001b[39m\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mmacro\u001b[39m\u001b[38;5;124m\"\u001b[39m)\n\u001b[1;32m     36\u001b[0m \u001b[38;5;28mprint\u001b[39m(\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mtest set macro F1: \u001b[39m\u001b[38;5;124m\"\u001b[39m \u001b[38;5;241m+\u001b[39m \u001b[38;5;28mstr\u001b[39m(np\u001b[38;5;241m.\u001b[39maround(val_macro_f1, \u001b[38;5;241m4\u001b[39m)))\n",
      "\u001b[0;31mNameError\u001b[0m: name 'sklearn' is not defined"
     ]
    }
   ],
   "source": [
    "for k in np.arange(NREPS):\n",
    "    for frac in SAMPLING_FRACS:\n",
    "        print(\"frac {}, rep {}\".format(frac, k))\n",
    "        if k==0 and frac==1.0:\n",
    "            continue\n",
    "        #downsample training set\n",
    "        sss = StratifiedShuffleSplit(n_splits=1, test_size=0.2, random_state=k) #same train/val set split for each frac in k\n",
    "        for index_train, index_val in sss.split(data, label):\n",
    "            np.random.seed(k)\n",
    "            index_train_small = np.random.choice(index_train, round(index_train.shape[0]*frac), replace=False)\n",
    "            X_train, y_train = data[index_train_small], label[index_train_small]\n",
    "            X_test, y_test = data[index_val], label[index_val]\n",
    "\n",
    "        print(\"Loaded data...\")\n",
    "\n",
    "        #train on train_dataset\n",
    "        \n",
    "        #hyperparameter tune using k-fold val on training data\n",
    "        cv_results = {}\n",
    "        for c in [1e-3, 1e-2, 1e-1, 1]:\n",
    "            print(\"c={}\".format(c))\n",
    "            lr = LogisticRegression(random_state=0, penalty=\"l1\", C=c, solver=\"liblinear\")\n",
    "            res = cross_validate(lr, X_train, y_train, scoring=['accuracy'])\n",
    "            cv_results[c] = np.mean(res['test_accuracy'])\n",
    "        print(cv_results)\n",
    "\n",
    "        #choose best c and calc performance on val_dataset\n",
    "        best_ind = np.argmax(list(cv_results.values()))\n",
    "        c = list(cv_results.keys())[best_ind]\n",
    "        print(\"best c={}\".format(c))\n",
    "        lr = LogisticRegression(random_state=0, penalty=\"l1\", C=c, solver=\"liblinear\")\n",
    "        lr.fit(X_train, y_train)\n",
    "        print(\"train set accuracy: \" + str(np.around(lr.score(X_train, y_train), 4)))\n",
    "        print(\"test set accuracy: \" + str(np.around(lr.score(X_test, y_test), 4)))\n",
    "        val_macro_f1 = f1_score(y_test, lr.predict(X_test), average=\"macro\")\n",
    "        print(\"test set macro F1: \" + str(np.around(val_macro_f1, 4)))\n",
    "        \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "frac 1.0, rep 0\n",
      "frac 0.75, rep 0\n",
      "frac 0.5, rep 0\n",
      "Loaded data...\n",
      "c=0.001\n",
      "c=0.01\n",
      "c=0.1\n",
      "c=1\n"
     ]
    }
   ],
   "source": [
    "nreps_record = []\n",
    "fracs_record = []\n",
    "c_record = []\n",
    "train_acc_record = []\n",
    "test_acc_record = []\n",
    "test_macrof1_record = []\n",
    "for k in np.arange(NREPS):\n",
    "    for frac in SAMPLING_FRACS:\n",
    "        nreps_record.append(k)\n",
    "        print(\"frac {}, rep {}\".format(frac, k))\n",
    "        fracs_record.append(frac)\n",
    "        if k==0 and frac==1.0:\n",
    "            continue\n",
    "        if k==0 and frac==0.75:\n",
    "            continue\n",
    "        #downsample training set\n",
    "        sss = StratifiedShuffleSplit(n_splits=1, test_size=0.2, random_state=k) #same train/val set split for each frac in k\n",
    "        for index_train, index_val in sss.split(data, label):\n",
    "            np.random.seed(k)\n",
    "            index_train_small = np.random.choice(index_train, round(index_train.shape[0]*frac), replace=False)\n",
    "            X_train, y_train = data[index_train_small], label[index_train_small]\n",
    "            X_test, y_test = data[index_val], label[index_val]\n",
    "\n",
    "        print(\"Loaded data...\")\n",
    "\n",
    "        #train on train_dataset\n",
    "        \n",
    "        #hyperparameter tune using k-fold val on training data\n",
    "        cv_results = {}\n",
    "        for c in [1e-3, 1e-2, 1e-1, 1]:\n",
    "            print(\"c={}\".format(c))\n",
    "            lr = LogisticRegression(random_state=0, penalty=\"l1\", C=c, solver=\"liblinear\")\n",
    "            res = cross_validate(lr, X_train, y_train, scoring=['accuracy'])\n",
    "            cv_results[c] = np.mean(res['test_accuracy'])\n",
    "        print(cv_results)\n",
    "\n",
    "        #choose best c and calc performance on val_dataset\n",
    "        best_ind = np.argmax(list(cv_results.values()))\n",
    "        c = list(cv_results.keys())[best_ind]\n",
    "        c_record.append(c)\n",
    "        print(\"best c={}\".format(c))\n",
    "        lr = LogisticRegression(random_state=0, penalty=\"l1\", C=c, solver=\"liblinear\")\n",
    "        lr.fit(X_train, y_train)\n",
    "        print(\"train set accuracy: \" + str(np.around(lr.score(X_train, y_train), 4)))\n",
    "        print(\"test set accuracy: \" + str(np.around(lr.score(X_test, y_test), 4)))\n",
    "        val_macro_f1 = f1_score(y_test, lr.predict(X_test), average=\"macro\")\n",
    "        print(\"test set macro F1: \" + str(np.around(val_macro_f1, 4)))\n",
    "        train_acc_record.append(np.around(lr.score(X_train, y_train), 4))\n",
    "        test_acc_record.append(np.around(lr.score(X_test, y_test), 4))\n",
    "        test_macrof1_record.append(np.around(val_macro_f1, 4))\n",
    "        \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 63,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2\n",
      "0.5\n",
      "0.1\n",
      "0.7831263696128561\n"
     ]
    }
   ],
   "source": [
    "print(k)\n",
    "print(frac)\n",
    "print(c)\n",
    "print(lr.score(X_test, y_test))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 65,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "frac 0.25, rep 0\n",
      "Loaded data...\n",
      "c=0.001\n",
      "c=0.01\n",
      "c=0.1\n",
      "c=1\n",
      "{0.001: 0.6259313367421475, 0.01: 0.7357195032870709, 0.1: 0.7432432432432432, 1: 0.7035792549306062}\n",
      "best c=0.1\n",
      "train set accuracy: 0.9736\n",
      "test set accuracy: 0.7528\n",
      "test set macro F1: 0.6355\n",
      "frac 0.1, rep 0\n",
      "Loaded data...\n",
      "c=0.001\n",
      "c=0.01\n",
      "c=0.1\n",
      "c=1\n",
      "{0.001: 0.5135149818351499, 0.01: 0.6961260540612605, 0.1: 0.7094595540445955, 1: 0.6765881745158817}\n",
      "best c=0.1\n",
      "train set accuracy: 0.9792\n",
      "test set accuracy: 0.7192\n",
      "test set macro F1: 0.5992\n"
     ]
    }
   ],
   "source": [
    "nreps_record = []\n",
    "fracs_record = []\n",
    "c_record = []\n",
    "train_acc_record = []\n",
    "test_acc_record = []\n",
    "test_macrof1_record = []\n",
    "\n",
    "NREPS=1\n",
    "SAMPLING_FRACS=[0.25,0.1]\n",
    "\n",
    "for k in np.arange(NREPS):\n",
    "    for frac in SAMPLING_FRACS:\n",
    "        nreps_record.append(k)\n",
    "        print(\"frac {}, rep {}\".format(frac, k))\n",
    "        fracs_record.append(frac)\n",
    "        if k==0 and frac==1.0:\n",
    "            continue\n",
    "        if k==0 and frac==0.75:\n",
    "            continue\n",
    "        #downsample training set\n",
    "        sss = StratifiedShuffleSplit(n_splits=1, test_size=0.2, random_state=k) #same train/val set split for each frac in k\n",
    "        for index_train, index_val in sss.split(data, label):\n",
    "            np.random.seed(k)\n",
    "            index_train_small = np.random.choice(index_train, round(index_train.shape[0]*frac), replace=False)\n",
    "            X_train, y_train = data[index_train_small], label[index_train_small]\n",
    "            X_test, y_test = data[index_val], label[index_val]\n",
    "\n",
    "        print(\"Loaded data...\")\n",
    "\n",
    "        #train on train_dataset\n",
    "        \n",
    "        #hyperparameter tune using k-fold val on training data\n",
    "        cv_results = {}\n",
    "        for c in [1e-3, 1e-2, 1e-1, 1]:\n",
    "            print(\"c={}\".format(c))\n",
    "            lr = LogisticRegression(random_state=0, penalty=\"l1\", C=c, solver=\"liblinear\")\n",
    "            res = cross_validate(lr, X_train, y_train, scoring=['accuracy'])\n",
    "            cv_results[c] = np.mean(res['test_accuracy'])\n",
    "        print(cv_results)\n",
    "\n",
    "        #choose best c and calc performance on val_dataset\n",
    "        best_ind = np.argmax(list(cv_results.values()))\n",
    "        c = list(cv_results.keys())[best_ind]\n",
    "        c_record.append(c)\n",
    "        print(\"best c={}\".format(c))\n",
    "        lr = LogisticRegression(random_state=0, penalty=\"l1\", C=c, solver=\"liblinear\")\n",
    "        lr.fit(X_train, y_train)\n",
    "        print(\"train set accuracy: \" + str(np.around(lr.score(X_train, y_train), 4)))\n",
    "        print(\"test set accuracy: \" + str(np.around(lr.score(X_test, y_test), 4)))\n",
    "        val_macro_f1 = f1_score(y_test, lr.predict(X_test), average=\"macro\")\n",
    "        print(\"test set macro F1: \" + str(np.around(val_macro_f1, 4)))\n",
    "        train_acc_record.append(np.around(lr.score(X_train, y_train), 4))\n",
    "        test_acc_record.append(np.around(lr.score(X_test, y_test), 4))\n",
    "        test_macrof1_record.append(np.around(val_macro_f1, 4))\n",
    "        \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# logistic regression model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "cv_results = {}\n",
    "for c in [1e-6, 1e-4, 1e-2, 1]:\n",
    "    cv_results[str(c)] = {}\n",
    "    lr = LogisticRegression(random_state=0, penalty=\"l1\", C=c, solver=\"liblinear\")\n",
    "    res = cross_validate(lr, X_train, y_train, scoring=['accuracy','f1_macro'])\n",
    "    cv_results[c]['test_accuracy'] = np.mean(res['test_accuracy'])\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'1e-06': {'test_accuracy': 0.04158144631117604,\n",
       "  'test_f1_macro': 0.007258445949778784},\n",
       " '0.0001': {'test_accuracy': 0.5174762600438276,\n",
       "  'test_f1_macro': 0.20680632503188856},\n",
       " '0.01': {'test_accuracy': 0.7841672753834917,\n",
       "  'test_f1_macro': 0.6449686172460548},\n",
       " '1': {'test_accuracy': 0.7628195763330898,\n",
       "  'test_f1_macro': 0.6593903887882162}}"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "cv_results"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'0.001': {'test_accuracy': 0.6888787436084733,\n",
       "  'test_f1_macro': 0.48298639855688597},\n",
       " '0.1': {'test_accuracy': 0.800018261504748,\n",
       "  'test_f1_macro': 0.6919380706270787}}"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "cv_results = {}\n",
    "for c in [1e-3, 1e-1]:\n",
    "    cv_results[str(c)] = {}\n",
    "    lr = LogisticRegression(random_state=0, penalty=\"l1\", C=c, solver=\"liblinear\")\n",
    "    res = cross_validate(lr, X_train, y_train, scoring=['accuracy','f1_macro'])\n",
    "    cv_results[str(c)]['test_accuracy'] = np.mean(res['test_accuracy'])\n",
    "    cv_results[str(c)]['test_f1_macro'] = np.mean(res['test_f1_macro'])\n",
    "cv_results"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.8116143170197224"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#based on cv results\n",
    "c = 0.1\n",
    "lr = LogisticRegression(random_state=0, penalty=\"l1\", C=c, solver=\"liblinear\")\n",
    "lr.fit(X_train, y_train)\n",
    "lr.score(X_test, y_test)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.7074023856802742"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sklearn.metrics.f1_score(y_test, lr.predict(X_test), average=\"macro\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "CD14+ Monocyte: 0.8982456140350877\n",
      "CD19+ B: 0.8382978723404255\n",
      "CD34+: 0.8125\n",
      "CD4+ T Helper2: 0.05\n",
      "CD4+/CD25 T Reg: 0.6814874696847211\n",
      "CD4+/CD45RA+/CD25- Naive T: 0.4732620320855615\n",
      "CD4+/CD45RO+ Memory: 0.4297385620915033\n",
      "CD56+ NK: 0.9264957264957265\n",
      "CD8+ Cytotoxic T: 0.8075162611418936\n",
      "CD8+/CD45RA+ Naive Cytotoxic: 0.8993691799339141\n",
      "Dendritic: 0.7613365155131265\n"
     ]
    }
   ],
   "source": [
    "#check accuracy per tissue class\n",
    "for ct in np.unique(zheng_data.obs.celltype):\n",
    "    print(ct+\": {}\".format(lr.score(X_test[y_test==ct], y_test[y_test==ct])))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "rna-pretrain",
   "language": "python",
   "name": "rna-pretrain"
  },
  "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.9.12"
  },
  "vscode": {
   "interpreter": {
    "hash": "ea531a0049080d5658e3e61cf7db54ce403b028279b105adc411fb7eae04af57"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}