{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.10.12","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"gpu","dataSources":[{"sourceId":59094,"databundleVersionId":7010844,"sourceType":"competition"}],"dockerImageVersionId":30587,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# What is about \n\nModeling for OP2.\n\nHere are \"last hour\" modeling - use PyBoost to directly predict genes (not tsvd-components + tsvd.inverse_transform). Just direct prediction of 18211 targets.\nTo predict all at once - we need more RAM than 16G Kaggle GPU RAM. \nSo we predict by 1000 targets at once and do it 18 times.\n\nWe had not enough time to tune params for that model. First trials give score 0.594 and quite diverse from tsvd-Pyboost. \nSo it might be promosing component for blend.\n\nReport by U900 team is here: \nhttps://www.kaggle.com/competitions/open-problems-single-cell-perturbations/discussion/460858\n\n\n","metadata":{}},{"cell_type":"markdown","source":"## Preliminaries","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-29T19:47:56.308985Z","iopub.execute_input":"2023-11-29T19:47:56.309400Z","iopub.status.idle":"2023-11-29T19:47:56.319768Z","shell.execute_reply.started":"2023-11-29T19:47:56.309370Z","shell.execute_reply":"2023-11-29T19:47:56.318798Z"},"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 )\n\n\nfeatures_columns = ['cell_type','sm_name'] # \n# selected_features_columns = ['cell_type','sm_name']  # [ 'sm_name'] # What features to consider - others not used \n#     # For OP2 task we have two features - cell-type and compound - both categorical\n#     # Some simple models (like Ridge) better work when one uses only one feature - 'sm_name' (=\"compound\") - so you can try that option. \n#     # But seems pyboost is powerful enough to proper work with the both features.\n#     # Table with results for Ridge can be found here: https://www.kaggle.com/code/alexandervc/op2-target-encoders?scriptVersionId=148576642&cellId=1\n#     # Or here: https://docs.google.com/spreadsheets/d/1APN63PMaWZygVjYimK9Ivt0RvifdAU5JRYkxiDn4szw/edit?usp=sharing (sheet \"Encoders\") and discussed here:\n    \nX_submit_categorical = df_id_map[features_columns] # selected_features_columns] # pd.DataFrame(df_id_map, columns= features_columns )\nX_full_categorical = df_de_train[features_columns]# selected_features_columns ]\n\nY_full = df_de_train.iloc[:,5:].values\nprint('X_submit_categorical.shape, X_full_categorical.shape ,  Y_full.shape', X_submit_categorical.shape, X_full_categorical.shape ,  Y_full.shape)","metadata":{"execution":{"iopub.status.busy":"2023-11-29T19:51:09.527838Z","iopub.execute_input":"2023-11-29T19:51:09.528186Z","iopub.status.idle":"2023-11-29T19:52:16.416093Z","shell.execute_reply.started":"2023-11-29T19:51:09.528153Z","shell.execute_reply":"2023-11-29T19:52:16.414960Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Custom CV-schemes\n\nSee the notebook: https://www.kaggle.com/code/alexandervc/op2-class-for-custom-cv-schemes#Ridge-CV-scoring.-Simple-example\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\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":"%%time\nimport random \nfrom sklearn.model_selection import KFold\n\n# 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\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","metadata":{"execution":{"iopub.status.busy":"2023-11-29T19:50:20.590516Z","iopub.execute_input":"2023-11-29T19:50:20.591198Z","iopub.status.idle":"2023-11-29T19:50:20.660041Z","shell.execute_reply.started":"2023-11-29T19:50:20.591158Z","shell.execute_reply":"2023-11-29T19:50:20.659199Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Specify Model and params","metadata":{}},{"cell_type":"code","source":"from sklearn.linear_model import Ridge\nfrom sklearn.decomposition import TruncatedSVD\nfrom sklearn.metrics import r2_score\nfrom sklearn.metrics import mean_absolute_error\nimport category_encoders as ce\n\nverbose = 10\n\n\nstr_inf_cfg = ''\n\n\nalpha = 2e5 # 2e5 seems optimum for  LOO  CT&Drug     (seems for various tsvd n_comps - same is true)\nalpha = 5e4 # 5e4 seems optimum for  LOO  drug only   (seems for various tsvd n_comps - same is true)\nmodel = Ridge(alpha=alpha)\nstr_model_id = 'Ridge'+str(alpha)\nstr_inf_cfg += str_model_id\n\nn_components = 300 # 300 seems to be optimum for encoding both CT and Drug\nn_components = 400 # 400 seems to be optimum for encoding both drug only\n# V33 Ridge50000.0_tsvd400_Dr_LOO_SB1_Tonya_SB0      0.960772 0.67219  0.252819 0.343247 \n# V31 Ridge200000.0_tsvd250_Dr_CT_LOO_SB1_Tonya_SB0  0.956498 0.668837 0.262345 0.363502\n\nreducer = TruncatedSVD(n_components=n_components, n_iter=7, random_state=42)\nstr_inf_cfg += '_tsvd'+str(n_components)\n\n\nlist_features_to_encode = ['cell_type','sm_name']# 'cell_type', ['sm_name']  # [ 'sm_name'] # What features to consider - others not used \n#     # For OP2 task we have two features - cell-type and compound - both categorical\n#     # Some simple models (like Ridge) better work when one uses only one feature - 'sm_name' (=\"compound\") - so you can try that option. \n#     # But seems pyboost is powerful enough to proper work with the both features.\n#     # Table with results for Ridge can be found here: https://www.kaggle.com/code/alexandervc/op2-target-encoders?scriptVersionId=148576642&cellId=1\n#     # Or here: https://docs.google.com/spreadsheets/d/1APN63PMaWZygVjYimK9Ivt0RvifdAU5JRYkxiDn4szw/edit?usp=sharing (sheet \"Encoders\") and discussed here:\n\nif 'sm_name' in list_features_to_encode: str_inf_cfg += '_Dr'\nif 'cell_type' in list_features_to_encode: str_inf_cfg += '_CT'\n\n\n# smoothing = 10; enc = ce.TargetEncoder(smoothing = smoothing )\nenc = ce.LeaveOneOutEncoder()\nstr_inf_cfg += '_LOO'\n\nn_selfblend = 1\nstr_inf_cfg += '_SB'+str(n_selfblend)\n\n\nCV_scheme = 'Tonya' #  'MT'\nvalid_size = 0 \nkf = KFold_custom(CV_scheme, valid_size = valid_size)   \nprint('CV_scheme:', CV_scheme)\nstr_inf_cfg += '_'+CV_scheme\n\nprint(str_inf_cfg)\n","metadata":{"execution":{"iopub.status.busy":"2023-11-29T19:50:20.661095Z","iopub.execute_input":"2023-11-29T19:50:20.661362Z","iopub.status.idle":"2023-11-29T19:50:20.671601Z","shell.execute_reply.started":"2023-11-29T19:50:20.661338Z","shell.execute_reply":"2023-11-29T19:50:20.670679Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time \n\n# that works for CPU also\n!pip install py-boost\n\nimport os\n# Optional: set the device to run\nos.environ[\"CUDA_DEVICE_ORDER\"] = \"PCI_BUS_ID\"\nos.environ[\"CUDA_VISIBLE_DEVICES\"] = \"0\"\n\nos.makedirs('../data', exist_ok=True)\n\nimport joblib\nfrom sklearn.datasets import make_regression\nimport numpy as np\n\n# simple case - just one class is used\nfrom py_boost import GradientBoosting, TLPredictor, TLCompiledPredictor\nfrom py_boost.cv import CrossValidation\n\n\nfrom sklearn.linear_model import Ridge\nfrom sklearn.decomposition import TruncatedSVD\nfrom sklearn.metrics import r2_score\nfrom sklearn.metrics import mean_absolute_error\nimport category_encoders as ce\n\n\nstr_inf_cfg = ''\n\n\nntrees = 50\nmax_depth = 5#10\nsubsample = 1 # From 0 to 1\ncolsample =  0.35  # From 0 to 1 \nlr = 0.01 #  From 0 to 1 \nmodel = GradientBoosting('mse', ntrees=ntrees, lr=lr,  max_depth=max_depth  , subsample=subsample, colsample=colsample, \n                         min_data_in_leaf=1, min_gain_to_split=0, verbose=100  )           \nstr_model_id = 'Pyboost'\nstr_inf_cfg += str_model_id\n\nn_components = 70\nreducer = TruncatedSVD(n_components=n_components, n_iter=7, random_state=42)\nstr_inf_cfg += '_tsvd'+str(n_components)\n\n\n\nenc = ce.QuantileEncoder(quantile =.8)\nstr_inf_cfg += '_QE80'\n# enc = ce.LeaveOneOutEncoder()\n# str_inf_cfg += '_LOO'\n\n\n\n\nlist_features_to_encode = ['cell_type','sm_name']# ['sm_name']  # [ 'sm_name'] # What features to consider - others not used \n#     # For OP2 task we have two features - cell-type and compound - both categorical\n#     # Some simple models (like Ridge) better work when one uses only one feature - 'sm_name' (=\"compound\") - so you can try that option. \n#     # But seems pyboost is powerful enough to proper work with the both features.\n#     # Table with results for Ridge can be found here: https://www.kaggle.com/code/alexandervc/op2-target-encoders?scriptVersionId=148576642&cellId=1\n#     # Or here: https://docs.google.com/spreadsheets/d/1APN63PMaWZygVjYimK9Ivt0RvifdAU5JRYkxiDn4szw/edit?usp=sharing (sheet \"Encoders\") and discussed here:\n\nif 'sm_name' in list_features_to_encode: str_inf_cfg += '_Dr'\nif 'cell_type' in list_features_to_encode: str_inf_cfg += '_CT'\n\n\n\nn_selfblend = 1\nstr_inf_cfg += '_SB'+str(n_selfblend)\n\n\nCV_scheme = 'Tonya' #  'Random' #   'MT'\nvalid_size = 0 \nkf = KFold_custom(CV_scheme, valid_size = valid_size)   \nprint('CV_scheme:', CV_scheme)\nstr_inf_cfg += '_'+CV_scheme\n\nprint(str_inf_cfg)\n\n    \n    ","metadata":{"execution":{"iopub.status.busy":"2023-11-29T19:50:20.673602Z","iopub.execute_input":"2023-11-29T19:50:20.673903Z","iopub.status.idle":"2023-11-29T19:50:32.139457Z","shell.execute_reply.started":"2023-11-29T19:50:20.673879Z","shell.execute_reply":"2023-11-29T19:50:32.138277Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Modeling","metadata":{}},{"cell_type":"code","source":"import gc","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nprint(str_inf_cfg)\n\ndf_stat = pd.DataFrame()\nY_pred_oof_blend = np.zeros( (614, 18211) ); \nY_submit = np.zeros( (255, 18211) ); i_blend_submit = 0\nfor i_selfblend in range(n_selfblend ):\n    t0_start_folds_loop = time.time()\n    list_mrrmse = [];list_r2 = []; list_mae = [];  list_corr_row = []; list_corr_col = [];  list_n_train_sample = []\n    for i_fold,(IX_train,IX_test)   in enumerate( kf.split()) :\n        if i_fold > 1: break ####################  \n            \n        t0_fold = time.time()\n        # \n        # Get train, test:\n        if hasattr(X_full_categorical ,'values'):\n            X_train_categorical,X_test_categorical = X_full_categorical[list_features_to_encode].values[IX_train], X_full_categorical[list_features_to_encode].values[IX_test]\n        else:\n            X_train_categorical,X_test_categorical = X_full_categorical[list_features_to_encode][IX_train], X_full_categorical[list_features_to_encode][IX_test]\n        if hasattr(X_submit_categorical,'values'):\n            X_submit_categorical_loc = X_submit_categorical[list_features_to_encode].values\n        else:\n            X_submit_categorical_loc = X_submit_categorical[list_features_to_encode]\n\n        # Reduce the targets \n        Y_train = Y_full[IX_train]\n        Y_train_red = reducer.fit_transform(Y_full[IX_train] ) \n\n        if verbose >= 100:\n            print('fold:', i_fold, 'X_train_categorical.shape, X_test_categorical.shape, Y_train_red.shape', X_train_categorical.shape, X_test_categorical.shape, Y_train_red.shape)\n\n        ### -------------- Target Encoding -------------------------------------------------------------------\n        for i_target in range(Y_train_red.shape[1]):\n            if i_target == 0:\n                X_train_encoded = enc.fit_transform(X_train_categorical, Y_train_red[:,i_target])\n                X_test_encoded = enc.transform(X_test_categorical)\n                X_submit_encoded = enc.transform(X_submit_categorical_loc)\n            else:\n                X_encoded_tmp = enc.fit_transform(X_train_categorical, Y_train_red[:,i_target])\n                X_train_encoded = np.concatenate( [X_train_encoded, X_encoded_tmp], axis = 1)\n                X_encoded_tmp = enc.transform(X_test_categorical)\n                X_test_encoded = np.concatenate( [X_test_encoded, X_encoded_tmp], axis = 1)\n                X_encoded_tmp = enc.transform(X_submit_categorical_loc)\n                X_submit_encoded = np.concatenate( [X_submit_encoded, X_encoded_tmp], axis = 1)\n\n        list_n_train_sample.append( X_train_encoded.shape[0] ) # Save to inf \n        \n        n_targ = 500\n        for i_targ in range(n_targ,19000,n_targ):\n            t0_target_bucket = time.time()\n            \n            ntrees = 500\n            max_depth = 6\n            subsample = 1 # From 0 to 1\n            colsample =  0.35  # From 0 to 1 \n            lr = 0.01 #  From 0 to 1 \n            model = GradientBoosting('mse', ntrees = ntrees, lr=lr,  max_depth=max_depth  , subsample=subsample, colsample=colsample, \n                                     min_data_in_leaf=1, min_gain_to_split=0, verbose=100  )  \n            \n            # Train model\n            model.fit(X_train_encoded, Y_train[:, i_targ-n_targ : i_targ] )\n            Y_test_pred = model.predict( X_test_encoded )\n            #Y_test_pred = reducer.inverse_transform( Y_test_pred_red )\n\n            # Prepare OOF: \n            #Y_pred_oof_blend[ IX_test ] =  Y_test_pred\n            Y_pred_oof_blend[ IX_test , i_targ-n_targ : i_targ ] = ( Y_pred_oof_blend[ IX_test  , i_targ-n_targ : i_targ ] *  i_selfblend + Y_test_pred ) / ( i_selfblend + 1)\n            # Prapare submit\n            Y_submit_pred = model.predict( X_submit_encoded )\n            #Y_submit_pred = reducer.inverse_transform( Y_submit_pred_red )\n            Y_submit[:, i_targ-n_targ : i_targ] = ( Y_submit[:, i_targ-n_targ : i_targ] *  i_blend_submit  +  Y_submit_pred ) / ( i_blend_submit  + 1  )\n            \n            print('Target buckets i_targ', i_targ, 'time:', np.round( time.time() - t0_target_bucket ) )\n            gc.collect()\n\n        # Scoring  ( score blended OOF)\n        Y_test = Y_full[IX_test]\n        Y_test_pred = Y_pred_oof_blend[ IX_test ]\n        mrrmse_one_fold = np.sqrt(np.square(Y_test - Y_test_pred).mean(axis=1)).mean();  r2_one_fold = r2_score( Y_test , Y_test_pred ) \n        mae_one_fold = mean_absolute_error( Y_test , Y_test_pred  )\n        list_mrrmse.append(mrrmse_one_fold); list_r2.append(r2_one_fold)\n        list_mae.append( mae_one_fold )\n        list_tmp =[ np.corrcoef( Y_test[k,:] , Y_test_pred[k,:]  )[0,1] for k in range(Y_test.shape[0]) ]\n        list_corr_row.append( np.mean(list_tmp ))\n        list_tmp_col =[ np.corrcoef( Y_test[:,k] , Y_test_pred[:,k] )[0,1] for k in range(Y_test.shape[1]) ]\n        list_corr_col.append( np.mean(list_tmp_col ))\n        if verbose >= 10:\n            print('fold:', i_fold, 'mrrmse:', np.round(mrrmse_one_fold,4), 'r2:', np.round(r2_one_fold,4), 'mae:', np.round(mae_one_fold,4), \n                  'corr row:', np.round(list_corr_row[-1],4), 'corr col:', np.round(list_corr_col[-1],4), \n                  'Time', np.round(time.time()- t0_fold,1) )\n\n\n\n    IX = len(df_stat)+1\n    df_stat.loc[IX,'model'] = str_inf_cfg +'_SB'+str(i_selfblend) #  str_model_id # 'NLPr_SB'+str(i_selfblend+1)+ '_MUL'+str(n_multiplex_train) + '_Epo'+str(epochs) + '_Noi'+str(prm_GaussianNoise) +\\\n    #     '_RBW' +str(int(restore_best_weights))\n    df_stat.loc[IX,'mrrmse'] = np.mean(list_mrrmse)\n    df_stat.loc[IX,'mae'] = np.mean(list_mae)\n    df_stat.loc[IX,'corr row'] = np.mean(list_corr_row)\n    df_stat.loc[IX,'corr col'] = np.mean(list_corr_col)\n    df_stat.loc[IX,'corr row med'] = np.median(list_corr_row)\n    df_stat.loc[IX,'corr col med'] = np.median(list_corr_col)\n    df_stat.loc[IX,'r2'] = np.mean(list_r2)\n    df_stat.loc[IX,'time'] = np.round(  time.time()- t0_start_folds_loop  ,1 )\n    for kk in range( len( list_mrrmse ) ):\n        df_stat.loc[IX,'mrrmse ' + str(kk)] = list_mrrmse[kk]\n    for kk in range( len( list_mae ) ):\n        df_stat.loc[IX,'mae ' + str(kk)] = list_mae[kk]\n    for kk in range( len( list_corr_row ) ):\n        df_stat.loc[IX,'corr row ' + str(kk)] = list_corr_row[kk]\n    for kk in range( len( list_corr_col ) ):\n        df_stat.loc[IX,'corr col ' + str(kk)] = list_corr_col[kk]    \n    df_stat.loc[IX,'n_feat'] = X_train_encoded.shape[1]\n    df_stat.loc[IX,'n_train'] = np.median( list_n_train_sample )\n\n\n    if verbose >= 1:\n        print('Average mrrmse:', np.round(np.mean(list_mrrmse ),4 ),  'mae:', np.round(np.mean(list_mae ),4 ), \n              'corr row:',  np.round(np.mean(list_corr_row),4),  'corr col:',  np.round(np.mean(list_corr_col),4) )\n        print('Average r2:', np.round(np.mean(list_r2 ),4 ), 'Time:', np.round(time.time() - t0_start_folds_loop, 1))   \n        print()\n\ndf_stat.round(6).to_csv('df_stat.csv')\npd.set_option('display.max_columns', None)\npd.set_option('display.max_rows', None)\ndisplay(df_stat.round(3) )","metadata":{"execution":{"iopub.status.busy":"2023-11-29T19:50:32.141393Z","iopub.execute_input":"2023-11-29T19:50:32.141757Z","iopub.status.idle":"2023-11-29T19:51:09.490310Z","shell.execute_reply.started":"2023-11-29T19:50:32.141712Z","shell.execute_reply":"2023-11-29T19:51:09.489528Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Show df_stat","metadata":{}},{"cell_type":"code","source":"%%time\npd.set_option('display.max_columns', None)\npd.set_option('display.max_rows', None)\ndisplay(df_stat  )\n# display(df_stat.sort_values(df_stat.columns[0]  ) )","metadata":{"execution":{"iopub.status.busy":"2023-11-29T19:51:09.491586Z","iopub.execute_input":"2023-11-29T19:51:09.491966Z","iopub.status.idle":"2023-11-29T19:51:09.500564Z","shell.execute_reply.started":"2023-11-29T19:51:09.491929Z","shell.execute_reply":"2023-11-29T19:51:09.499606Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Save OOF and submit ","metadata":{}},{"cell_type":"code","source":"%%time\nif 1:\n    str_inf4file = str_inf_cfg\n    fn4save = 'Y_oof_Index_from1_'+str_inf4file \n    fn4save\n    np.save( fn4save+'.npy', Y_pred_oof_blend )\n    df_oof_for_save =  pd.DataFrame( Y_pred_oof_blend )\n    df_oof_for_save.index = range(1,615)\n    df_oof_for_save.index.name = 'id'\n    df_oof_for_save.to_csv(  fn4save + '.csv')\n    print( df_oof_for_save.shape )\n\n\n    df_submit = pd.DataFrame(Y_submit, columns = df_de_train.columns[5:])\n    df_submit.index.name = 'id'\n    print( df_submit.shape )\n    display(df_submit)\n    fn4save = 'Y_submit_'+str_inf4file \n    print( fn4save )\n    np.save( fn4save+'.npy', Y_submit )\n    df_submit.to_csv(  fn4save + '.csv')\n","metadata":{"execution":{"iopub.status.busy":"2023-11-29T19:51:09.502227Z","iopub.execute_input":"2023-11-29T19:51:09.502557Z","iopub.status.idle":"2023-11-29T19:51:09.512311Z","shell.execute_reply.started":"2023-11-29T19:51:09.502524Z","shell.execute_reply":"2023-11-29T19:51:09.511455Z"},"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-29T19:51:09.513728Z","iopub.execute_input":"2023-11-29T19:51:09.514284Z","iopub.status.idle":"2023-11-29T19:51:09.525186Z","shell.execute_reply.started":"2023-11-29T19:51:09.514251Z","shell.execute_reply":"2023-11-29T19:51:09.524375Z"},"trusted":true},"execution_count":null,"outputs":[]}]}