{"metadata":{"kernelspec":{"display_name":"Python 3","language":"python","name":"python3"},"language_info":{"codemirror_mode":{"name":"ipython","version":3},"file_extension":".py","mimetype":"text/x-python","name":"python","nbconvert_exporter":"python","pygments_lexer":"ipython3","version":"3.10.13"},"kaggle":{"accelerator":"nvidiaTeslaT4","dataSources":[{"sourceId":70367,"databundleVersionId":9188054,"sourceType":"competition"}],"dockerImageVersionId":30746,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":true},"papermill":{"default_parameters":{},"duration":1087.424982,"end_time":"2024-08-03T13:00:18.164826","environment_variables":{},"exception":null,"input_path":"__notebook__.ipynb","output_path":"__notebook__.ipynb","parameters":{},"start_time":"2024-08-03T12:42:10.739844","version":"2.5.0"}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# Ariel Data Challenge 2024: Introductory model: training\n\nIn this notebook, we show how to train and cross-validate a model. At the end, we save the model so that it can be used for inference in the separate notebook [ADC24 Intro inference](https://www.kaggle.com/code/ambrosm/adc24-intro-inference).\n","metadata":{"_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","papermill":{"duration":0.009657,"end_time":"2024-08-03T12:42:14.074227","exception":false,"start_time":"2024-08-03T12:42:14.06457","status":"completed"},"tags":[]}},{"cell_type":"code","source":"!pip install lightautoml","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-11-30T10:35:21.513433Z","iopub.execute_input":"2024-11-30T10:35:21.514088Z","iopub.status.idle":"2024-11-30T10:37:18.6065Z","shell.execute_reply.started":"2024-11-30T10:35:21.51405Z","shell.execute_reply":"2024-11-30T10:37:18.60496Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"!pip install dask[dataframe]","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-11-30T09:48:05.3687Z","iopub.execute_input":"2024-11-30T09:48:05.369368Z","iopub.status.idle":"2024-11-30T09:48:14.975199Z","shell.execute_reply.started":"2024-11-30T09:48:05.369328Z","shell.execute_reply":"2024-11-30T09:48:14.974365Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"from lightautoml.automl.presets.tabular_presets import TabularAutoML\nfrom lightautoml.tasks import Task\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-11-30T10:38:18.846803Z","iopub.execute_input":"2024-11-30T10:38:18.847139Z","iopub.status.idle":"2024-11-30T10:38:18.851481Z","shell.execute_reply.started":"2024-11-30T10:38:18.847113Z","shell.execute_reply":"2024-11-30T10:38:18.850536Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import pandas as pd\nimport polars as pl\nimport matplotlib.pyplot as plt\nimport numpy as np\nimport seaborn as sns\nimport scipy.stats\nfrom tqdm import tqdm\nimport pickle\n\nfrom sklearn.model_selection import cross_val_predict\nfrom sklearn.linear_model import Ridge\nfrom sklearn.metrics import r2_score, mean_squared_error","metadata":{"_kg_hide-input":true,"papermill":{"duration":3.023083,"end_time":"2024-08-03T12:42:17.107326","exception":false,"start_time":"2024-08-03T12:42:14.084243","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2024-11-30T09:48:14.976659Z","iopub.execute_input":"2024-11-30T09:48:14.976949Z","iopub.status.idle":"2024-11-30T09:48:16.669095Z","shell.execute_reply.started":"2024-11-30T09:48:14.97692Z","shell.execute_reply":"2024-11-30T09:48:16.668459Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# A look at the data\n\nWe start by reading the metadata:","metadata":{"papermill":{"duration":0.008428,"end_time":"2024-08-03T12:42:17.180297","exception":false,"start_time":"2024-08-03T12:42:17.171869","status":"completed"},"tags":[]}},{"cell_type":"code","source":"train_adc_info = pd.read_csv('/kaggle/input/ariel-data-challenge-2024/train_adc_info.csv',\n                           index_col='planet_id')\n# test_adc_info = pd.read_csv('/kaggle/input/ariel-data-challenge-2024/test_adc_info.csv',\n#                            index_col='planet_id')\ntrain_labels = pd.read_csv('/kaggle/input/ariel-data-challenge-2024/train_labels.csv',\n                           index_col='planet_id')\nwavelengths = pd.read_csv('/kaggle/input/ariel-data-challenge-2024/wavelengths.csv')\naxis_info = pd.read_parquet('/kaggle/input/ariel-data-challenge-2024/axis_info.parquet')\n","metadata":{"papermill":{"duration":0.188106,"end_time":"2024-08-03T12:42:17.37713","exception":false,"start_time":"2024-08-03T12:42:17.189024","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2024-11-30T09:48:16.670715Z","iopub.execute_input":"2024-11-30T09:48:16.671036Z","iopub.status.idle":"2024-11-30T09:48:16.988186Z","shell.execute_reply.started":"2024-11-30T09:48:16.671014Z","shell.execute_reply":"2024-11-30T09:48:16.987506Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"Some facts about the data:\n- We have 673 planets for training. These planets belong to two different stars.\n- There will be roughly 800 planets for testing (but the test data is hidden).\n- The competition is a multi-output regression task with 283 targets to predict.","metadata":{"papermill":{"duration":0.008486,"end_time":"2024-08-03T12:42:17.394634","exception":false,"start_time":"2024-08-03T12:42:17.386148","status":"completed"},"tags":[]}},{"cell_type":"markdown","source":"## The FGS1 data\n\nHaving read the metadata, we'll tackle the FGS1 data (Fine Guidance System). The FGS1 measurements consist of one file per planet (673 files for 673 planets for training). For now, we ignore the calibration files.\n\nEach file contains 135,000 rows of images taken at 0.1 second time steps. Each row is a 32\\*32 image at a single wavelength.\n\nWe read a sample file:","metadata":{"papermill":{"duration":0.008446,"end_time":"2024-08-03T12:42:17.411955","exception":false,"start_time":"2024-08-03T12:42:17.403509","status":"completed"},"tags":[]}},{"cell_type":"code","source":"planet_id = 14485303\nf_signal = pd.read_parquet(f'/kaggle/input/ariel-data-challenge-2024/train/{planet_id}/FGS1_signal.parquet')\nf_signal","metadata":{"papermill":{"duration":1.759363,"end_time":"2024-08-03T12:42:19.180097","exception":false,"start_time":"2024-08-03T12:42:17.420734","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2024-11-30T09:48:16.989162Z","iopub.execute_input":"2024-11-30T09:48:16.989509Z","iopub.status.idle":"2024-11-30T09:48:18.601703Z","shell.execute_reply.started":"2024-11-30T09:48:16.989477Z","shell.execute_reply":"2024-11-30T09:48:18.600702Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"Every row of the file corresponds to an image of a star. The images come in pairs, and the second image is lighter than the first one:","metadata":{"papermill":{"duration":0.009711,"end_time":"2024-08-03T12:42:19.199727","exception":false,"start_time":"2024-08-03T12:42:19.190016","status":"completed"},"tags":[]}},{"cell_type":"code","source":"_, (ax1, ax2) = plt.subplots(1, 2, figsize=(12, 4))\nsns.heatmap(f_signal.iloc[0].values.reshape(32, 32), ax=ax1, vmin=0, vmax=52000)\nax1.set_aspect('equal')\nsns.heatmap(f_signal.iloc[1].values.reshape(32, 32), ax=ax2, vmin=0, vmax=52000)\nax2.set_aspect('equal')\nplt.suptitle('A pair of FGS1 images')\nplt.show()","metadata":{"_kg_hide-input":true,"papermill":{"duration":0.559852,"end_time":"2024-08-03T12:42:19.769794","exception":false,"start_time":"2024-08-03T12:42:19.209942","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2024-11-30T09:48:18.602763Z","iopub.execute_input":"2024-11-30T09:48:18.603023Z","iopub.status.idle":"2024-11-30T09:48:19.304996Z","shell.execute_reply.started":"2024-11-30T09:48:18.603Z","shell.execute_reply":"2024-11-30T09:48:19.30412Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"To see the time series, we first have to compute the difference between the even and the odd frames to get the net signal (67500 time steps). We then take the mean over all 1024 pixels. The net signal is very noisy, and we smoothen it by computing a moving average. The plot of the smoothened signal clearly shows that the signal intensity is reduced (i.e., the image gets darker) while the planet passes in front of the star (between time steps 23500 and 44000).\n\nThe left diagram shows a planet with a strong reduction of the signal intensity, the right diagram shows a planet with a weak reduction:","metadata":{"papermill":{"duration":0.011386,"end_time":"2024-08-03T12:42:19.791847","exception":false,"start_time":"2024-08-03T12:42:19.780461","status":"completed"},"tags":[]}},{"cell_type":"code","source":"_, ((ax1, ax2), (ax3, ax4)) = plt.subplots(2, 2, sharex=True, figsize=(12, 4))\n\nplanet_id = 14485303\nf_signal = pd.read_parquet(f'/kaggle/input/ariel-data-challenge-2024/train/{planet_id}/FGS1_signal.parquet')\n\nmean_signal = f_signal.values.mean(axis=1)\nnet_signal = mean_signal[1::2] - mean_signal[0::2]\ncum_signal = net_signal.cumsum()\nwindow=800\nsmooth_signal = (cum_signal[window:] - cum_signal[:-window]) / window\n\nax1.set_title('FGS1: time series of planet with strong signal')\nax1.plot(net_signal, label='raw signal')\nax1.legend()\nax3.plot(smooth_signal, color='c', label='smoothened signal')\nax3.legend()\nax3.set_xlabel('time step')\nfor time_step in [20500, 23500, 44000, 47000]:\n    ax3.axvline(time_step, color='gray')\n\nplanet_id = 4249337798\nf_signal = pd.read_parquet(f'/kaggle/input/ariel-data-challenge-2024/train/{planet_id}/FGS1_signal.parquet')\n\nmean_signal = f_signal.values.mean(axis=1)\nnet_signal = mean_signal[1::2] - mean_signal[0::2]\ncum_signal = net_signal.cumsum()\nwindow=800\nsmooth_signal = (cum_signal[window:] - cum_signal[:-window]) / window\n\nax2.set_title('FGS1: time series of planet with weak signal')\nax2.plot(net_signal, label='raw signal')\nax2.legend()\nax4.plot(smooth_signal, color='c', label='smoothened signal')\nax4.legend()\nax4.set_xlabel('time step')\nfor time_step in [20500, 23500, 44000, 47000]:\n    ax4.axvline(time_step, color='gray')\n\n# plt.suptitle('FGS1 time series', y=0.96)\nplt.show()","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2024-11-30T09:48:19.305984Z","iopub.execute_input":"2024-11-30T09:48:19.306208Z","iopub.status.idle":"2024-11-30T09:48:22.99438Z","shell.execute_reply.started":"2024-11-30T09:48:19.306189Z","shell.execute_reply":"2024-11-30T09:48:22.993493Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## The AIRS data\n\nAIRS is the other sensor of the satellite. It produces one file per planet as well. Each file contains 11,250 rows of images captured at constant time steps. Each 32 x 356 image has been flattened into 11392 columns.","metadata":{"papermill":{"duration":0.077595,"end_time":"2024-08-03T13:00:12.362301","exception":false,"start_time":"2024-08-03T13:00:12.284706","status":"completed"},"tags":[]}},{"cell_type":"code","source":"planet_id = 4249337798\na_signal = pd.read_parquet(f'/kaggle/input/ariel-data-challenge-2024/train/{planet_id}/AIRS-CH0_signal.parquet')\na_signal","metadata":{"papermill":{"duration":2.313418,"end_time":"2024-08-03T13:00:14.755198","exception":false,"start_time":"2024-08-03T13:00:12.44178","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2024-11-30T09:48:22.995472Z","iopub.execute_input":"2024-11-30T09:48:22.995778Z","iopub.status.idle":"2024-11-30T09:48:24.859225Z","shell.execute_reply.started":"2024-11-30T09:48:22.995753Z","shell.execute_reply":"2024-11-30T09:48:24.857716Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"a_signal = a_signal.values.reshape(11250, 32, 356)\n\nplt.figure(figsize=(10, 3))\nsns.heatmap(a_signal[1])\nplt.ylabel('spatial dimension')\nplt.xlabel('wavelength dimension')\nplt.show()","metadata":{"_kg_hide-input":true,"papermill":{"duration":0.807236,"end_time":"2024-08-03T13:00:15.817499","exception":false,"start_time":"2024-08-03T13:00:15.010263","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2024-11-30T09:48:24.860529Z","iopub.execute_input":"2024-11-30T09:48:24.860889Z","iopub.status.idle":"2024-11-30T09:48:25.392549Z","shell.execute_reply.started":"2024-11-30T09:48:24.860856Z","shell.execute_reply":"2024-11-30T09:48:25.391727Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"The data again is a time series, and we can see how the star is obscured while the planet is passing in front of it.","metadata":{"papermill":{"duration":0.08019,"end_time":"2024-08-03T13:00:15.979865","exception":false,"start_time":"2024-08-03T13:00:15.899675","status":"completed"},"tags":[]}},{"cell_type":"code","source":"mean_signal = a_signal.mean(axis=2).mean(axis=1)\nnet_signal = mean_signal[1::2] - mean_signal[0::2]\ncum_signal = net_signal.cumsum()\nwindow=80\nsmooth_signal = (cum_signal[window:] - cum_signal[:-window]) / window\n\n_, (ax1, ax2) = plt.subplots(2, 1, sharex=True)\nax1.plot(net_signal, label='raw net signal')\nax1.legend()\nax2.plot(smooth_signal, color='c', label='smoothened net signal')\nax2.legend()\nax2.set_xlabel('time')\nfor time_step in [20500, 23500, 44000, 47000]:\n    ax2.axvline(time_step * 11250 // 135000, color='gray')\nplt.suptitle('AIRS-CH0 time series', y=0.96)\nplt.show()\n","metadata":{"_kg_hide-input":true,"papermill":{"duration":0.727242,"end_time":"2024-08-03T13:00:16.788841","exception":false,"start_time":"2024-08-03T13:00:16.061599","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2024-11-30T09:48:25.39512Z","iopub.execute_input":"2024-11-30T09:48:25.395383Z","iopub.status.idle":"2024-11-30T09:48:25.892668Z","shell.execute_reply.started":"2024-11-30T09:48:25.395362Z","shell.execute_reply":"2024-11-30T09:48:25.891931Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Reading the data\n\nWe now read the FGS1 data and the AIRS-CH0 data for all 673 training planets. As the dataset doesn't fit into RAM completely, we keep only two one-dimensional time series for every planet. At the end we'll have\n1. A time series with 67500 steps per planet taken from the FGS1 data, and\n2. A time series with 5625 steps per planet taken from the AIRS-CH0 data.\n\nWe use the Jupyter `%%writefile` cell magic to save the function code to a file. This ensures that the inference notebook will process the test data in exactly the same way as this notebook processes the training data.","metadata":{}},{"cell_type":"code","source":"%%writefile f_read_and_preprocess.py\n\ndef f_read_and_preprocess(dataset, adc_info, planet_ids):\n    \"\"\"Read the FGS1 files for all planet_ids and extract the time series.\n    \n    Parameters\n    dataset: 'train' or 'test'\n    adc_info: metadata dataframe, either train_adc_info or test_adc_info\n    planet_ids: list of planet ids\n    \n    Returns\n    dataframe with one row per planet_id and 67500 values per row\n    \n    \"\"\"\n    f_raw_train = np.full((len(planet_ids), 67500), np.nan, dtype=np.float32)\n    for i, planet_id in tqdm(list(enumerate(planet_ids))):\n        f_signal = pl.read_parquet(f'/kaggle/input/ariel-data-challenge-2024/{dataset}/{planet_id}/FGS1_signal.parquet')\n        mean_signal = f_signal.cast(pl.Int32).sum_horizontal().cast(pl.Float32).to_numpy() / 1024 # mean over the 32*32 pixels\n        net_signal = mean_signal[1::2] - mean_signal[0::2]\n        f_raw_train[i] = net_signal\n    return f_raw_train\n    ","metadata":{"papermill":{"duration":1067.258124,"end_time":"2024-08-03T13:00:08.197814","exception":false,"start_time":"2024-08-03T12:42:20.93969","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2024-11-30T09:48:25.893841Z","iopub.execute_input":"2024-11-30T09:48:25.894223Z","iopub.status.idle":"2024-11-30T09:48:25.900233Z","shell.execute_reply.started":"2024-11-30T09:48:25.894189Z","shell.execute_reply":"2024-11-30T09:48:25.899438Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"%%time\nexec(open('f_read_and_preprocess.py', 'r').read())\nf_raw_train = f_read_and_preprocess('train', train_adc_info, train_labels.index)\nwith open('f_raw_train.pickle', 'wb') as f:\n    pickle.dump(f_raw_train, f)\n","metadata":{"papermill":{"duration":1067.258124,"end_time":"2024-08-03T13:00:08.197814","exception":false,"start_time":"2024-08-03T12:42:20.93969","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2024-11-30T09:48:25.901171Z","iopub.execute_input":"2024-11-30T09:48:25.901392Z","iopub.status.idle":"2024-11-30T09:58:36.611659Z","shell.execute_reply.started":"2024-11-30T09:48:25.901373Z","shell.execute_reply":"2024-11-30T09:58:36.610753Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"%%writefile a_read_and_preprocess.py\ndef a_read_and_preprocess(dataset, adc_info, planet_ids):\n    \"\"\"Read the AIRS-CH0 files for all planet_ids and extract the time series.\n    \n    Parameters\n    dataset: 'train' or 'test'\n    adc_info: metadata dataframe, either train_adc_info or test_adc_info\n    planet_ids: list of planet ids\n    \n    Returns\n    dataframe with one row per planet_id and 5625 values per row\n    \n    \"\"\"\n    a_raw_train = np.full((len(planet_ids), 5625), np.nan, dtype=np.float32)\n    for i, planet_id in tqdm(list(enumerate(planet_ids))):\n        signal = pl.read_parquet(f'/kaggle/input/ariel-data-challenge-2024/{dataset}/{planet_id}/AIRS-CH0_signal.parquet')\n        mean_signal = signal.cast(pl.Int32).sum_horizontal().cast(pl.Float32).to_numpy() / (32*356) # mean over the 32*356 pixels\n        net_signal = mean_signal[1::2] - mean_signal[0::2]\n        a_raw_train[i] = net_signal\n    return a_raw_train\n    ","metadata":{"execution":{"iopub.status.busy":"2024-11-30T09:58:36.612971Z","iopub.execute_input":"2024-11-30T09:58:36.613608Z","iopub.status.idle":"2024-11-30T09:58:36.618987Z","shell.execute_reply.started":"2024-11-30T09:58:36.613574Z","shell.execute_reply":"2024-11-30T09:58:36.618181Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"%%time\nexec(open('a_read_and_preprocess.py', 'r').read())\na_raw_train = a_read_and_preprocess('train', train_adc_info, train_labels.index)\nwith open('a_raw_train.pickle', 'wb') as f:\n    pickle.dump(a_raw_train, f)\n","metadata":{"execution":{"iopub.status.busy":"2024-11-30T10:15:59.102289Z","iopub.execute_input":"2024-11-30T10:15:59.103094Z","iopub.status.idle":"2024-11-30T10:30:35.664458Z","shell.execute_reply.started":"2024-11-30T10:15:59.103068Z","shell.execute_reply":"2024-11-30T10:30:35.663563Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"As a plausibility check, we plot the means of all time series:","metadata":{}},{"cell_type":"code","source":"plt.figure(figsize=(6, 2))\nplt.plot(f_raw_train.mean(axis=0))\nfor time_step in [20500, 23500, 44000, 47000]:\n    plt.axvline(time_step, color='gray')\nplt.xlabel('time step')\nplt.title('FGS1: Overall mean')\nplt.show()\n\nplt.figure(figsize=(6, 2))\nplt.plot(a_raw_train.mean(axis=0))\nfor time_step in [20500, 23500, 44000, 47000]:\n    plt.axvline(time_step * 11250 // 135000, color='gray')\nplt.xlabel('time step')\nplt.title('AIRS-CH0: Overall mean')\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-11-30T10:30:52.862743Z","iopub.execute_input":"2024-11-30T10:30:52.863317Z","iopub.status.idle":"2024-11-30T10:30:53.204111Z","shell.execute_reply.started":"2024-11-30T10:30:52.863288Z","shell.execute_reply":"2024-11-30T10:30:53.203225Z"},"_kg_hide-input":true,"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Feature engineering\n\nWe want to know how much darker the images get when the planet obscures the star. The time series diagrams above show that the planets reduce the brightness of the stars (on average) by 0.2 % (from 228.2 to 227.6 or from 1371 to 1368).","metadata":{}},{"cell_type":"code","source":"%%writefile feature_engineering.py\n\ndef feature_engineering(f_raw, a_raw):\n    \"\"\"Create a dataframe with two features from the raw data.\n    \n    Parameters:\n    f_raw: ndarray of shape (n_planets, 67500)\n    a_raw: ndarray of shape (n_planets, 5625)\n    \n    Return value:\n    df: DataFrame of shape (n_planets, 2)\n    \"\"\"\n    obscured = f_raw[:, 23500:44000].mean(axis=1)\n    unobscured = (f_raw[:, :20500].mean(axis=1) + f_raw[:, 47000:].mean(axis=1)) / 2\n    f_relative_reduction = (unobscured - obscured) / unobscured\n    obscured = a_raw[:, 1958:3666].mean(axis=1)\n    unobscured = (a_raw[:, :1708].mean(axis=1) + a_raw[:, 3916:].mean(axis=1)) / 2\n    a_relative_reduction = (unobscured - obscured) / unobscured\n\n    df = pd.DataFrame({'a_relative_reduction': a_relative_reduction,\n                       'f_relative_reduction': f_relative_reduction})\n    \n    return df\n","metadata":{"execution":{"iopub.status.busy":"2024-11-30T10:30:55.775838Z","iopub.execute_input":"2024-11-30T10:30:55.776498Z","iopub.status.idle":"2024-11-30T10:30:55.782242Z","shell.execute_reply.started":"2024-11-30T10:30:55.776469Z","shell.execute_reply":"2024-11-30T10:30:55.781273Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"exec(open('feature_engineering.py', 'r').read())\n\ntrain = feature_engineering(f_raw_train, a_raw_train)","metadata":{"execution":{"iopub.status.busy":"2024-11-30T10:30:58.917478Z","iopub.execute_input":"2024-11-30T10:30:58.918213Z","iopub.status.idle":"2024-11-30T10:30:58.939969Z","shell.execute_reply.started":"2024-11-30T10:30:58.918185Z","shell.execute_reply":"2024-11-30T10:30:58.939109Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"The following scatterplot shows a strong correlation between the signal reduction when the planet is in front of the star and one of the targets we want to predict. The other targets have a similarly high correlation.","metadata":{"papermill":{"duration":0.073746,"end_time":"2024-08-03T13:00:08.349433","exception":false,"start_time":"2024-08-03T13:00:08.275687","status":"completed"},"tags":[]}},{"cell_type":"code","source":"\ncolor_array = np.array(plt.rcParams['axes.prop_cycle'].by_key()['color'])\nplt.scatter(train.a_relative_reduction, train_labels.wl_1, s=15, alpha=0.5,\n            c=color_array[train_adc_info.star])\nplt.xlabel('relative signal reduction when planet is in front')\nplt.ylabel('target')\nplt.title('Correlation between relative signal reduction and target')\nplt.gca().set_aspect('equal')\npoints = [plt.Line2D([0], [0], label=f'star {i}', marker='o', markersize=3,\n         markeredgecolor=color_array[i], markerfacecolor=color_array[i], linestyle='') for i in range(2)]\n\nplt.legend(handles=points)\nplt.show()","metadata":{"_kg_hide-input":true,"papermill":{"duration":0.478644,"end_time":"2024-08-03T13:00:08.902093","exception":false,"start_time":"2024-08-03T13:00:08.423449","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2024-11-30T10:31:00.884963Z","iopub.execute_input":"2024-11-30T10:31:00.885786Z","iopub.status.idle":"2024-11-30T10:31:01.186325Z","shell.execute_reply.started":"2024-11-30T10:31:00.885754Z","shell.execute_reply":"2024-11-30T10:31:01.185539Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"train","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-11-30T10:31:03.271765Z","iopub.execute_input":"2024-11-30T10:31:03.272377Z","iopub.status.idle":"2024-11-30T10:31:03.281869Z","shell.execute_reply.started":"2024-11-30T10:31:03.272352Z","shell.execute_reply":"2024-11-30T10:31:03.280863Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"train_labels","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-11-30T10:32:07.554922Z","iopub.execute_input":"2024-11-30T10:32:07.555261Z","iopub.status.idle":"2024-11-30T10:32:07.580647Z","shell.execute_reply.started":"2024-11-30T10:32:07.555236Z","shell.execute_reply":"2024-11-30T10:32:07.58Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# The model and the cross-validation\n\nTo keep things simple, we predict the targets with ridge regression.\n\nWe are interested in three cross-validation metrics:\n1. The R2 score is above 0.9, which confirms the correlation we've seen in the scatterplot.\n2. The root mean squared error will be the predicted uncertainty.\n3. The competition metric gives an indication of the leaderboard score. Unfortunately the competition metric depends on the value of `sigma_true`, which I don't know.","metadata":{"papermill":{"duration":0.076791,"end_time":"2024-08-03T13:00:09.056736","exception":false,"start_time":"2024-08-03T13:00:08.979945","status":"completed"},"tags":[]}},{"cell_type":"markdown","source":"## 최후의 수단 AutoML 사용","metadata":{}},{"cell_type":"code","source":"train_data = train.copy()\ntrain_labels_df = train_labels.reset_index(drop=True)\ntrain_data = pd.concat([train_data, train_labels_df], axis=1)\ntrain_data","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-11-30T10:40:08.193359Z","iopub.execute_input":"2024-11-30T10:40:08.194369Z","iopub.status.idle":"2024-11-30T10:40:08.219353Z","shell.execute_reply.started":"2024-11-30T10:40:08.194329Z","shell.execute_reply":"2024-11-30T10:40:08.218536Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"task = Task('multi:reg', loss='mse')\nautoml = TabularAutoML(task=task, timeout=600)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-11-30T11:03:08.233744Z","iopub.execute_input":"2024-11-30T11:03:08.23459Z","iopub.status.idle":"2024-11-30T11:03:08.266753Z","shell.execute_reply.started":"2024-11-30T11:03:08.234564Z","shell.execute_reply":"2024-11-30T11:03:08.265793Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"roles = {'target': train_labels.columns.tolist()}\noof_pred = automl.fit_predict(train_data, roles=roles)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-11-30T11:03:12.495643Z","iopub.execute_input":"2024-11-30T11:03:12.496319Z","iopub.status.idle":"2024-11-30T11:09:31.237533Z","shell.execute_reply.started":"2024-11-30T11:03:12.496289Z","shell.execute_reply":"2024-11-30T11:09:31.236719Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"print(automl.create_model_str_desc())","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-11-30T11:21:19.183142Z","iopub.execute_input":"2024-11-30T11:21:19.18352Z","iopub.status.idle":"2024-11-30T11:21:19.188533Z","shell.execute_reply.started":"2024-11-30T11:21:19.183489Z","shell.execute_reply":"2024-11-30T11:21:19.187634Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"y_true = train_labels.values\ny_pred = oof_pred.data","metadata":{"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"r2_scores = [r2_score(y_true[:, i], y_pred[:, i]) for i in range(y_true.shape[1])]\nrmse_scores = [mean_squared_error(y_true[:, i], y_pred[:, i], squared=False) for i in range(y_true.shape[1])]\n\nprint(f\"# Average R2 score: {np.mean(r2_scores):.3f}\")\nprint(f\"# Average RMSE: {np.mean(rmse_scores):.6f}\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-11-30T11:09:31.247756Z","iopub.execute_input":"2024-11-30T11:09:31.248059Z","iopub.status.idle":"2024-11-30T11:09:31.521786Z","shell.execute_reply.started":"2024-11-30T11:09:31.248031Z","shell.execute_reply":"2024-11-30T11:09:31.520753Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"print(f\"# R2 score: {r2_score(train_labels, oof_pred):.3f}\")\nsigma_pred = mean_squared_error(train_labels, oof_pred, squared=False)\nprint(f\"# Root mean squared error: {sigma_pred:.6f}\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-11-30T11:41:14.759156Z","iopub.execute_input":"2024-11-30T11:41:14.759488Z","iopub.status.idle":"2024-11-30T11:41:14.785751Z","shell.execute_reply.started":"2024-11-30T11:41:14.759462Z","shell.execute_reply":"2024-11-30T11:41:14.784573Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"model = Ridge(alpha=1e-12)\n\noof_pred = cross_val_predict(model, train, train_labels)\n\nprint(f\"# R2 score: {r2_score(train_labels, oof_pred):.3f}\")\nsigma_pred = mean_squared_error(train_labels, oof_pred, squared=False)\nprint(f\"# Root mean squared error: {sigma_pred:.6f}\")\n# R2 score: 0.971\n# Root mean squared error: 0.000293","metadata":{"papermill":{"duration":0.700113,"end_time":"2024-08-03T13:00:09.837053","exception":false,"start_time":"2024-08-03T13:00:09.13694","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2024-11-30T11:41:18.269752Z","iopub.execute_input":"2024-11-30T11:41:18.270445Z","iopub.status.idle":"2024-11-30T11:41:18.335539Z","shell.execute_reply.started":"2024-11-30T11:41:18.270416Z","shell.execute_reply":"2024-11-30T11:41:18.334722Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"oof_pred_array = oof_pred.data\noof_pred = pd.DataFrame(oof_pred_array, index=train_labels.index, columns=train_labels.columns)\noof_pred","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-11-30T11:31:37.413718Z","iopub.execute_input":"2024-11-30T11:31:37.414508Z","iopub.status.idle":"2024-11-30T11:31:37.436562Z","shell.execute_reply.started":"2024-11-30T11:31:37.414481Z","shell.execute_reply":"2024-11-30T11:31:37.435742Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"col = 1  \n\ny_true = train_labels.iloc[:, col].values  # Convert to NumPy array\ny_pred = oof_pred_array[:, col]\n\n# Plot\nplt.scatter(y_true, y_pred, s=15, c='lightgreen')\nplt.gca().set_aspect('equal')\nplt.xlabel('True Values')\nplt.ylabel('Predicted Values')\nplt.title('Predicted vs. True Values')\nplt.show()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2024-11-30T11:38:50.034616Z","iopub.execute_input":"2024-11-30T11:38:50.035443Z","iopub.status.idle":"2024-11-30T11:38:50.292198Z","shell.execute_reply.started":"2024-11-30T11:38:50.035416Z","shell.execute_reply":"2024-11-30T11:38:50.291521Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# col = 1\n# plt.scatter(oof_pred[:,col], train_labels.iloc[:,col], s=15, c='lightgreen')\n# plt.gca().set_aspect('equal')\n# plt.xlabel('y_pred')\n# plt.ylabel('y_true')\n# plt.title('Comparing y_true and y_pred')\n# plt.show()","metadata":{"papermill":{"duration":0.700113,"end_time":"2024-08-03T13:00:09.837053","exception":false,"start_time":"2024-08-03T13:00:09.13694","status":"completed"},"tags":[],"_kg_hide-input":true,"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"%%writefile competition_score.py\n# Adapted from https://www.kaggle.com/code/metric/ariel-gaussian-log-likelihood\nclass ParticipantVisibleError(Exception):\n    pass\n\ndef competition_score(\n        solution: pd.DataFrame,\n        submission: pd.DataFrame,\n        naive_mean: float,\n        naive_sigma: float,\n        sigma_true: float,\n        row_id_column_name='planet_id',\n    ) -> float:\n    '''\n    This is a Gaussian Log Likelihood based metric. For a submission, which contains the predicted mean (x_hat) and variance (x_hat_std),\n    we calculate the Gaussian Log-likelihood (GLL) value to the provided ground truth (x). We treat each pair of x_hat,\n    x_hat_std as a 1D gaussian, meaning there will be 283 1D gaussian distributions, hence 283 values for each test spectrum,\n    the GLL value for one spectrum is the sum of all of them.\n\n    Inputs:\n        - solution: Ground Truth spectra (from test set)\n            - shape: (nsamples, n_wavelengths)\n        - submission: Predicted spectra and errors (from participants)\n            - shape: (nsamples, n_wavelengths*2)\n        naive_mean: (float) mean from the train set.\n        naive_sigma: (float) standard deviation from the train set.\n        sigma_true: (float) essentially sets the scale of the outputs.\n    '''\n\n    del solution[row_id_column_name]\n    del submission[row_id_column_name]\n\n    if submission.min().min() < 0:\n        raise ParticipantVisibleError('Negative values in the submission')\n    for col in submission.columns:\n        if not pd.api.types.is_numeric_dtype(submission[col]):\n            raise ParticipantVisibleError(f'Submission column {col} must be a number')\n\n    n_wavelengths = len(solution.columns)\n    if len(submission.columns) != n_wavelengths*2:\n        raise ParticipantVisibleError('Wrong number of columns in the submission')\n\n    y_pred = submission.iloc[:, :n_wavelengths].values\n    # Set a non-zero minimum sigma pred to prevent division by zero errors.\n    sigma_pred = np.clip(submission.iloc[:, n_wavelengths:].values, a_min=10**-15, a_max=None)\n    y_true = solution.values\n\n    GLL_pred = np.sum(scipy.stats.norm.logpdf(y_true, loc=y_pred, scale=sigma_pred))\n    GLL_true = np.sum(scipy.stats.norm.logpdf(y_true, loc=y_true, scale=sigma_true * np.ones_like(y_true)))\n    GLL_mean = np.sum(scipy.stats.norm.logpdf(y_true, loc=naive_mean * np.ones_like(y_true), scale=naive_sigma * np.ones_like(y_true)))\n\n    submit_score = (GLL_pred - GLL_mean)/(GLL_true - GLL_mean)\n    return float(np.clip(submit_score, 0.0, 1.0))","metadata":{"_kg_hide-input":true,"papermill":{"duration":0.029243,"end_time":"2024-08-03T12:42:17.162838","exception":false,"start_time":"2024-08-03T12:42:17.133595","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2024-11-30T11:41:32.576376Z","iopub.execute_input":"2024-11-30T11:41:32.576756Z","iopub.status.idle":"2024-11-30T11:41:32.582744Z","shell.execute_reply.started":"2024-11-30T11:41:32.576725Z","shell.execute_reply":"2024-11-30T11:41:32.581962Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"%%writefile postprocessing.py\n\ndef postprocessing(pred_array, index, sigma_pred):\n    \"\"\"Create a submission dataframe from its components\n    \n    Parameters:\n    pred_array: ndarray of shape (n_samples, 283)\n    index: pandas.Index of length n_samples with name 'planet_id'\n    sigma_pred: float\n    \n    Return value:\n    df: DataFrame of shape (n_samples, 566) with planet_id as index\n    \"\"\"\n    return pd.concat([pd.DataFrame(pred_array.clip(0, None), index=index, columns=wavelengths.columns),\n                      pd.DataFrame(sigma_pred, index=index, columns=[f\"sigma_{i}\" for i in range(1, 284)])],\n                     axis=1)\n\n","metadata":{"execution":{"iopub.status.busy":"2024-11-30T11:41:33.386651Z","iopub.execute_input":"2024-11-30T11:41:33.387351Z","iopub.status.idle":"2024-11-30T11:41:33.39276Z","shell.execute_reply.started":"2024-11-30T11:41:33.387322Z","shell.execute_reply":"2024-11-30T11:41:33.391816Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"exec(open('competition_score.py', 'r').read())\nexec(open('postprocessing.py', 'r').read())\n\noof_df = postprocessing(oof_pred, train_adc_info.index, sigma_pred)\ndisplay(oof_df)\n\ngll_score = competition_score(train_labels.copy().reset_index(),\n                              oof_df.copy().reset_index(),\n                              naive_mean=train_labels.values.mean(),\n                              naive_sigma=train_labels.values.std(),\n                              sigma_true=0.000003)\nprint(f\"# Estimated competition score: {gll_score:.3f}\")\n# Estimated competition score: 0.259","metadata":{"papermill":{"duration":0.231627,"end_time":"2024-08-03T13:00:10.149224","exception":false,"start_time":"2024-08-03T13:00:09.917597","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2024-11-30T11:41:34.248163Z","iopub.execute_input":"2024-11-30T11:41:34.248828Z","iopub.status.idle":"2024-11-30T11:41:34.743638Z","shell.execute_reply.started":"2024-11-30T11:41:34.248801Z","shell.execute_reply":"2024-11-30T11:41:34.742713Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Refitting and saving the model","metadata":{}},{"cell_type":"code","source":"# Refit the model to the full dataset\nmodel.fit(train, train_labels)\nwith open('model.pickle', 'wb') as f:\n    pickle.dump(model, f)\nwith open('sigma_pred.pickle', 'wb') as f:\n    pickle.dump(sigma_pred, f)\n","metadata":{"execution":{"iopub.status.busy":"2024-11-30T11:41:37.798083Z","iopub.execute_input":"2024-11-30T11:41:37.798506Z","iopub.status.idle":"2024-11-30T11:41:37.808956Z","shell.execute_reply.started":"2024-11-30T11:41:37.798474Z","shell.execute_reply":"2024-11-30T11:41:37.808026Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Submission","metadata":{"papermill":{"duration":0.077479,"end_time":"2024-08-03T13:00:10.305108","exception":false,"start_time":"2024-08-03T13:00:10.227629","status":"completed"},"tags":[]}},{"cell_type":"code","source":"# Load the data\ntest_adc_info = pd.read_csv('/kaggle/input/ariel-data-challenge-2024/test_adc_info.csv',\n                           index_col='planet_id')\nsample_submission = pd.read_csv('/kaggle/input/ariel-data-challenge-2024/sample_submission.csv',\n                                index_col='planet_id')\nf_raw_test = f_read_and_preprocess('test', test_adc_info, sample_submission.index)\na_raw_test = a_read_and_preprocess('test', test_adc_info, sample_submission.index)\ntest = feature_engineering(f_raw_test, a_raw_test)\n\n# Load the model\nwith open('model.pickle', 'rb') as f:\n    model = pickle.load(f)\nwith open('sigma_pred.pickle', 'rb') as f:\n    sigma_pred = pickle.load(f)\n\n# Predict\ntest_pred = model.predict(test)\n\n# Package into submission file\nsub_df = postprocessing(test_pred, sample_submission.index, sigma_pred)\ndisplay(sub_df)\nsub_df.to_csv('submission.csv')\n#!head submission.csv","metadata":{"papermill":{"duration":1.822768,"end_time":"2024-08-03T13:00:12.204995","exception":false,"start_time":"2024-08-03T13:00:10.382227","status":"completed"},"tags":[],"execution":{"iopub.status.busy":"2024-11-30T11:41:39.895352Z","iopub.execute_input":"2024-11-30T11:41:39.895688Z","iopub.status.idle":"2024-11-30T11:41:43.078731Z","shell.execute_reply.started":"2024-11-30T11:41:39.895651Z","shell.execute_reply":"2024-11-30T11:41:43.07744Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null}]}