{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.11.13","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"none","dataSources":[{"sourceId":101849,"databundleVersionId":13093295,"isSourceIdPinned":false,"sourceType":"competition"},{"sourceId":13302394,"sourceType":"datasetVersion","datasetId":8076891}],"dockerImageVersionId":31153,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"This notebook demonstrates the training of gradient boostings for `dmu_airs`, `dsigma_fgs`, `dmu_fgs`, `dsigma_fgs`. The attached dataset contains saved features from fitted models (the exact features used for boosting inference in main `ARIEL-2025-inference`).","metadata":{}},{"cell_type":"code","source":"%%writefile score.py\n\nimport numpy as np\nimport pandas as pd\nimport pandas.api.types\nimport scipy.stats\n\n\nclass ParticipantVisibleError(Exception):\n    pass\n\n\ndef competition_score(\n    solution: pd.DataFrame,\n    submission: pd.DataFrame,\n    row_id_column_name: str,\n    naive_mean: float,\n    naive_sigma: float,\n    fsg_sigma_true: float = 1e-6,\n    airs_sigma_true: float = 1e-5,\n    fgs_weight: float = 1,\n):\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        fsg_sigma_true: (float) standard deviation from the FSG1 instrument for the test set.\n        airs_sigma_true: (float) standard deviation from the AIRS instrument for the test set.\n        fgs_weight: (float) relative weight of the fgs channel\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 pandas.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    sigma_true = np.append(\n        np.array(\n            [\n                fsg_sigma_true,\n            ]\n        ),\n        np.ones(n_wavelengths - 1) * airs_sigma_true,\n    )\n    y_true = solution.values\n\n    GLL_pred = scipy.stats.norm.logpdf(y_true, loc=y_pred, scale=sigma_pred)\n    GLL_true = scipy.stats.norm.logpdf(y_true, loc=y_true, scale=sigma_true * np.ones_like(y_true))\n    GLL_mean = scipy.stats.norm.logpdf(y_true, loc=naive_mean * np.ones_like(y_true), scale=naive_sigma * np.ones_like(y_true))\n\n    # normalise the score, right now it becomes a matrix instead of a scalar.\n    ind_scores = (GLL_pred - GLL_mean) / (GLL_true - GLL_mean)\n\n    # print(GLL_mean.shape)\n\n    # print(GLL_pred, GLL_true)\n\n    # ind_scores[:, 0] = 0.5\n\n    weights = np.append(np.array([fgs_weight]), np.ones(len(solution.columns) - 1))\n    weights = weights * np.ones_like(ind_scores)\n    submit_score = np.average(ind_scores, weights=weights)\n    return float(np.clip(submit_score, 0.0, 1.0)), ind_scores, np.average(ind_scores, weights=np.append(np.array([fgs_weight]), np.ones(len(solution.columns) - 1)), axis=1)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-10-08T14:38:14.044065Z","iopub.execute_input":"2025-10-08T14:38:14.044413Z","iopub.status.idle":"2025-10-08T14:38:14.052057Z","shell.execute_reply.started":"2025-10-08T14:38:14.044365Z","shell.execute_reply":"2025-10-08T14:38:14.051222Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import pandas as pd\nimport numpy as np\nimport matplotlib.pyplot as plt\nfrom catboost import CatBoostRegressor\nfrom sklearn.model_selection import train_test_split\nfrom sklearn.model_selection import KFold\nimport xgboost as xgb\nimport lightgbm as lgb\nimport tqdm\n\n\nDATA_ROOT = \"/kaggle/input/ariel-data-challenge-2025\"\nnpy_root = \"/kaggle/input/calibrated/features/\"\n\nmus = np.load(f\"{npy_root}/mus.npy\")\nsigmas = np.load(f\"{npy_root}/sigmas.npy\")\nfeatures = np.load(f\"{npy_root}/features_fgs.npy\")\ncleans = np.load(f\"{npy_root}/cleans_fgs.npy\") & np.load(f\"{npy_root}/cleans_airs.npy\")\n\ngt_fgs = pd.read_csv(f'{DATA_ROOT}/train.csv', index_col='planet_id').values[:, :1]\ngt_airs = pd.read_csv(f'{DATA_ROOT}/train.csv', index_col='planet_id').values[:, 1:]\n\ngt_dmu = (gt_fgs - mus[:, :1])\ngt_dsigma = np.abs(gt_fgs[:, :1] - mus[:, :1]) - sigmas[:, :1]\nfeatures = features[:, 1:]\n\ncleans[289] = False\ncleans[496] = False\nsigmas[289] = 0.01\nsigmas[496] = 0.01\n\nnew_mus = np.copy(mus)\nnew_sigmas = np.copy(sigmas)\n\nn_splits = 5\nkfold = KFold(n_splits=n_splits, shuffle=False)\n\n# for i, (train_index, test_index) in enumerate(kfold.split(features)):\nfor i, (train_index, test_index) in tqdm.tqdm(enumerate(kfold.split(features)), total=n_splits):\n    test = features[test_index]\n    clean = cleans[test_index]\n    X_train, y_train_dmu, y_train_dsigma = features[train_index], gt_dmu[train_index], gt_dsigma[train_index]\n    X_test, y_test_dmu, y_test_dsigma = features[test_index], gt_dmu[test_index], gt_dsigma[test_index]\n\n    X_train = X_train[sigmas[train_index].mean(axis=1) < 0.0003]\n    y_train_dmu = y_train_dmu[sigmas[train_index].mean(axis=1) < 0.0003]\n    y_train_dsigma = y_train_dsigma[sigmas[train_index].mean(axis=1) < 0.0003]\n\n    model_dmu = CatBoostRegressor(\n        # objective=\"MultiRMSE\",\n        iterations=1000,\n        learning_rate=0.015,\n        depth=4,\n        verbose=False,\n        random_state=42,\n        # task_type=\"GPU\",\n    )\n    model_dmu.fit(X_train, y_train_dmu)\n    model_dmu.save_model(f\"boost_dmu_fgs_v4.cbm\")\n    y_pred_dmu = model_dmu.predict(X_test)[:, None]\n    # y_pred_dmu = np.repeat(y_pred_dmu, 47, axis=1)\n\n    model_sigma = CatBoostRegressor(\n        iterations=1000,\n        learning_rate=0.02,\n        depth=4,\n        verbose=False,\n        random_state=42,\n        # task_type=\"GPU\",\n    )\n    model_sigma.fit(X_train, y_train_dsigma)\n    model_sigma.save_model(f\"boost_dsigma_fgs_v4.cbm\")\n    y_pred_dsigma = model_sigma.predict(X_test)[:, None]\n\n    mu_old = mus[test_index].copy()\n    sigma_old = sigmas[test_index].copy()\n\n    mask_dmu = ((np.abs(y_pred_dmu) / mu_old[:, :1] < 0.03) & clean[:, None]).squeeze()\n    mask_sigma = ((0.7 < (y_pred_dsigma + sigma_old[:, :1]) / sigma_old[:, :1]) & clean[:, None]).squeeze()\n\n    mu_old[mask_dmu, :1] = mu_old[mask_dmu, :1] + y_pred_dmu[mask_dmu, :]\n    sigma_old[mask_sigma, :1] = sigma_old[mask_sigma, :1] + y_pred_dsigma[mask_sigma, :]\n    sigma_old[mask_dmu, :1] = sigma_old[mask_dmu, :1] * 0.9\n\n    new_mus[test_index] = mu_old\n    new_sigmas[test_index] = sigma_old\n\n#######################################################################3\n\ndef format_mu(planet_codes, mu, wavelengths):\n    return pd.DataFrame(mu.clip(0, None), index=planet_codes, columns=wavelengths.columns)\n\ndef format_sigma(planet_codes, sigma):\n    return pd.DataFrame(sigma, index=planet_codes, columns=[f\"sigma_{i}\" for i in range(1, 284)])\n\ndef format_to_submission(planet_codes, mus, sigmas):\n    mus = pd.concat(mus)\n    sigmas = pd.concat(sigmas)\n    df_index = pd.DataFrame({'planet_id': list(planet_codes)}, index=planet_codes)\n    df_submission = df_index.join(mus).join(sigmas)\n    return df_submission.reset_index(drop=True)\n\nfrom score import competition_score\nwavelengths = pd.read_csv(f'{DATA_ROOT}/wavelengths.csv')\ntrain_labels = pd.read_csv(f'{DATA_ROOT}/train.csv', index_col='planet_id')\nstar_info = pd.read_csv(f\"{DATA_ROOT}/train_star_info.csv\")\nplanet_ids = star_info[\"planet_id\"].tolist()\n\nmu_final = format_mu(planet_ids, new_mus, wavelengths)\nsigma_final = format_sigma(planet_ids, new_sigmas)\n\nsubmission_df = format_to_submission(train_labels.index, mus=[mu_final], sigmas=[sigma_final])\n\nprint(f\"MSE: {np.square(train_labels.values[:, 1:] - mu_final.loc[train_labels.index].values[:, 1:]).mean()}\")\n\ngll_score_avg, sep_scores, ind_avg = competition_score(\n    train_labels.copy().reset_index(),\n    submission_df.copy(),\n    row_id_column_name=\"planet_id\",\n    naive_mean=train_labels.values.mean(),\n    naive_sigma=train_labels.values.std(),\n    fsg_sigma_true=1e-6,\n    airs_sigma_true=1e-5,\n    fgs_weight=57.846,\n)\nprint('SCORE:', gll_score_avg)","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true,"execution":{"iopub.status.busy":"2025-10-08T14:39:18.237692Z","iopub.execute_input":"2025-10-08T14:39:18.238249Z","iopub.status.idle":"2025-10-08T14:39:32.156003Z","shell.execute_reply.started":"2025-10-08T14:39:18.238222Z","shell.execute_reply":"2025-10-08T14:39:32.15505Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import pandas as pd\nimport numpy as np\nimport matplotlib.pyplot as plt\nfrom catboost import CatBoostRegressor\nfrom sklearn.model_selection import train_test_split\nfrom sklearn.model_selection import KFold\n# import xgboost as xgb\n# import lightgbm as lgb\nimport tqdm\n\n\nDATA_ROOT = \"/kaggle/input/ariel-data-challenge-2025\"\nnpy_root = \"/kaggle/input/calibrated/features/\"\n\nmus = np.load(f\"{npy_root}/mus.npy\")\nsigmas = np.load(f\"{npy_root}/sigmas.npy\")\nfeatures = np.load(f\"{npy_root}/features_airs.npy\")\ncleans = np.load(f\"{npy_root}/cleans_airs.npy\")\n\ngt_fgs = pd.read_csv(f'{DATA_ROOT}/train.csv', index_col='planet_id').values[:, 0]\ngt_airs = pd.read_csv(f'{DATA_ROOT}/train.csv', index_col='planet_id').values[:, 1:]\ngt_dmu = gt_airs[:, 1:].mean(axis=1) - mus[:, 1:].mean(axis=1)\ngt_sigma = features[:, 0] - sigmas.mean(axis=1)\nfeatures = features[:, 1:]\n\ncleans[289] = False\ncleans[496] = False\nsigmas[289] = 0.01\nsigmas[496] = 0.01\n\nnew_mus = np.copy(mus)\nnew_sigmas = np.copy(sigmas)\n\nn_splits = 5\nkfold = KFold(n_splits=n_splits, shuffle=False)\n\n# for i, (train_index, test_index) in enumerate(kfold.split(features)):\nfor i, (train_index, test_index) in tqdm.tqdm(enumerate(kfold.split(features)), total=n_splits):\n    test = features[test_index]\n    X_test, y_test_dmu, y_test_sigma = features[test_index], gt_dmu[test_index], gt_sigma[test_index]\n\n    clean = cleans[test_index]\n\n    X_train, y_train_dmu, y_train_sigma = features[train_index], gt_dmu[train_index], gt_sigma[train_index]\n    X_train = X_train[sigmas[train_index].mean(axis=1) < 0.0003]\n    y_train_dmu = y_train_dmu[sigmas[train_index].mean(axis=1) < 0.0003]\n    y_train_sigma = y_train_sigma[sigmas[train_index].mean(axis=1) < 0.0003]\n\n    model_dmu = CatBoostRegressor(\n        iterations=1000,\n        # learning_rate=0.015,\n        learning_rate=0.015,\n        depth=4,\n        verbose=False,\n        random_state=42,\n    )\n    model_dmu.fit(X_train, y_train_dmu)\n    model_dmu.save_model(f\"boost_dmu_airs_v4.cbm\")\n    y_pred_dmu = model_dmu.predict(X_test)\n\n    X_train = X_train[:, :-250]\n\n    model_sigma = CatBoostRegressor(\n        # objective=\"MultiRMSE\",\n        iterations=1500,\n        learning_rate=0.015,\n        depth=2,\n        # depth=4,\n        verbose=False,\n        random_state=42,\n        # task_type=\"GPU\",\n    )\n    model_sigma.fit(X_train, y_train_sigma)\n    model_sigma.save_model(f\"boost_dsigma_airs_v4.cbm\")\n    y_pred_sigma = model_sigma.predict(X_test)\n\n    mu_old = mus[test_index].copy()\n    sigma_old = sigmas[test_index].copy()\n\n    mask_dmu = (np.abs(y_pred_dmu) / mu_old[:, 1:].mean(axis=1) < 0.02) & clean\n    mask_sigma = (0.3 < (y_pred_sigma + sigma_old.mean(axis=1)) / sigma_old.mean(axis=1)) & clean\n\n    mu_old[mask_dmu, 1:] = mu_old[mask_dmu, 1:] + y_pred_dmu[mask_dmu, None]\n    # mu_old[mask_dmu, :] = mu_old[mask_dmu, :] + y_pred_dmu[mask_dmu, None]\n    # sigma_old[mask_sigma, 1:] = y_pred_sigma[mask_sigma, None]\n    sigma_old[mask_sigma, 1:] = sigma_old[mask_sigma, 1:] + y_pred_sigma[mask_sigma, None]\n    sigma_old[mask_dmu, 1:] = sigma_old[mask_dmu, 1:] * 0.9\n\n    new_mus[test_index] = mu_old\n    new_sigmas[test_index] = sigma_old\n\n#######################################################################3\n\ndef format_mu(planet_codes, mu, wavelengths):\n    return pd.DataFrame(mu.clip(0, None), index=planet_codes, columns=wavelengths.columns)\n\ndef format_sigma(planet_codes, sigma):\n    return pd.DataFrame(sigma, index=planet_codes, columns=[f\"sigma_{i}\" for i in range(1, 284)])\n\ndef format_to_submission(planet_codes, mus, sigmas):\n    mus = pd.concat(mus)\n    sigmas = pd.concat(sigmas)\n    df_index = pd.DataFrame({'planet_id': list(planet_codes)}, index=planet_codes)\n    df_submission = df_index.join(mus).join(sigmas)\n    return df_submission.reset_index(drop=True)\n\nerrors = np.abs(gt_airs - mus[:, 1:])\n\nfrom score import competition_score\nwavelengths = pd.read_csv(f'{DATA_ROOT}/wavelengths.csv')\ntrain_labels = pd.read_csv(f'{DATA_ROOT}/train.csv', index_col='planet_id')\nstar_info = pd.read_csv(f\"{DATA_ROOT}/train_star_info.csv\")\nplanet_ids = star_info[\"planet_id\"].tolist()\n\nmu_final = format_mu(planet_ids, new_mus, wavelengths)\nsigma_final = format_sigma(planet_ids, new_sigmas)\n\nsubmission_df = format_to_submission(train_labels.index, mus=[mu_final], sigmas=[sigma_final])\n\nprint(f\"MSE: {np.square(train_labels.values[:, 1:] - mu_final.loc[train_labels.index].values[:, 1:]).mean()}\")\n\ngll_score_avg, sep_scores, ind_avg = competition_score(\n    train_labels.copy().reset_index(),\n    submission_df.copy(),\n    row_id_column_name=\"planet_id\",\n    naive_mean=train_labels.values.mean(),\n    naive_sigma=train_labels.values.std(),\n    fsg_sigma_true=1e-6,\n    airs_sigma_true=1e-5,\n    fgs_weight=57.846,\n)\nprint('SCORE:', gll_score_avg)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-10-08T14:40:48.080259Z","iopub.execute_input":"2025-10-08T14:40:48.08138Z","iopub.status.idle":"2025-10-08T14:42:07.283629Z","shell.execute_reply.started":"2025-10-08T14:40:48.081325Z","shell.execute_reply":"2025-10-08T14:42:07.282734Z"}},"outputs":[],"execution_count":null}]}