{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# What is about ?\n\nHere is a custom class cross-validation class which supports proposed validation schemes AmbrosM, MT, etc  as well as standard sklearn random \"Kfolds\". It can be used in the same interface as sklearn \"Kfolds\":\n\n    for i, (train_index, test_index) in enumerate(kf.split(X)):\n        X_train,X_test = X_full[IX_train], X_full[IX_test]\n        \nIt can also support returning 3 sub-folds: train, validation, test . Where validation is used for a kind of \"early stopping\" for neural networks and boostings. \n\nCV-schemes by AmbrosM and MT are built-in:\n\n    kf = KFold_custom(CV_scheme = 'MT')\n    kf = KFold_custom(CV_scheme = 'AmbrosM') \n    kf = KFold_custom(CV_scheme = 'Random') # wrapper for sklearn Kfolds with shuffle = True\n    kf = KFold_custom(CV_scheme = CV_scheme, random_state = 42 ) # fixing random_state ensures reprodicibility of folds. Fixed by 42 by default\n\nTo work with triples: \n\n    kf = KFold_custom(CV_scheme = CV_scheme, valid_size = 0.2 )\n    for i_fold,(IX_train ,IX_valid,IX_test)   in enumerate( kf.split()) :\n        X_train,X_valid,X_test = X_full[IX_train], X_full[IX_valid], X_full[IX_test]\n        \nHere the train part of underlying CV-scheme will be split into two subparts - real-train and validation, where validation can be used for e.g. early stopping or choosing epochs number. While \"test\" part will be the same as in the underlying scheme ! The size of validation part is controlled by valid_size : it takes \"pre-train\" and cuts from it \"valid_size\" part to create validation. As usually valid_size is from 0 to 1. Zero (default) implies absense of validation subfold and train and test will be returned by kf.split as usual sklearn does. \n\n\nLater reconsider the simple Ridge model considered by Ambros and MT. And detailze the comparaison of the validation schemes.\nThe optimial value is around alpha = 1 for both schemes ( 2 - for AmbrosM, 0.5 for MT). The model gives LB0.613, and blended with \"priors\" LB 0.605. \n\nMore analysis of the CV schemes can be found in slides: https://docs.google.com/presentation/d/1wiz0Wmt4D54pqMMsIOyJHuQYMZ3hTBZQQnjbLzwoGYY/edit?usp=sharing and sheet: https://docs.google.com/spreadsheets/d/1APN63PMaWZygVjYimK9Ivt0RvifdAU5JRYkxiDn4szw/edit?usp=sharing\n\n\nPS \n\nSee related notebook with NN analysis: \nhttps://www.kaggle.com/code/alexandervc/op2-kishan-s-nn-streamlined-and-blended\n\nCV schemes (see discussion https://www.kaggle.com/competitions/open-problems-single-cell-perturbations/discussion/444494 ):\n\nAmbrosM: https://www.kaggle.com/code/ambrosm/scp-quickstart\n\nMT: https://www.kaggle.com/code/masato114/scp-quickstart-another-cv-strategy\n\nKishan: https://www.kaggle.com/code/liudacheldieva/neural-network-regression/notebook#Predicting-on-test-data \n","metadata":{}},{"cell_type":"code","source":"# This Python 3 environment comes with many helpful analytics libraries installed\n# It is defined by the kaggle/python Docker image: https://github.com/kaggle/docker-python\n# For example, here's several helpful packages to load\n\nimport time\nt0start = time.time() \n\nimport numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\n\nimport pandas as pd\nimport tensorflow as tf\nimport numpy as np\nimport matplotlib.pyplot as plt\nimport seaborn as sns\n\n# Input data files are available in the read-only \"../input/\" directory\n# For example, running this (by clicking run or pressing Shift+Enter) will list all files under the input directory\n\nimport os\nfor dirname, _, filenames in os.walk('/kaggle/input'):\n    for filename in filenames:\n        print(os.path.join(dirname, filename))\n\n# You can write up to 20GB to the current directory (/kaggle/working/) that gets preserved as output when you create a version using \"Save & Run All\" \n# You can also write temporary files to /kaggle/temp/, but they won't be saved outside of the current session","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2023-11-11T21:33:00.520232Z","iopub.execute_input":"2023-11-11T21:33:00.520704Z","iopub.status.idle":"2023-11-11T21:33:09.704482Z","shell.execute_reply.started":"2023-11-11T21:33:00.520668Z","shell.execute_reply":"2023-11-11T21:33:09.703143Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Load data","metadata":{}},{"cell_type":"code","source":"%%time\nfn = '/kaggle/input/open-problems-single-cell-perturbations/de_train.parquet'\ndf_de_train = pd.read_parquet(fn)# , index_col = 0)\nprint(df_de_train.shape)\ndisplay(df_de_train )\n\nplt.figure(figsize = (20,4) )\nv = df_de_train.iloc[:,5:].max(axis = 0 ).sort_values(ascending = False, key = abs )\nplt.plot(v.values,'*-')\nplt.title('Max-abs DE for genes',fontsize = 20 )\nplt.grid()\nplt.show()\ndisplay(v.head(15))\n\n# %%time\nfn = '/kaggle/input/open-problems-single-cell-perturbations/id_map.csv'\ndf_id_map = pd.read_csv(fn,index_col = 0)\nprint(df_id_map.shape)\ndisplay(df_id_map)\nfn = '/kaggle/input/open-problems-single-cell-perturbations/sample_submission.csv'\ndf_sample_submit = pd.read_csv(fn, index_col = 0)\nprint(df_sample_submit.shape)\ndisplay( df_sample_submit )","metadata":{"execution":{"iopub.status.busy":"2023-11-11T21:33:09.707172Z","iopub.execute_input":"2023-11-11T21:33:09.708091Z","iopub.status.idle":"2023-11-11T21:33:18.254409Z","shell.execute_reply.started":"2023-11-11T21:33:09.708040Z","shell.execute_reply":"2023-11-11T21:33:18.252929Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Custom CV (AmbrosM, MT, etc )","metadata":{}},{"cell_type":"markdown","source":"## Preprations for CV class\n\nPrepare indices for AmbrosM and MT scheme\n\nCV schemes (see discussion https://www.kaggle.com/competitions/open-problems-single-cell-perturbations/discussion/444494 ):\n\nAmbrosM: https://www.kaggle.com/code/ambrosm/scp-quickstart\n\nMT: https://www.kaggle.com/code/masato114/scp-quickstart-another-cv-strategy\n\nKishan's notebook: https://www.kaggle.com/code/liudacheldieva/neural-network-regression/notebook#Predicting-on-test-data ","metadata":{}},{"cell_type":"code","source":"# Prepare indexes for AmbrosM scheme:\n# For AmbrosM scheme: \nfolds_index_data_AmbrosM = [ ]\ntrain_sm_names = ['Idelalisib', 'Crizotinib', 'Linagliptin', 'Palbociclib', 'Dabrafenib', 'Alvocidib', 'LDN 193189', 'R428', 'Porcn Inhibitor III', \n  'Belinostat', 'Foretinib', 'MLN 2238', 'Penfluridol', 'Dactolisib', 'O-Demethylated Adapalene', 'Oprozomib (ONX 0912)', 'CHIR-99021']\nlist_fold_ids =  ['NK cells', 'T cells CD4+', 'T cells CD8+', 'T regulatory cells'] \nfor fold_id in list_fold_ids:\n    mask_va = (df_de_train.cell_type == fold_id) & ~df_de_train.sm_name.isin(train_sm_names)\n    mask_tr = ~mask_va # 485 or 487 training rows\n    IX_train = np.where( mask_tr > 0 )[0]\n    IX_test = np.where( mask_va > 0 )[0]\n    #print(fold_id,  len(IX_test), type(IX_test), IX_test[:3], len(IX_train), type(IX_train), IX_train[:3]  )\n    folds_index_data_AmbrosM.append( [IX_train, IX_test ])\n\n# For MT CV scheme:\nfolds_index_data_MT = []\nfold_to_compounds = {0: ['Alvocidib', 'Belinostat', 'Foretinib', 'LDN 193189',  'Linagliptin', 'O-Demethylated Adapalene'],\n 1: ['Dabrafenib', 'Dactolisib', 'Idelalisib', 'MLN 2238', 'Palbociclib', 'Porcn Inhibitor III'],\n 2: ['CHIR-99021', 'Crizotinib', 'Oprozomib (ONX 0912)', 'Penfluridol',  'R428']}\nfor fold_id in [0,1,2]:\n    mask_va = df_de_train['cell_type'].isin(['Myeloid cells', 'B cells']) & df_de_train['sm_name'].isin(fold_to_compounds[fold_id])\n    mask_tr = ~mask_va    \n    IX_train = np.where( mask_tr > 0 )[0]\n    IX_test = np.where( mask_va > 0 )[0]\n    #print(fold_id,  len(IX_test), type(IX_test), IX_test[:3], len(IX_train), type(IX_train), IX_train[:3]  )\n    folds_index_data_MT.append( [IX_train, IX_test ])    \n    \n    \ndict_folds_Antonina = {0: [4, 15, 18, 41, 47, 53, 56, 57, 59, 61, 62, 66, 78, 83, 86, 93, 100, 110, 115, 116, 117, 120, 123, 131, 148, 152, 161, 185, 189, 193, 197, 198, 205, 207, 209, 210, 213, 218, 219, 225, 228, 230, 250, 252, 257, 263, 286, 292, 294, 304, 306, 310, 325, 327, 336, 338, 339, 350, 351, 356, 357, 360, 362, 369, 379, 385, 395, 419, 424, 430, 432, 434, 438, 441, 443, 445, 448, 452, 456, 465, 467, 471, 472, 474, 498, 508, 523, 534, 542, 544, 574, 581, 584, 589, 601, 607], 1: [1, 3, 5, 14, 24, 38, 40, 49, 50, 54, 55, 60, 79, 82, 87, 102, 112, 114, 118, 119, 127, 132, 146, 149, 166, 179, 183, 184, 186, 190, 192, 201, 203, 204, 216, 227, 242, 253, 261, 264, 271, 289, 295, 296, 299, 300, 309, 323, 324, 329, 333, 337, 354, 387, 390, 393, 406, 408, 415, 416, 420, 421, 422, 436, 442, 449, 453, 458, 462, 463, 469, 484, 487, 489, 490, 492, 497, 501, 505, 506, 509, 513, 515, 516, 546, 548, 566, 569, 585, 587, 588, 592, 593, 595, 600, 605], 2: [0, 16, 19, 44, 48, 51, 52, 58, 63, 65, 70, 81, 85, 88, 89, 92, 111, 113, 122, 125, 135, 137, 139, 150, 164, 175, 191, 202, 206, 211, 217, 220, 221, 229, 239, 243, 247, 248, 251, 258, 259, 260, 267, 270, 274, 281, 288, 290, 291, 303, 322, 326, 328, 332, 347, 359, 364, 370, 383, 384, 389, 392, 401, 423, 428, 454, 455, 461, 464, 494, 496, 504, 510, 511, 524, 525, 526, 530, 540, 541, 551, 565, 568, 570, 572, 573, 575, 577, 579, 580, 597, 598, 603, 606, 610, 611, 613], 3: [21, 22, 26, 36, 42, 45, 46, 64, 67, 69, 71, 80, 84, 90, 101, 103, 134, 136, 140, 147, 151, 160, 163, 165, 167, 174, 176, 177, 180, 181, 182, 187, 188, 195, 196, 199, 200, 212, 215, 226, 231, 240, 241, 249, 255, 256, 262, 265, 269, 285, 287, 297, 298, 307, 312, 319, 321, 330, 340, 342, 344, 345, 348, 358, 366, 368, 382, 391, 400, 403, 404, 405, 425, 427, 431, 433, 435, 444, 450, 451, 457, 459, 460, 473, 485, 499, 500, 502, 503, 507, 512, 514, 527, 528, 529, 531, 543, 547, 553, 554, 564, 571, 576, 583, 586, 591, 596, 602, 604, 608, 609], 4: [2, 6, 7, 17, 20, 23, 25, 27, 28, 29, 37, 39, 43, 68, 91, 121, 124, 126, 128, 129, 130, 133, 138, 153, 162, 178, 194, 208, 214, 222, 223, 224, 232, 244, 245, 246, 254, 266, 268, 272, 273, 282, 283, 284, 293, 301, 302, 305, 308, 311, 320, 331, 334, 335, 341, 343, 346, 349, 352, 353, 355, 361, 363, 365, 367, 377, 378, 380, 381, 386, 388, 394, 396, 397, 398, 399, 402, 407, 417, 418, 426, 429, 437, 439, 440, 446, 447, 466, 468, 470, 481, 482, 483, 486, 488, 491, 493, 495, 532, 533, 545, 549, 550, 552, 555, 562, 563, 567, 578, 582, 590, 594, 599, 612], 5: [8, 9, 10, 11, 12, 13, 30, 31, 32, 33, 34, 35, 72, 73, 74, 75, 76, 77, 94, 95, 96, 97, 98, 99, 104, 105, 106, 107, 108, 109, 141, 142, 143, 144, 145, 154, 155, 156, 157, 158, 159, 168, 169, 170, 171, 172, 173, 233, 234, 235, 236, 237, 238, 275, 276, 277, 278, 279, 280, 313, 314, 315, 316, 317, 318, 371, 372, 373, 374, 375, 376, 409, 410, 411, 412, 413, 414, 475, 476, 477, 478, 479, 480, 517, 518, 519, 520, 521, 522, 535, 536, 537, 538, 539, 556, 557, 558, 559, 560, 561]}\n\nfolds_index_data_Antonina = []    \nfor i in range(5):  \n    l_valid = np.array( dict_folds_Antonina[i] )\n    l_train = np.array( [k for k in range(614) if k not in l_valid ] )\n    folds_index_data_Antonina.append( [l_train, l_valid ])    \nprint( len(folds_index_data_Antonina ) )\n    ","metadata":{"execution":{"iopub.status.busy":"2023-11-11T21:33:18.256433Z","iopub.execute_input":"2023-11-11T21:33:18.256787Z","iopub.status.idle":"2023-11-11T21:33:18.324325Z","shell.execute_reply.started":"2023-11-11T21:33:18.256756Z","shell.execute_reply":"2023-11-11T21:33:18.322904Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## CV-custom class","metadata":{}},{"cell_type":"code","source":"import random \nfrom sklearn.model_selection import KFold\n\nclass KFold_custom:\n    '''\n    Class with similar to sklearn \"KFold\" class, interface.\n    Support of CV schemes proposed by AmbrosM, MT, Kishan and standard sklearn KFold\n    Supports options to return IX_train,IX_valid,IX_test - triplet - if valid_size > 0\n    Examples:\n    kf = KFold_custom('AmbrosM'):     \n    for i_fold,(IX_train,IX_test)   in enumerate( kf.split()) :\n        pass\n    kf = KFold_custom('AmbrosM', valid_size = 0.1 )     \n    for i_fold,(IX_train ,IX_valid,IX_test)   in enumerate( kf.split()) :\n        pass\n        \n    kf = KFold_custom(CV_scheme = 'Random') # wrapper for sklearn Kfolds with shuffle = True\n    kf = KFold_custom(CV_scheme = CV_scheme, random_state = 42 ) # fixing random_state ensures reprodicibility of folds \n    '''\n    def __init__(self, CV_scheme,  valid_size = 0, n_splits=5, random_state = 42, verbose = 0 ):\n        self.CV_scheme = CV_scheme\n        self.CV_scheme_inf =  CV_scheme\n        self.valid_size = np.clip(valid_size,0,1)\n        self.random_state = random_state\n\n        self.folds_index_data = [ [np.arange(0), np.arange(0), np.arange(0)] ] # Example data\n        self.n_splits = 1 # Example data\n        if CV_scheme == 'Kishan1':\n            # One fold scheme which contains train, valid, test parts \n            IX_train_Kishan1 = np.arange(429)\n            IX_val_Kishan1 = np.arange(429,521)\n            IX_test_Kishan1 = np.arange(521,614)\n            self.n_splits = 1\n            self.folds_index_data =  [   [IX_train_Kishan1, IX_val_Kishan1, IX_test_Kishan1]  ]\n        elif CV_scheme == 'Tonya':\n            self.n_splits = 5\n            self.folds_index_data = folds_index_data_Antonina\n        elif CV_scheme == 'AmbrosM':\n            self.n_splits = 4\n            self.folds_index_data = folds_index_data_AmbrosM\n        elif CV_scheme == 'MT':\n            self.n_splits = 3\n            self.folds_index_data = folds_index_data_MT\n        elif 'Random'.lower() in CV_scheme.lower():\n            self.n_splits = n_splits\n            self.CV_scheme_inf = 'Random_'+str(n_splits) +'_'+str(random_state )\n            kf = KFold(n_splits=n_splits, random_state = random_state, shuffle=True )#,\n            self.folds_index_data = list( kf.split( np.arange(614) ) )\n        elif 'Full'.lower() in CV_scheme.lower():\n            # Just return the full set as both train and test - it useful to re-train the model on the entire data\n            IX_train_full = np.arange(614)\n            IX_test_full = np.arange(614)\n            self.n_splits = 1\n            self.folds_index_data = [   [IX_train_full, IX_test_full]  ]\n        else:\n            s = 'Uncrecognized CV_scheme ' + str(CV_scheme) \n            raise ValueError(s)\n\n        self.verbose = verbose \n        if verbose >= 10:\n            print(self.CV_scheme, self.n_splits, len(self.folds_index_data) )\n            #print(self.folds_index_data)\n            \n    def split(self, X=None):\n        '''\n        X - NOT used, just for compatibility with sklearn \n        '''\n        for item in self.folds_index_data:\n            if self.CV_scheme in ['Kishan1']:\n                if self.valid_size > 0:\n                    yield item[0],item[1],item[2]\n                else :\n                    yield np.array( list(set(item[0])|set(item[1]))  ) , item[2] # Return train, test only \n            elif self.CV_scheme in ['Tonya', 'AmbrosM','MT', 'Random', 'Full']:\n                if self.valid_size == 0:\n                    yield item[0],item[1]\n                else:\n                    # Split \"full-train\"->( real-train, valid ) \n                    index_for_real_train_part = int(  (1-self.valid_size) * len(item[0]) )\n                    if self.random_state is None:\n                        IX_train =  item[0][:index_for_real_train_part]\n                        IX_valid =  item[0][index_for_real_train_part:]\n                    elif self.random_state == -1:\n                        p = np.random.permutation(len(item[0])) \n                        IX_train =  item[0][p][:index_for_real_train_part]\n                        IX_valid =  item[0][p][index_for_real_train_part:]\n                    else:\n                        np.random.seed(self.random_state) # Set temporary random seed\n                        p = np.random.permutation(len(item[0])) \n                        np.random.seed(random.randint(0,30000)) # Randomize seed again - use Python random, not numpy       \n                        IX_train =  item[0][p][:index_for_real_train_part]\n                        IX_valid =  item[0][p][index_for_real_train_part:]\n                    yield IX_train, IX_valid, item[1]\n                    \n        \n    def get_list_main_CV_schemes(self):\n        return ['AmbrosM','MT','Kishan1','Random']\n    def get_n_splits(self, X=np.array([])):\n        return self.n_splits\n        \n    \nprint();print();    \nprint('--------------------------------Examples-----------------------------------------------')\nprint();print();    \n        \nCV_scheme = 'Kishan1'\nprint(CV_scheme,  )\nkf = KFold_custom('Kishan1', valid_size = 0.18)     \n#list( kf.split() )\nprint(\" i_fold, len(IX_train),len(IX_valid),len(IX_test), IX_train[:3],IX_valid[:3], IX_test[:3],  type(IX_train),type(IX_valid),type(IX_test)  \" )\nfor i_fold,(IX_train ,IX_valid,IX_test)   in enumerate( kf.split(None)) :\n    print(i_fold, len(IX_train),len(IX_valid),len(IX_test), IX_train[:3],IX_valid[:3], IX_test[:3],  type(IX_train),type(IX_valid),type(IX_test) )\nprint()\n\nfor CV_scheme in ['Tonya','AmbrosM','MT','Random']:\n    valid_size = 0\n    print(CV_scheme, 'valid_size', valid_size )\n    kf = KFold_custom(CV_scheme, valid_size = valid_size)     \n    #list( kf.split() )\n    print(\" i_fold, len(IX_train),len(IX_test), IX_train[:3], IX_test[:3],  type(IX_train),type(IX_test)  \" )\n    for i_fold,(IX_train,IX_test)   in enumerate( kf.split(None)) :\n        print(i_fold, len(IX_train), len(IX_test), IX_train[:3], IX_test[:3],  type(IX_train), type(IX_test) )\n    \n    print()\n    valid_size = 0.2 \n    print(CV_scheme, 'valid_size', valid_size )\n    kf = KFold_custom(CV_scheme, valid_size = valid_size )     \n    #list( kf.split() )\n    print(\" i_fold, len(IX_train),len(IX_valid),len(IX_test), IX_train[:3],IX_valid[:3], IX_test[:3],  type(IX_train),type(IX_valid),type(IX_test)  \" )\n    for i_fold,(IX_train ,IX_valid,IX_test)   in enumerate( kf.split(None)) :\n        print(i_fold, len(IX_train),len(IX_valid),len(IX_test), IX_train[:3],IX_valid[:3], IX_test[:3],  type(IX_train),type(IX_valid),type(IX_test) )\n    print()","metadata":{"execution":{"iopub.status.busy":"2023-11-11T21:33:18.327953Z","iopub.execute_input":"2023-11-11T21:33:18.328807Z","iopub.status.idle":"2023-11-11T21:33:18.565983Z","shell.execute_reply.started":"2023-11-11T21:33:18.328770Z","shell.execute_reply":"2023-11-11T21:33:18.564633Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Encoding ","metadata":{}},{"cell_type":"code","source":"%%time\n\nfrom sklearn.preprocessing import OneHotEncoder\n\nfeatures_columns = ['cell_type','sm_name']\nselected_features_columns = [ 'sm_name']\nX_submit_categorical = df_id_map[selected_features_columns] # pd.DataFrame(df_id_map, columns= features_columns )\nX_full_categorical = df_de_train[ selected_features_columns ]\n\n\n# Create an instance of the encoder\nencoder = OneHotEncoder()\n\n\nX_full = encoder.fit_transform(X_full_categorical).toarray()\nX_submit = encoder.transform(X_submit_categorical).toarray()\n\nY_full = df_de_train.iloc[:,5:].values\n\nprint('X_full.shape, X_submit.shape ,  Y_full.shape', X_full.shape, X_submit.shape ,  Y_full.shape)\n","metadata":{"execution":{"iopub.status.busy":"2023-11-11T21:33:18.567424Z","iopub.execute_input":"2023-11-11T21:33:18.567757Z","iopub.status.idle":"2023-11-11T21:33:18.626760Z","shell.execute_reply.started":"2023-11-11T21:33:18.567728Z","shell.execute_reply":"2023-11-11T21:33:18.625431Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Ridge CV scoring. Simple example ","metadata":{}},{"cell_type":"markdown","source":"## MT CV-scheme","metadata":{}},{"cell_type":"code","source":"%%time \nfrom sklearn.linear_model import Ridge\nfrom sklearn.decomposition import TruncatedSVD\nfrom sklearn.metrics import r2_score\n\nverbose = 1\n\nn_components = 30\nreducer = TruncatedSVD(n_components=n_components, n_iter=7, random_state=42)\n\nalpha = 1\nmodel = Ridge(alpha=alpha)\n\nCV_scheme = 'MT'\nprint('CV_scheme:', CV_scheme)\nvalid_size = 0 \nkf = KFold_custom(CV_scheme, valid_size = valid_size)   \n\nlist_mrrmse = [];list_r2 = []\nfor i_fold,(IX_train,IX_test)   in enumerate( kf.split()) :\n    X_train,X_test = X_full[IX_train], X_full[IX_test]\n    \n    Y_train_red = reducer.fit_transform(Y_full[IX_train] ) \n    if verbose >= 10:\n        print('fold:', i_fold, 'X_train.shape, X_test.shape, Y_train_red.shape', X_train.shape, X_test.shape, Y_train_red.shape)\n    \n    model.fit(X_train, Y_train_red)\n    Y_test_pred_red = model.predict( X_test )\n    Y_test_pred = reducer.inverse_transform( Y_test_pred_red )\n    \n    Y_test = Y_full[IX_test]\n    mrrmse = np.sqrt(np.square(Y_test - Y_test_pred).mean(axis=1)).mean();  r2 = r2_score( Y_test , Y_test_pred ) \n    print('fold:', i_fold, 'mrrmse:', np.round(mrrmse,4), 'r2:', np.round(r2,4),  )\n    list_mrrmse.append(mrrmse); list_r2.append(r2)\nprint()\nprint('Average mrrmse:', np.round(np.mean(list_mrrmse ),4 ))\nprint('Average r2:', np.round(np.mean(list_r2 ),4 ))","metadata":{"execution":{"iopub.status.busy":"2023-11-11T21:33:18.627919Z","iopub.execute_input":"2023-11-11T21:33:18.628261Z","iopub.status.idle":"2023-11-11T21:33:24.369898Z","shell.execute_reply.started":"2023-11-11T21:33:18.628234Z","shell.execute_reply":"2023-11-11T21:33:24.368272Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## AmbrosM CV scheme","metadata":{}},{"cell_type":"code","source":"%%time \nfrom sklearn.linear_model import Ridge\nfrom sklearn.decomposition import TruncatedSVD\nfrom sklearn.metrics import r2_score\n\nverbose = 1\n\nn_components = 30\nreducer = TruncatedSVD(n_components=n_components, n_iter=7, random_state=42)\n\nalpha = 1\nmodel = Ridge(alpha=alpha)\n\nCV_scheme = 'AmbrosM'\nprint('CV_scheme:', CV_scheme)\nvalid_size = 0 \nkf = KFold_custom(CV_scheme, valid_size = valid_size)   \n\nlist_mrrmse = [];list_r2 = []\nfor i_fold,(IX_train,IX_test)   in enumerate( kf.split()) :\n    X_train,X_test = X_full[IX_train], X_full[IX_test]\n    \n    Y_train_red = reducer.fit_transform(Y_full[IX_train] ) \n    if verbose >= 10:\n        print('fold:', i_fold, 'X_train.shape, X_test.shape, Y_train_red.shape', X_train.shape, X_test.shape, Y_train_red.shape)\n    \n    model.fit(X_train, Y_train_red)\n    Y_test_pred_red = model.predict( X_test )\n    Y_test_pred = reducer.inverse_transform( Y_test_pred_red )\n    \n    Y_test = Y_full[IX_test]\n    mrrmse = np.sqrt(np.square(Y_test - Y_test_pred).mean(axis=1)).mean();  r2 = r2_score( Y_test , Y_test_pred ) \n    print('fold:', i_fold, 'mrrmse:', np.round(mrrmse,4), 'r2:', np.round(r2,4),  )\n    list_mrrmse.append(mrrmse); list_r2.append(r2)\nprint()\nprint('Average mrrmse:', np.round(np.mean(list_mrrmse ),4 ))\nprint('Average r2:', np.round(np.mean(list_r2 ),4 ))","metadata":{"execution":{"iopub.status.busy":"2023-11-11T21:33:24.378026Z","iopub.execute_input":"2023-11-11T21:33:24.379115Z","iopub.status.idle":"2023-11-11T21:33:31.705555Z","shell.execute_reply.started":"2023-11-11T21:33:24.379052Z","shell.execute_reply":"2023-11-11T21:33:31.704325Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Compare several CV. Computation (7min)","metadata":{}},{"cell_type":"code","source":"%%time \nfrom sklearn.linear_model import Ridge\nfrom sklearn.decomposition import TruncatedSVD\nfrom sklearn.metrics import r2_score\n\nverbose = 1000\n\nn_components = 30\nreducer = TruncatedSVD(n_components=n_components, n_iter=7, random_state=42)\n\nalpha = 1\nmodel = Ridge(alpha=alpha)\n\n# CV_scheme = 'MT'\n# print('CV_scheme:', CV_scheme)\nvalid_size = 0 \ndf_stat = pd.DataFrame();\n\nlist_alpha = [0.0001,0.0002,0.0005,0.001, 0.002, 0.005, 0.01,0.02, 0.05,0.1,0.2,0.5,1,2,5,10,20,50 ,1e2,2e2,5e2] # 1e3,1e4,1e5,1e6]:\n# list_alpha = [0.001,1,2] # 0.0001,0.0002,0.0005,0.001, 0.002, 0.005, 0.01,0.02, 0.05,0.1,0.2,0.5,1,2,5,10,20,50 ,1e2,2e2,5e2] # 1e3,1e4,1e5,1e6]:\nlist_CV_scheme =  ['Tonya', 'AmbrosM','MT','Random']\n\nfor CV_scheme in  list_CV_scheme:\n    print();print('CV_scheme', CV_scheme)\n    list_mrrmse = []\n    IX_stat = 0;\n    for alpha in list_alpha:\n        if verbose > 0:\n            print('alpha:',alpha)\n        model = Ridge(alpha=alpha)\n        list_mrrmse_folds = []; list_r2_folds = []\n        kf = KFold_custom(CV_scheme, valid_size = valid_size)   \n        for i_fold,(IX_train,IX_test)   in enumerate( kf.split()) :\n            X_train,X_test = X_full[IX_train], X_full[IX_test]\n\n            Y_train_red = reducer.fit_transform(Y_full[IX_train] ) \n            if verbose >= 10:\n                print('fold:', i_fold, 'X_train.shape, X_test.shape, Y_train_red.shape', X_train.shape, X_test.shape, Y_train_red.shape)\n\n            model.fit(X_train, Y_train_red)\n            Y_test_pred_red = model.predict( X_test )\n            Y_test_pred = reducer.inverse_transform( Y_test_pred_red )\n\n            Y_test = Y_full[IX_test]\n            mrrmse = np.sqrt(np.square(Y_test - Y_test_pred).mean(axis=1)).mean();  r2 = r2_score( Y_test , Y_test_pred ) \n            print('fold:', i_fold, 'mrrmse:', np.round(mrrmse,4), 'r2:', np.round(r2,4),  )\n            list_mrrmse_folds.append( mrrmse ); list_r2_folds.append(r2)\n\n            df_stat.loc[IX_stat,'alpha' ] = alpha\n            df_stat.loc[IX_stat,'mrrmse fold '+str(i_fold) +  ' ' + CV_scheme] = mrrmse\n            df_stat.loc[IX_stat,'r2 fold '+str(i_fold) + ' ' + CV_scheme] = r2\n        IX_stat += 1\n\n        list_mrrmse.append(np.mean( list_mrrmse_folds ))\n        print('Average:', 'mrrmse:', np.round(np.mean(list_mrrmse_folds),4), 'r2:',np.round(np.mean(list_r2_folds),4)  )\ndf_stat    ","metadata":{"execution":{"iopub.status.busy":"2023-11-11T21:33:31.707417Z","iopub.execute_input":"2023-11-11T21:33:31.707881Z","iopub.status.idle":"2023-11-11T21:43:59.085559Z","shell.execute_reply.started":"2023-11-11T21:33:31.707811Z","shell.execute_reply":"2023-11-11T21:43:59.084646Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Scores correlations for different CV-schemes","metadata":{}},{"cell_type":"code","source":"%%time\nfor str_score in ['mrrmse', 'r2']:\n    df_stat_averages_over_folds = pd.DataFrame()\n    for CV_scheme in  list_CV_scheme:\n        list_selected_cols = [t for t in df_stat.columns if (CV_scheme in t ) and (str_score  in t ) ]\n        #print(list_selected_cols)\n        df_stat_averages_over_folds[CV_scheme] = df_stat[list_selected_cols].mean(axis = 1)\n    print(str_score, 'correlations for different CV-schemes:')\n    display(df_stat_averages_over_folds.corr().round(2))\n","metadata":{"execution":{"iopub.status.busy":"2023-11-11T21:52:11.320143Z","iopub.execute_input":"2023-11-11T21:52:11.320624Z","iopub.status.idle":"2023-11-11T21:52:11.370752Z","shell.execute_reply.started":"2023-11-11T21:52:11.320588Z","shell.execute_reply":"2023-11-11T21:52:11.369603Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_stat.to_csv('df_stat.csv')\ndf_stat_save = df_stat.copy()","metadata":{"execution":{"iopub.status.busy":"2023-11-11T21:43:59.145864Z","iopub.execute_input":"2023-11-11T21:43:59.146339Z","iopub.status.idle":"2023-11-11T21:43:59.160226Z","shell.execute_reply.started":"2023-11-11T21:43:59.146280Z","shell.execute_reply":"2023-11-11T21:43:59.158550Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"pd.set_option('display.max_columns', 500)\npd.set_option('display.max_rows', 500)\n\nprint( df_stat.shape )\ndisplay( df_stat )\n","metadata":{"execution":{"iopub.status.busy":"2023-11-11T21:43:59.161739Z","iopub.execute_input":"2023-11-11T21:43:59.162134Z","iopub.status.idle":"2023-11-11T21:43:59.235136Z","shell.execute_reply.started":"2023-11-11T21:43:59.162100Z","shell.execute_reply":"2023-11-11T21:43:59.233942Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Plots CV compare","metadata":{"execution":{"iopub.status.busy":"2023-10-24T16:50:49.399441Z","iopub.execute_input":"2023-10-24T16:50:49.399881Z","iopub.status.idle":"2023-10-24T16:50:49.405087Z","shell.execute_reply.started":"2023-10-24T16:50:49.399844Z","shell.execute_reply":"2023-10-24T16:50:49.404164Z"}}},{"cell_type":"code","source":"df_stat = df_stat_save.copy()\nencoding = 'onehot'\nstr_model_inf = 'Ridge; '+ encoding + ' for compound only'\nfor score_name in ['mrrmse','r2']:\n    print();print(); print('Score:', score_name ); print()\n    for CV_scheme in list_CV_scheme:\n        list_selected_cols =[ col for col in df_stat.columns if (score_name in col) and (CV_scheme in col) ]\n        col_average_name = score_name + ' ' + CV_scheme   + ' ' +  'average'\n        df_stat[col_average_name] = df_stat[ list_selected_cols ].mean(axis = 1)\n        list_selected_cols = [ col_average_name ] + list_selected_cols\n        plt.figure(figsize = (20,4))\n        \n        str_inf = 'Best alphas: '\n        for col in list_selected_cols:\n            v = df_stat[col]\n            plt.plot(np.log10(list_alpha), df_stat[col] ,'*-' , label = col)\n            #print('Best alpha:', df_stat.sort_values(col)['alpha'].iat[0], 'for ', col )\n            if score_name in ['mrrmse']:\n                str_score = str( np.round(df_stat.sort_values(col)[col].iat[0],4) )\n                str_inf += str(df_stat.sort_values(col)['alpha'].iat[0]) + ' for ' + col + ' ' + str_score +' ; '\n            else:\n                str_score = str( np.round(df_stat.sort_values(col,ascending = False)[col].iat[0],4) )\n                str_inf += str(df_stat.sort_values(col,ascending = False)['alpha'].iat[0]) + ' for ' + col  + ' ' + str_score + ' ; '\n                \n        print(str_inf)\n        plt.xlabel('LOG10 alpha',fontsize = 20 )\n        plt.legend(fontsize = 12 )\n        plt.title('scores: '  + score_name + ', CV: ' + CV_scheme + ', Model: ' + str_model_inf, fontsize = 20 )\n        plt.grid()\n        plt.show()\n\n\n# df_stat","metadata":{"execution":{"iopub.status.busy":"2023-11-11T21:43:59.236967Z","iopub.execute_input":"2023-11-11T21:43:59.237423Z","iopub.status.idle":"2023-11-11T21:44:02.649099Z","shell.execute_reply.started":"2023-11-11T21:43:59.237379Z","shell.execute_reply":"2023-11-11T21:44:02.647822Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\npd.set_option('display.max_columns', 500)\npd.set_option('display.max_rows', 500)\n\nprint( df_stat.shape )\nlist_cols1 = [t for t in df_stat.columns if ('average' in t) or ('alpha' == t) ]\nlist_cols2 = [t for t in df_stat.columns if 'average' not in t ]\ndf_stat = df_stat[list_cols1+list_cols2]\ndisplay( df_stat[list_cols1+list_cols2] )\ndf_stat.to_csv('df_stat.csv')","metadata":{"execution":{"iopub.status.busy":"2023-11-11T21:44:02.650711Z","iopub.execute_input":"2023-11-11T21:44:02.651073Z","iopub.status.idle":"2023-11-11T21:44:02.752298Z","shell.execute_reply.started":"2023-11-11T21:44:02.651042Z","shell.execute_reply":"2023-11-11T21:44:02.751135Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"m = (df_stat['alpha'].iloc[:,0] == 0.1) | (df_stat['alpha'].iloc[:,0] == 0.2) | ( (df_stat['alpha'].iloc[:,0] == 1) )  | ( (df_stat['alpha'].iloc[:,0] == .5) )   | ( (df_stat['alpha'].iloc[:,0] == 1) )\ndf_stat[m].round(4)","metadata":{"execution":{"iopub.status.busy":"2023-11-11T21:44:02.753947Z","iopub.execute_input":"2023-11-11T21:44:02.754264Z","iopub.status.idle":"2023-11-11T21:44:02.814030Z","shell.execute_reply.started":"2023-11-11T21:44:02.754237Z","shell.execute_reply":"2023-11-11T21:44:02.813023Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"pd.set_option('display.max_columns', 50)\npd.set_option('display.max_rows', 100)\n","metadata":{"execution":{"iopub.status.busy":"2023-11-11T21:44:02.815330Z","iopub.execute_input":"2023-11-11T21:44:02.815652Z","iopub.status.idle":"2023-11-11T21:44:02.821206Z","shell.execute_reply.started":"2023-11-11T21:44:02.815624Z","shell.execute_reply":"2023-11-11T21:44:02.819922Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Prepare submit\n\n\nAlpha values around 1 are optimal for mrrmse diffent  CV. Prepare submit with that param.","metadata":{}},{"cell_type":"code","source":"%%time \nfrom sklearn.linear_model import Ridge\nfrom sklearn.decomposition import TruncatedSVD\nfrom sklearn.metrics import r2_score\n\nverbose = 1\n\nn_components = 30\nreducer = TruncatedSVD(n_components=n_components, n_iter=7, random_state=42)\n\nalpha = 1\nmodel = Ridge(alpha=alpha)\n\nY_full_red = reducer.fit_transform(Y_full ) \nmodel.fit(X_full, Y_full_red)\nY_submit_pred_red = model.predict( X_submit )\nY_submit_pred = reducer.inverse_transform( Y_submit_pred_red )\n    \n    \ndf_submit = pd.DataFrame(Y_submit_pred, columns = df_de_train.columns[5:])\ndf_submit.index.name = 'id'\nprint( df_submit.shape )\ndisplay(df_submit)\nfn_save = 'submission_' + 'RidgeAlpha'+str(alpha)+'_tsvd'+str(n_components)+'.csv' \nprint(fn_save)\ndf_submit.to_csv(fn_save)\n","metadata":{"execution":{"iopub.status.busy":"2023-11-11T21:44:02.822329Z","iopub.execute_input":"2023-11-11T21:44:02.822713Z","iopub.status.idle":"2023-11-11T21:44:17.838825Z","shell.execute_reply.started":"2023-11-11T21:44:02.822681Z","shell.execute_reply":"2023-11-11T21:44:17.837923Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Blend with \"priors\"\n\nThe trick typically uplifts the score about 0.01\n\nTrick to add \"priors\" (aggregates over the cell type and compound - which typically improves scores by 0.01 ): ( df_submit = (1-w_compound - w_cell_type)(df_submit) + w_compound df_submit_aggr_compound + w_cell_type * df_submit_aggr_cell_type )\n\nLIUDA CHELDIEVA: https://www.kaggle.com/code/liudacheldieva/streamlined-baseline-approach\n\nZXMKCD : https://www.kaggle.com/code/zmcxjt/streamlined-baseline-approach","metadata":{}},{"cell_type":"code","source":"%%time \n#group by drug and take the mean\ndf_tmp = df_de_train.iloc[:, [1] + list(range(5, df_de_train.shape[1]))] # Take only numeric columns and \"sm_name\"\ndf_aggr = df_tmp.groupby('sm_name').mean().reset_index()\nprint(df_aggr.shape)\ndf_submit_aggr_compound = pd.merge( df_id_map.reset_index(),  df_aggr, on='sm_name', how = 'left' ).sort_values('id').drop(columns = ['cell_type', 'sm_name']).set_index('id')\nprint(df_submit_aggr_compound.shape)\ndisplay(df_submit_aggr_compound.head(3))\n\nprint( )\n\ndf_tmp = df_de_train.iloc[:, [0] + list(range(5, df_de_train.shape[1]))] # Take only numeric columns and \"cell_type\"\ndf_aggr = df_tmp.groupby('cell_type').mean().reset_index()\nprint(df_aggr.shape)\ndf_submit_aggr_cell_type = pd.merge( df_id_map.reset_index(),  df_aggr, on='cell_type', how = 'left' ).sort_values('id').drop(columns = ['cell_type', 'sm_name']).set_index('id')\nprint(df_submit_aggr_cell_type.shape)\ndisplay(df_submit_aggr_cell_type.head(3))","metadata":{"execution":{"iopub.status.busy":"2023-11-11T21:44:17.840175Z","iopub.execute_input":"2023-11-11T21:44:17.840729Z","iopub.status.idle":"2023-11-11T21:44:18.649671Z","shell.execute_reply.started":"2023-11-11T21:44:17.840695Z","shell.execute_reply":"2023-11-11T21:44:18.648492Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nw_cell_type = 0.1\nw_compound = 0.35\ndf_submit_2 = (1-w_compound - w_cell_type)*(df_submit) + w_compound * df_submit_aggr_compound +  w_cell_type * df_submit_aggr_cell_type\nfn_save2 = 'submission_' + 'RidgeAlpha'+str(alpha)+'_tsvd'+str(n_components)+'_blendWithPriors_'+str(w_compound)+'_'+str( w_cell_type )+'.csv' \ndf_submit_2.to_csv(fn_save2)\nprint(fn_save2)\nprint(df_submit_2.shape)\ndisplay(df_submit_2.head(3) )","metadata":{"execution":{"iopub.status.busy":"2023-11-11T21:44:18.651245Z","iopub.execute_input":"2023-11-11T21:44:18.651584Z","iopub.status.idle":"2023-11-11T21:44:32.108220Z","shell.execute_reply.started":"2023-11-11T21:44:18.651552Z","shell.execute_reply":"2023-11-11T21:44:32.107207Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Final timing","metadata":{}},{"cell_type":"code","source":"print('%.1f seconds passed total '%(time.time()-t0start) )\nprint('%.1f minutes passed total '%( (time.time()-t0start)/60)  )\nprint('%.2f hours passed total '%( (time.time()-t0start)/3600)  )","metadata":{"execution":{"iopub.status.busy":"2023-11-11T21:44:32.109795Z","iopub.execute_input":"2023-11-11T21:44:32.110168Z","iopub.status.idle":"2023-11-11T21:44:32.117049Z","shell.execute_reply.started":"2023-11-11T21:44:32.110136Z","shell.execute_reply":"2023-11-11T21:44:32.115906Z"},"trusted":true},"execution_count":null,"outputs":[]}]}