{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.10.14","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"tpu1vmV38","dataSources":[{"sourceId":56537,"databundleVersionId":8015876,"sourceType":"competition"},{"sourceId":8543596,"sourceType":"datasetVersion","datasetId":5102115},{"sourceId":8543979,"sourceType":"datasetVersion","datasetId":5016801},{"sourceId":8544497,"sourceType":"datasetVersion","datasetId":5016941},{"sourceId":8553406,"sourceType":"datasetVersion","datasetId":5051423},{"sourceId":8562207,"sourceType":"datasetVersion","datasetId":5071761},{"sourceId":8573447,"sourceType":"datasetVersion","datasetId":5079224},{"sourceId":8579041,"sourceType":"datasetVersion","datasetId":5130383},{"sourceId":8580145,"sourceType":"datasetVersion","datasetId":5131018},{"sourceId":8581023,"sourceType":"datasetVersion","datasetId":5131829}],"isInternetEnabled":false,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"import numpy as np\nimport matplotlib.pyplot as plt\nimport pandas as pd\nimport time\nimport gc\n\nprint('done')","metadata":{"execution":{"iopub.status.busy":"2024-06-05T16:30:27.320361Z","iopub.execute_input":"2024-06-05T16:30:27.321068Z","iopub.status.idle":"2024-06-05T16:30:30.357251Z","shell.execute_reply.started":"2024-06-05T16:30:27.321032Z","shell.execute_reply":"2024-06-05T16:30:30.35648Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import tensorflow as tf\nfrom tensorflow import keras\nfrom tensorflow.keras import layers\n\nprint('done')","metadata":{"execution":{"iopub.status.busy":"2024-06-05T16:31:17.191862Z","iopub.execute_input":"2024-06-05T16:31:17.192511Z","iopub.status.idle":"2024-06-05T16:31:47.391153Z","shell.execute_reply.started":"2024-06-05T16:31:17.19247Z","shell.execute_reply":"2024-06-05T16:31:47.390184Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# # detect and init the TPU\ntpu = tf.distribute.cluster_resolver.TPUClusterResolver()\n\n# # instantiate a distribution strategy\ntf.tpu.experimental.initialize_tpu_system(tpu)\ntpu_strategy = tf.distribute.TPUStrategy(tpu)\n\nprint('done')","metadata":{"execution":{"iopub.status.busy":"2024-06-05T16:32:03.193878Z","iopub.execute_input":"2024-06-05T16:32:03.194282Z","iopub.status.idle":"2024-06-05T16:32:11.644282Z","shell.execute_reply.started":"2024-06-05T16:32:03.19425Z","shell.execute_reply":"2024-06-05T16:32:11.643513Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print(tpu_strategy)","metadata":{"execution":{"iopub.status.busy":"2024-06-05T16:32:25.499902Z","iopub.execute_input":"2024-06-05T16:32:25.500337Z","iopub.status.idle":"2024-06-05T16:32:25.505156Z","shell.execute_reply.started":"2024-06-05T16:32:25.500302Z","shell.execute_reply":"2024-06-05T16:32:25.504338Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Data","metadata":{}},{"cell_type":"markdown","source":"**Prepare data**","metadata":{}},{"cell_type":"code","source":"sx = np.r_[ \nnp.arange(0,60), np.arange(60,120), np.arange(120,180), \nnp.arange(191,240), np.arange(240,300), np.arange(300,360), \nnp.array([360,361,362,363,364,365,366,367,368,369,370,371,372,373,374,375]), \nnp.arange(376,436), np.arange(436,463), np.arange(496,523) \n]\nnp.arange(180,192)\nsy = np.r_[ \nnp.arange(0,60), np.arange(60,120), np.arange(148,180), np.arange(180,240),\n#np.arange(180,192), np.arange(193,240),\nnp.arange(240,300), np.arange(300,360), \nnp.array([360, 361, 362, 363, 364, 365, 366, 367]) \n]\n\nsx.shape, sy.shape","metadata":{"execution":{"iopub.status.busy":"2024-06-01T19:33:40.875646Z","iopub.execute_input":"2024-06-01T19:33:40.876058Z","iopub.status.idle":"2024-06-01T19:33:40.887317Z","shell.execute_reply.started":"2024-06-01T19:33:40.876028Z","shell.execute_reply":"2024-06-01T19:33:40.886712Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# columns\n\nq = pd.read_csv('/kaggle/input/leap-atmospheric-physics-ai-climsim/train.csv', nrows=2)\n\ncols = q.columns.values.astype(str)\ncols = cols[1:] # remove ids column\n\n# columns, x/input (selection based on std/mad, uniqe, etc)\ncolsx = cols[sx]\n\n# columns, y/target\ncolsyy0 = cols[556:] # all y are from 556\n\n# selection based on std/mad, uniqe, etc\ncolsy0 = colsyy0[sy] # selection (with 0)\n\n# remove 0-stuff\nq = pd.read_csv('/kaggle/input/leap-atmospheric-physics-ai-climsim/sample_submission.csv', nrows=2)\nq = q.iloc[1].values[1:] # y scale sample submission\nsy2 = np.where(q[sy] != 0)[0] # select sy, and then where scaling not 0\ncolsy = colsy0[sy2]\n\nprint(colsx.shape, '\\t', colsyy0.shape, colsy0.shape, colsy.shape) ","metadata":{"execution":{"iopub.status.busy":"2024-06-01T19:33:43.895535Z","iopub.execute_input":"2024-06-01T19:33:43.895873Z","iopub.status.idle":"2024-06-01T19:33:43.955136Z","shell.execute_reply.started":"2024-06-01T19:33:43.895845Z","shell.execute_reply":"2024-06-01T19:33:43.954448Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"n = 500000\nmx = 479\nmy = 292\n\n#n, mx = dx1n.shape\n#_, my = dy1n.shape\nprint(n, mx, my)","metadata":{"execution":{"iopub.status.busy":"2024-06-01T19:33:47.567558Z","iopub.execute_input":"2024-06-01T19:33:47.567884Z","iopub.status.idle":"2024-06-01T19:33:47.572026Z","shell.execute_reply.started":"2024-06-01T19:33:47.567857Z","shell.execute_reply":"2024-06-01T19:33:47.571375Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# free mem\n\nq=1\n\ngc.collect()","metadata":{"execution":{"iopub.status.busy":"2024-06-01T19:33:49.736376Z","iopub.execute_input":"2024-06-01T19:33:49.737231Z","iopub.status.idle":"2024-06-01T19:33:49.784905Z","shell.execute_reply.started":"2024-06-01T19:33:49.737197Z","shell.execute_reply":"2024-06-01T19:33:49.78428Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Neural network","metadata":{}},{"cell_type":"code","source":"# # detect and init the TPU\n# tpu = tf.distribute.cluster_resolver.TPUClusterResolver()\n\n# # instantiate a distribution strategy\n# tf.tpu.experimental.initialize_tpu_system(tpu)\n# tpu_strategy = tf.distribute.TPUStrategy(tpu)\n\n# print('done')","metadata":{"execution":{"iopub.status.busy":"2024-06-01T19:33:52.789425Z","iopub.execute_input":"2024-06-01T19:33:52.790125Z","iopub.status.idle":"2024-06-01T19:33:53.051873Z","shell.execute_reply.started":"2024-06-01T19:33:52.790089Z","shell.execute_reply":"2024-06-01T19:33:53.050922Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"tf.keras.backend.clear_session()\ntf.random.set_seed(1)\nfrom keras.regularizers import L2\n\nreg=5e-6\n\ninp = keras.Input(shape=(mx,))\ndense = layers.Dense( units=8*64, activation=\"relu\", use_bias=1, kernel_regularizer=L2(reg))\nx = dense(inp)\ndense = layers.Dense( units=6*64, activation=\"relu\", use_bias=1, kernel_regularizer=L2(reg))\nx = dense(x)\ndense = layers.Dense( units=4*64, activation=\"tanh\", use_bias=1, kernel_regularizer=L2(reg))\nx = dense(x)\nout = layers.Dense( units=my, activation=\"linear\", use_bias=0, kernel_regularizer=L2(reg))(x)\n\n# instantiating the model in the strategy scope creates the model on the TPU\n#with tpu_strategy.scope():\nmodela = keras.Model(inputs=inp, outputs=out)\n\nmodela.summary()\n#keras.utils.plot_model(modela, \"arch.png\", 1, 1, dpi=86, show_layer_activations=1)\n\n# compile\n\nmodela.compile(\noptimizer=tf.keras.optimizers.Adam(learning_rate=1e-3),\n\n#loss = clos,\nloss = tf.keras.losses.MeanSquaredError(),\n\n#metrics=[cmet]\n#metrics=[tf.keras.metrics.MeanSquaredError()]\nmetrics = [tf.keras.metrics.R2Score()]\n)\n    \nprint('done')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"class trn_generator(keras.utils.PyDataset):  \n\n    def __init__(self):\n        bla=3\n        #self.x, self.y = x_path, y_path\n        #self.batch_size = batch_size\n    \n    def __len__(self):\n        # Return number of batches.\n        return 80\n    \n    def __getitem__(self, idx):\n        \n        a = np.arange(500000)\n        l = int(500000/10)\n        \n        i = idx//8\n        s = a[i*l:(i+1)*l]    \n\n        j = idx%8\n        dx = np.load('/kaggle/input/2024-04-leap-d'+str(j+1)+'/dx1_'+str(j+1)+'.npy')[s,:]\n        dy = np.load('/kaggle/input/2024-04-leap-d'+str(j+1)+'/dy1_'+str(j+1)+'.npy')[s,:]\n\n        #print(dx.shape, dy.shape, idx, i, j)\n        #time.sleep(2)\n\n        return dx, dy\n    \nprint('done')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"class val_generator(keras.utils.PyDataset):  \n\n    def __init__(self):\n        bla=2\n        \n    def __len__(self):\n        # Return number of batches.\n        return 80\n    \n    def __getitem__(self, idx):\n\n        a = np.arange(250000)\n        l = int(250000/10)\n\n        i = idx//8\n        s = a[i*l:(i+1)*l]    \n    \n        j = idx%8\n        dx = np.load('/kaggle/input/2024-04-leap-d'+str(j+1)+'/dx2_'+str(j+1)+'.npy')[s,:]\n        dy = np.load('/kaggle/input/2024-04-leap-d'+str(j+1)+'/dy2_'+str(j+1)+'.npy')[s,:]\n    \n        #print(dx.shape, dy.shape, idx, i, j)\n        #time.sleep(2)\n\n        return dx, dy\n    \n\nprint('done')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# run\n\nepi = 200 # number of epochs\n\nt1 = time.time()\nh = modela.fit(\nx=trn_generator(), validation_data=val_generator(), \nepochs=epi, verbose=1, steps_per_epoch=80+1, validation_freq=10\n)\nt2 = time.time()\nprint( '%.1f' %((t2-t1)/60), 'min' )","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# plot training history\n\nepval = np.linspace(0, h.epoch[-1], len(h.history['val_loss']))\nmet = list(h.history.keys()) # 4 metrics\n\nfig, ax = plt.subplots(2, figsize=(6,8))\n\n# plot 1st (loss)\n\nax[0].plot(h.epoch, h.history[met[0]], label='train');\nax[0].plot(epval, h.history[met[2]], label='val');\nax[0].set_xlabel('epochs'); ax[0].set_ylabel('loss:  ' + h.model.loss.name);\n#ax[0].set_ylim([0.25,0.35]); ax[0].set_yscale('linear');\nax[0].grid(); ax[0].legend();\n\n#% plot 2nd (metric)\n\nax[1].plot(h.epoch, h.history[met[1]], label='train')\nax[1].plot(epval, h.history[met[3]], label='val')\nax[1].set_xlabel('epochs'); ax[1].set_ylabel('metric:  '+met[1]);\nax[1].grid(); ax[1].legend()\nax[1].set_ylim(-0.1,None)\n\nprint(met)\nprint('done')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**Save model**","metadata":{}},{"cell_type":"code","source":"modela.save('modela.keras')\n\nprint('done')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Results","metadata":{}},{"cell_type":"code","source":"gc.collect()","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"dx3n = np.load('/kaggle/input/2024-04-leap-d1/dx3_1.npy')\ndy3n = np.load('/kaggle/input/2024-04-leap-d1/dy3_1.npy')\n\nfor j in range(2,8):\n    dx3n = np.r_[dx3n, np.load('/kaggle/input/2024-04-leap-d'+str(j)+'/dx3_'+str(j)+'.npy')]\n    dy3n = np.r_[dy3n, np.load('/kaggle/input/2024-04-leap-d'+str(j)+'/dy3_'+str(j)+'.npy')]\n    \nprint(dx3n.shape, dy3n.shape)","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"dy3np = modela.predict(x=dx3n) # dy3n predicted (test set)\n\nprint(dy3np.shape)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"dx3n = 1\n\ngc.collect()","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# dy2np = modela.predict(x=dx2n) # dy2n predicted (validation set)\n\n# dy2np.shape, dy2n.shape","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# dy1np = modela.predict(x=dx1n) # dy1n predicted (train set)\n\n# dy1np.shape, dy1n.shape","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**MSE**","metadata":{}},{"cell_type":"code","source":"# # mean squared error (train)\n# print('train:\\n')\n\n# mse = np.mean((dy1n - dy1np)**2)\n# print('mse (train): ', '%.4f' %mse)\n# mses = np.mean((dy1n - dy1np)**2, axis=0)\n# plt.hist(mses, 50, label=str(my)+' mse(s)'); \n# plt.title('train'); plt.legend(); plt.xlabel('mse'); plt.show()\n# for i,j in zip(colsy, mses): print(i, '\\t', '%.4f' %j)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# # mean squared error (validation)\n# print('val:\\n')\n\n# mse = np.mean((dy2n - dy2np)**2)\n# print('mse (val): ', '%.4f' %mse)\n# mses = np.mean((dy2n - dy2np)**2, axis=0)\n# plt.hist(mses, 50, label=str(my)+' mse(s)'); \n# plt.title('validation'); plt.legend(); plt.xlabel('mse'); plt.show()\n# for i,j in zip(colsy, mses): print(i, '\\t', '%.4f' %j)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# mean squared error\nprint('test:\\n')\n\nmse = np.mean((dy3n - dy3np)**2)\nprint('mse: ', '%.4f' %mse)\n\n# mean squared error(s)\nmses = np.mean((dy3n - dy3np)**2, axis=0)\nprint( 'mse(s) (main): ', '%.4f' %(mses.mean()) )\nplt.hist(mses, 50, label=str(my)+' mse(s)'); \nplt.title('test'); plt.legend(); plt.xlabel('mse'); plt.show()\nfor i,j in zip(colsy, mses): print(i, '\\t', '%.4f' %j)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Results (R2)","metadata":{}},{"cell_type":"code","source":"from sklearn.metrics import r2_score\n\nprint('done')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# # R2 error (train)\n# print('train:\\n')\n\n# r2e = 1 - np.sum((dy1np - dy1n)**2) / np.sum((dy1n - dy1n.mean())**2)\n# print('r2e: ', '%.4f' %r2e)\n\n# r2es = np.zeros(my)\n# for i in range(my):\n#     r2es[i] = 1 - np.sum((dy1np[:,i] - dy1n[:,i])**2) / np.sum((dy1n[:,i] - dy1n[:,i].mean())**2)\n# print('r2e(s) (mean): ', '%.4f' %np.mean(r2es))\n# print('r2e (sklearn) :', '%.4f' %(r2_score(dy1n, dy1np)) )\n# plt.hist(r2es, 50, label=str(my)+' r2e(s)'); \n# plt.title('train'); plt.legend(); plt.xlabel('r2e(s)'); plt.show()\n# for i,j in zip(colsy, r2es): print(i, '\\t', '%.4f' %j)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# # R2 error (validation)\n# print('val:\\n')\n\n# r2e = 1 - np.sum((dy2np - dy2n)**2) / np.sum((dy2n - dy2n.mean())**2)\n# print('r2e: ', '%.4f' %r2e)\n\n# r2es = np.zeros(my)\n# for i in range(my):\n#     r2es[i] = 1 - np.sum((dy2np[:,i] - dy2n[:,i])**2) / np.sum((dy2n[:,i] - dy2n[:,i].mean())**2)\n# print('r2e(s) (mean): ', '%.4f' %np.mean(r2es))\n# print('r2e (sklearn) :', '%.4f' %(r2_score(dy2n, dy2np)) )\n# plt.hist(r2es, 50, label=str(my)+' r2e(s)'); \n# plt.title('validation'); plt.legend(); plt.xlabel('r2e(s)'); plt.show()\n# for i,j in zip(colsy, r2es): print(i, '\\t', '%.4f' %j)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# R2 error (test)\nprint('test:\\n')\n\nr2e = 1 - np.sum((dy3np - dy3n)**2) / np.sum((dy3n - dy3n.mean())**2)\nprint('r2e: ', '%.4f' %r2e)\n\nr2es = np.zeros(my)\nfor i in range(my):\n    r2es[i] = 1 - np.sum((dy3np[:,i] - dy3n[:,i])**2) / np.sum((dy3n[:,i] - dy3n[:,i].mean())**2)\nprint('r2e(s) (mean): ', '%.4f' %np.mean(r2es))\nprint('r2e (sklearn) :', '%.4f' %(r2_score(dy3n, dy3np)) )\nplt.hist(r2es, 50, label=str(my)+' r2e(s)'); \nplt.title('test'); plt.legend(); plt.xlabel('r2e(s)'); plt.show()\nfor i,j in zip(colsy, r2es): print(i, '\\t', '%.4f' %j)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**unnorm**","metadata":{}},{"cell_type":"code","source":"# normy = np.load('/kaggle/input/2024-04-leap-norm-x10/normy.npy')\n\n# dy3np = dy3np*normy[1,:] + normy[0,:]\n# dy3n = dy3n*normy[1,:] + normy[0,:]\n\n# # R2 error (test un)\n# print('test un:\\n')\n\n# r2e = 1 - np.sum((dy3np - dy3n)**2) / np.sum((dy3n - dy3n.mean())**2)\n# print('r2e: ', '%.4f' %r2e)\n\n# r2es = np.zeros(my)\n# for i in range(my):\n#     r2es[i] = 1 - np.sum((dy3np[:,i] - dy3n[:,i])**2) / np.sum((dy3n[:,i] - dy3n[:,i].mean())**2)\n# print('r2e(s) (mean): ', '%.4f' %np.mean(r2es))\n# print('r2e (sklearn) :', '%.4f' %(r2_score(dy3n, dy3np)) )\n# plt.hist(r2es, 50, label=str(my)+' r2e(s)'); \n# plt.title('test un'); plt.legend(); plt.xlabel('r2e(s)'); plt.show()\n# for i,j in zip(colsy, r2es): print(i, '\\t', '%.4f' %j)","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# The end","metadata":{}}]}