{"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":"none","dataSources":[{"sourceId":70367,"databundleVersionId":9188054,"sourceType":"competition"},{"sourceId":9747089,"sourceType":"datasetVersion","datasetId":5967131}],"dockerImageVersionId":30775,"isInternetEnabled":false,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"# ELEVENTH!\n# https://www.youtube.com/watch?v=HbDnxzrbxn4\n\nimport pandas as pd\nimport numpy as np\nimport itertools\nimport functools\nimport joblib\nimport gc\n\nimport pywt\nimport math\n\nfrom functools import reduce\n\nfrom sklearn.linear_model import LinearRegression\n\nfrom scipy.ndimage import median_filter\nfrom scipy.optimize import minimize\nfrom scipy.ndimage import gaussian_filter1d\nfrom scipy.signal import find_peaks\n\nimport matplotlib.pyplot as plt\n\nfrom tqdm import tqdm\n\nmode = 'test'\n# mode = 'train'","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def sigma_clip_custom(data, sigma=3, maxiters=5):\n    mask = np.zeros_like(data, dtype=bool)\n    \n    for _ in range(maxiters):\n        if np.all(~mask):\n            mean = data.mean()\n            std = data.std()\n        else:\n            mean = data[~mask].mean()\n            std = data[~mask].std()\n        \n        new_mask = (data > mean + sigma * std) | (data < mean - sigma * std)\n        updated_mask = mask | new_mask\n        \n        if np.array_equal(updated_mask, mask):\n            break\n        \n        mask = updated_mask\n    \n    return mask\n\ndef apply_linear_corr(linear_corr, clean_signal):\n    linear_corr = np.flip(linear_corr, axis=0)\n    for x, y in itertools.product(\n                range(clean_signal.shape[1]), range(clean_signal.shape[2])\n            ):\n        poli = np.poly1d(linear_corr[:, x, y])\n        clean_signal[:, x, y] = poli(clean_signal[:, x, y])\n    return clean_signal\n\ndef clean_dark(signal, dark, dt):\n    dark = np.tile(dark, (signal.shape[0], 1, 1))\n    signal -= dark* dt[:, np.newaxis, np.newaxis]\n    return signal\n\ndef preproc_airs(planet_id, adc_info, dataset='train'):\n    cut_inf, cut_sup = 39, 321\n    original_sensor_size = [11250, 32, 356]\n#     target_sensor_size = [1, 32, cut_sup-cut_inf]\n    target_sensor_size = [1, 32, 356]\n    sensor = 'AIRS-CH0'\n    linear_corr = (6, 32, 356)\n\n    \n    signal = pd.read_parquet(f'/kaggle/input/ariel-data-challenge-2024/{dataset}/{planet_id}/{sensor}_signal.parquet').to_numpy()\n    dark_frame = pd.read_parquet(f'/kaggle/input/ariel-data-challenge-2024/{dataset}/' + str(planet_id) + '/' + sensor + '_calibration/dark.parquet', engine='pyarrow').to_numpy()\n    dead_frame = pd.read_parquet(f'/kaggle/input/ariel-data-challenge-2024/{dataset}/' + str(planet_id) + '/' + sensor + '_calibration/dead.parquet', engine='pyarrow').to_numpy()\n    flat_frame = pd.read_parquet(f'/kaggle/input/ariel-data-challenge-2024/{dataset}/' + str(planet_id) + '/' + sensor + '_calibration/flat.parquet', engine='pyarrow').to_numpy()\n    linear_corr = pd.read_parquet(f'/kaggle/input/ariel-data-challenge-2024/{dataset}/' + str(planet_id) + '/' + sensor + '_calibration/linear_corr.parquet').values.astype(np.float64).reshape(linear_corr)\n\n    signal = signal.reshape(original_sensor_size) \n    \n    gain = adc_info.loc[planet_id][f'{sensor}_adc_gain']\n    offset = adc_info.loc[planet_id][f'{sensor}_adc_offset']\n    \n    signal = signal / gain + offset\n\n    hot = sigma_clip_custom(dark_frame, sigma=5, maxiters=5)\n#     hot = sigma_clip(dark_frame, sigma=5, maxiters=5).mask\n\n    dt = np.ones(len(signal))*0.1\n    dt[1::2] += 4.5 \n\n    signal = signal.clip(0)\n    linear_corr_signal = apply_linear_corr(linear_corr, signal)\n    signal = clean_dark(linear_corr_signal, dark_frame, dt)\n\n    flat = flat_frame.reshape(target_sensor_size)\n    flat[dead_frame.reshape(target_sensor_size)] = np.nan\n    flat[hot.reshape(target_sensor_size)] = np.nan\n    signal = signal / flat\n    \n    signal = signal[:,10:22,:]\n\n    mean_signal = np.nanmean(signal, axis=1) \n    cds_signal = (mean_signal[1::2] - mean_signal[0::2])\n\n    return cds_signal","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"wavelenghts_file =  f'/kaggle/input/ariel-data-challenge-2024/wavelengths.csv'\nadc_info_file =  f'/kaggle/input/ariel-data-challenge-2024/{mode}_adc_info.csv'\n\n# y_file =  f'/kaggle/input/ariel-data-challenge-2024/train_labels.csv'\n# y_df = pd.read_csv(y_file, index_col='planet_id')\nwavelenghts_df =  f'/kaggle/input/ariel-data-challenge-2024/wavelengths.csv'\n\nadc_info = pd.read_csv(adc_info_file, index_col='planet_id')\nplanet_ids = adc_info.index.to_list()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Utilities functions\ndef linear_transform(signal):\n    x = np.arange(len(signal)).reshape(-1, 1)\n    reg = LinearRegression().fit(x, signal)\n    trend_line = reg.predict(x).flatten()\n    signal_transformed = signal - trend_line    \n    return signal_transformed, trend_line\n\ndef extend_signal_with_polyfit(x, y, extend_points, deg):\n    p = np.polyfit(x, y, deg)\n    x_left = np.linspace(x[0] - extend_points, x[0], extend_points)\n    y_left = np.polyval(p, x_left)\n    x_right = np.linspace(x[-1], x[-1] + extend_points, extend_points)\n    y_right = np.polyval(p, x_right)\n    x_extended = np.concatenate([x_left, x, x_right])\n    y_extended = np.concatenate([y_left, y, y_right])\n\n    return x_extended, y_extended\n\ndef polynomial_transform2(signal, degree, K, L):\n    x_start = np.arange(K)\n    x_end = np.arange(L, len(signal))\n\n    y_start = signal[:K]\n    y_end = signal[L:]\n\n    x_concat = np.concatenate([x_start, x_end])\n    y_concat = np.concatenate([y_start, y_end])\n\n    extend_points = 300\n    x = np.arange(len(signal))\n    x, signal_extended = extend_signal_with_polyfit(x_concat, y_concat, extend_points, deg=degree-1)\n\n    coefs = np.polyfit(x, signal_extended, degree)\n    \n    x_full = np.arange(len(signal))\n    trend_poly = np.polyval(coefs, x_full)\n\n    signal_transformed = signal - trend_poly\n\n    return signal_transformed, trend_poly\n\n\ndef refit_with_polyfit(x, y, deg):\n    nn = len(x)\n    w = np.ones(nn)\n    # w[ : nn // 2] = np.linspace(1, 0.9, nn // 2)\n    # w[nn // 2 : ] = np.linspace(0.9, 1, len(w[nn // 2:]))\n    p = np.polyfit(x, y, deg, w=w)\n\n    return np.polyval(p, x)\n\ndef wavelet_denoising(signal, wavelet='db29', strength=1.0, mode='symmetric', level=None):\n    coeff = pywt.wavedec(signal, wavelet, mode=mode, level=level)\n    coeff[1:] = [pywt.threshold(c, strength * np.std(c)) for c in coeff[1:]]\n    return pywt.waverec(coeff, wavelet, mode=mode)\n\ndef polyfit_denoising(signal, extend_points=100, max_deg=4):\n    x = np.arange(len(signal))\n\n    if len(signal) < 2 * extend_points:\n        return signal\n    x_extended, signal_extended = extend_signal_with_polyfit(x, signal, extend_points, deg=1)\n\n    for deg in range(2, max_deg+1):\n        signal_refit = refit_with_polyfit(x_extended, signal_extended, deg)\n        signal_extended[:extend_points] = signal_refit[:extend_points]\n        signal_extended[-extend_points:] = signal_refit[-extend_points:]\n\n    signal_refit = refit_with_polyfit(x_extended, signal_extended, max_deg)\n\n    return signal_refit[extend_points:-extend_points]\n\ndef transit_formula(s):\n    return 1/s - 1","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Transformation tools\n\ndef dummy(**kwargs):\n    return kwargs\n\n\ndef bin_wavelengths(**kwargs):\n    fromwl = kwargs['from_wavelength']\n    towl = kwargs['to_wavelength']\n    signal = kwargs['signal']\n    kwargs['signal'] = signal[:, fromwl:towl].mean(axis=1)\n    return kwargs\n\ndef normalize_signal(**kwargs):\n    # signal = signal.copy()\n    # wavelength_means = signal[:, :].mean(axis=0)\n    # signal[:, :] /= wavelength_means\n    signal = kwargs['signal']\n    normalized_signal = signal.copy()\n    wavelength_means = signal.mean(axis=0)\n    normalized_signal /= wavelength_means\n    kwargs['signal'] = normalized_signal\n    return kwargs\n\ndef shrink_signal(ratio=12, **kwargs):\n    signal = kwargs['signal']\n    n = signal.shape[0]\n    signal_windowed = signal[:n - (n % ratio)].reshape(-1, ratio)\n    signal_means = signal_windowed.mean(axis=1)\n    signal_stds = signal_windowed.std(axis=1)\n    kwargs['signal'] = signal_means\n    kwargs['signal_stds'] = signal_stds\n\n    kwargs['margin'] //= ratio\n    if 'ingress_left' in kwargs:\n        kwargs['ingress_left'] //= ratio\n        kwargs['ingress_right'] //= ratio\n        kwargs['egress_left'] //= ratio\n        kwargs['egress_right'] //= ratio\n\n    return kwargs\n\ndef linear_detrend(**kwargs):\n    signal = kwargs['signal']\n    signal_transformed, trend_line = linear_transform(signal)\n\n    kwargs['signal'] = signal_transformed + trend_line.mean()\n    kwargs['trend_line'] = trend_line\n    kwargs['signal_before_detrend'] = signal\n\n    return kwargs\n\ndef polynomial_detrend(**kwargs):\n    signal = kwargs['signal']\n    degree = kwargs.get('degree', 5)\n    K = kwargs['ingress_left'] - 10\n    L = kwargs['egress_right'] + 10\n\n    signal_transformed, trend_poly = polynomial_transform2(signal, degree, K, L)\n\n    kwargs['signal'] = signal_transformed\n    kwargs['trend_poly'] = trend_poly\n\n    return kwargs\n\ndef get_transit_phases2(**kwargs):\n    signal = kwargs['signal'] + 1\n\n    grad1 = np.gradient(signal)\n    mn = np.argmin(grad1)\n    mx = np.argmax(grad1)\n    \n#     peaks = find_peaks(grad1, width=30)[0]\n#     ingress_left = peaks[peaks < mn][-1]\n#     ingress_right = peaks[peaks > mn][0]\n    \n#     dips = find_peaks(-grad1, width=20)[0]\n#     egress_left = dips[dips < mx][-1]\n#     egress_right = dips[dips > mx][0]\n\n    peaks = find_peaks(grad1, width=25)[0]\n    ingress_left = peaks[(peaks < mn) & (mn - peaks <= 250)][-1] if np.any((peaks < mn) & (mn - peaks <= 250)) else mn-110\n    ingress_right = peaks[(peaks > mn) & (peaks - mn <= 250)][0] if np.any((peaks > mn) & (peaks - mn <= 250)) else mn+110\n    \n    dips = find_peaks(-grad1, width=25)[0]\n    egress_left = dips[(dips < mx) & (mx - dips <= 250)][-1] if np.any((dips < mx) & (mx - dips <= 250)) else mx-110\n    egress_right = dips[(dips > mx) & (dips - mx <= 250)][0] if np.any((dips > mx) & (dips - mx <= 250)) else mx+110\n\n    kwargs['ingress_left'] = ingress_left\n    kwargs['ingress_right'] = ingress_right\n    kwargs['egress_left'] = egress_left\n    kwargs['egress_right'] = egress_right\n\n    return kwargs\n\n\ndef get_transit_phases(edge_margin=300, edge_margin2=250, sigma1=100, sigma2=10, **kwargs):\n    signal  = kwargs['signal']\n    singal_transformed, trend_line = linear_transform(signal)\n    signal_denoised = gaussian_filter1d(singal_transformed, sigma=sigma1)\n    grad = np.gradient(signal_denoised)\n    grad_transformed, _ = linear_transform(grad)\n    ingress_mean = np.argmin(grad_transformed[edge_margin:-edge_margin]) + edge_margin\n    egress_mean = np.argmax(grad_transformed[edge_margin:-edge_margin]) + edge_margin\n\n    grad2 = np.gradient(grad_transformed)\n    grad2_denoised = gaussian_filter1d(grad2, sigma=sigma2)\n    grad2_transformed, grad_trend_line = linear_transform(grad2_denoised)\n\n    ingress_left = np.argmin(grad2_transformed[ingress_mean-edge_margin2:ingress_mean]) + ingress_mean-edge_margin2\n    ingress_right = np.argmax(grad2_transformed[ingress_mean:ingress_mean+edge_margin2]) + ingress_mean\n\n    egress_left = np.argmax(grad2_transformed[egress_mean-edge_margin2:egress_mean]) + egress_mean-edge_margin2\n    egress_right = np.argmin(grad2_transformed[egress_mean:egress_mean+edge_margin2]) + egress_mean \n\n    kwargs['ingress_left'] = ingress_left\n    kwargs['ingress_right'] = ingress_right\n    kwargs['egress_left'] = egress_left\n    kwargs['egress_right'] = egress_right\n\n    return kwargs\n\n\ndef split_by_phases(**kwargs):\n    margin = kwargs['margin']\n    signal = kwargs['signal']\n    ingress_left = kwargs['ingress_left']\n    ingress_right = kwargs['ingress_right']\n    egress_left = kwargs['egress_left']\n    egress_right = kwargs['egress_right']\n\n    signal_before = signal[:ingress_left-margin]\n    signal_transit = signal[ingress_right+margin:egress_left-margin]\n    signal_after = signal[egress_right+margin:]\n\n    signals = [signal_before, signal_transit, signal_after]\n    kwargs['signals'] = signals\n\n    return kwargs\n\n\n# Wavelet denoising\ndef denoise_signal_with_wavelet(wavelet='db29', strength=1.0, mode='symmetric', level=None, **kwargs):\n    signal = kwargs['signal']\n    denoised_signal = wavelet_denoising(signal, wavelet=wavelet, strength=strength, mode=mode, level=level)\n    kwargs['signal'] = denoised_signal\n    return kwargs\n\n\ndef filter_signal_with_gaussian(sigma=10, **kwargs):\n    signal = kwargs['signal']\n    filtered_signal = gaussian_filter1d(signal, sigma=sigma)\n    kwargs['signal'] = filtered_signal\n    kwargs['signal_before_gaussian'] = signal\n    return kwargs\n\ndef filter_signal_with_polyfit(extend_points=100, max_deg=4, **kwargs):\n    signal = kwargs['signal']\n    filtered_signal = polyfit_denoising(signal, extend_points=extend_points, max_deg=max_deg)\n    kwargs['signal'] = filtered_signal\n    return kwargs\n\ndef filter_signal_with_median(window=5, **kwargs):\n    signal = kwargs['signal']\n    filtered_signal = median_filter(signal, size=window)\n    kwargs['signal'] = filtered_signal\n    return kwargs\n\ndef denoise_signals_with_wavelet(wavelet='db29', strength=1.0, mode='symmetric', level=None, **kwargs):\n    signals = kwargs['signals']\n    denoised_signals = []\n    for signal in signals:\n        denoised_signal = wavelet_denoising(signal, wavelet=wavelet, strength=strength, mode=mode, level=level)\n        denoised_signals.append(denoised_signal)\n    kwargs['signals'] = denoised_signals\n    return kwargs\n\ndef filter_signals_with_gaussian(sigma=10, **kwargs):\n    signals = kwargs['signals']\n    filtered_signals = []\n    for signal in signals:\n        filtered_signal = gaussian_filter1d(signal, sigma=sigma)\n        filtered_signals.append(filtered_signal)\n    kwargs['signals'] = filtered_signals\n    return kwargs\n\ndef filter_signals_with_polyfit(extend_points=100, max_deg=4, **kwargs):\n    signals = kwargs['signals']\n    filtered_signals = []\n    for signal in signals:\n        filtered_signal = polyfit_denoising(signal, extend_points=extend_points, max_deg=max_deg)\n        filtered_signals.append(filtered_signal)\n    kwargs['signals'] = filtered_signals\n    return kwargs\n\ndef reconstruct_signals(**kwargs):\n    margin = kwargs['margin']\n    signal_before, signal_transit, signal_after = kwargs['signals']\n    ingress_left = kwargs['ingress_left']\n    ingress_right = kwargs['ingress_right']\n    egress_left = kwargs['egress_left']\n    egress_right = kwargs['egress_right']\n\n    signal_reconstructed = np.concatenate([\n        signal_before, \n        np.ones(ingress_right-ingress_left + 2 * margin), \n        signal_transit, \n        np.ones(egress_right-egress_left + 2 * margin), \n        signal_after]\n    )\n\n    mask = np.ones_like(signal_reconstructed, dtype=bool)\n    mask[ingress_left-margin:ingress_right+margin] = False\n    mask[egress_left-margin:egress_right+margin] = False\n    x = np.arange(signal_reconstructed.shape[0])\n    x_reconstructed = x[mask]\n    y_reconstructed = signal_reconstructed[mask]\n\n    kwargs['x'] = x_reconstructed\n    kwargs['signal'] = y_reconstructed\n    return kwargs\n\n\ndef objective(adjustment, x, signal, ingress, egress, deg=4):\n    s = adjustment[0]\n    adjusted_signal = signal.copy()\n\n    adjusted_signal[ingress:egress] *= (1 + s)\n\n    coefficients = np.polyfit(x, adjusted_signal, deg=deg)\n    polyfit_values = np.polyval(coefficients, x)\n\n    mse = np.mean((adjusted_signal - polyfit_values) ** 2) * 1e3\n    return mse\n\ndef fit_transit(deg=4, **kwargs):\n    initial_guess = [0.0001]\n    bounds = [(0, 1)]\n\n    x = kwargs['x']\n    signal = 100 * kwargs['signal']\n    signals = kwargs['signals']\n    ingress = len(signals[0])\n    egress = len(signals[0]) + len(signals[1])\n    kwargs['ingress'] = ingress\n    kwargs['egress'] = egress\n    result = minimize(lambda adjustment: objective(adjustment, x, signal, ingress, egress, deg=deg), \n                          initial_guess, bounds=bounds, method='Nelder-Mead')\n    adjustment = result.x[-1]\n    kwargs['adjustment'] = adjustment\n    return kwargs\n\n\ndef get_features(**kwargs):\n    features = []\n    for feature in kwargs['features']:\n        feature = kwargs[feature]\n        if type(feature) == list:\n            features.extend(feature)\n        else:\n            features.append(feature)\n    kwargs['features'] = features\n    return kwargs\n\ndef fit_transit_poly(**kwargs):\n    sig = kwargs['signal'][kwargs['ingress_right']+0:kwargs['egress_left']-0] + 1\n    x_fitted = np.arange(len(sig))\n    extend_points=200\n    \n    x_fitted, sig2 = extend_signal_with_polyfit(x_fitted, sig, extend_points, deg=2)\n    x_fitted = np.arange(len(sig2)) \n    coefficients = np.polyfit(x_fitted, sig2, 3)\n    fitted_curve = np.polyval(coefficients, x_fitted)\n    fitted_curve = fitted_curve[extend_points:-extend_points]\n    x_fitted = np.arange(len(fitted_curve)) + result['ingress_right']\n    \n    kwargs['transit_poly_coeffs'] = coefficients\n    kwargs['transit_poly_signal'] = fitted_curve\n    return kwargs\n\ndef collect_features(**kwargs):\n    transit_start = kwargs['ingress_right']\n    transit_finish = kwargs['egress_left']\n    signal = kwargs['signal'] + 1\n    signal_transit = signal[transit_start:transit_finish]\n    transit_poly_signal = kwargs['transit_poly_signal']\n    all_s = [] \n    all_s.append(signal_transit.mean())\n    all_s.append(kwargs['transit_poly_coeffs'][-1])\n    all_s.append(transit_poly_signal.mean())\n    all_s.append(transit_poly_signal[len(transit_poly_signal) // 2])\n    all_s.append(signal_transit.min())\n    all_s.append(signal_transit.max())\n    percentiles = np.percentile(signal_transit, np.arange(35, 65, 5))\n    all_s.extend(percentiles)    \n    \n    td_features = []\n    for s in all_s:\n        x = transit_formula(s)\n        td_features.append(x)\n\n    td_features.append(signal_transit.std())\n    td_features.extend(kwargs['transit_poly_coeffs'].tolist())\n    kwargs['td_features'] = td_features\n\n    other_features = []\n    ingress_signal = kwargs['signal'][kwargs['ingress_left']:kwargs['ingress_right']] + 1\n    egress_signal = kwargs['signal'][kwargs['egress_left']:kwargs['egress_right']] + 1\n    other_features.extend([ingress_signal.max(), ingress_signal.min(), egress_signal.max(), egress_signal.min()])\n    other_features.append(signal_transit.mean())\n    other_features.append(kwargs['transit_poly_coeffs'][-1])\n    other_features.append(transit_poly_signal.mean())\n    other_features.append(transit_poly_signal[len(transit_poly_signal) // 2])\n    kwargs['other_features'] = other_features\n    \n    return kwargs\n\n\ndef collect_wavebin_features(**kwargs):\n    transit_start = kwargs['ingress_right']\n    transit_finish = kwargs['egress_left']\n    signal = kwargs['signal'] + 1\n    signal_transit = signal[transit_start:transit_finish]\n    transit_poly_signal = kwargs['transit_poly_signal']\n    all_s = [] \n    all_s.append(signal_transit.mean())\n    all_s.append(kwargs['transit_poly_coeffs'][-1])\n    all_s.append(transit_poly_signal.mean())\n    all_s.append(transit_poly_signal[len(transit_poly_signal) // 2])\n    features = []\n    for s in all_s:\n        x = transit_formula(s)\n        features.append(x)\n\n    kwargs['td_features'] = features\n    \n    return kwargs\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"basic_pipeline_phases = [\n    lambda  **kwargs: bin_wavelengths(**kwargs),\n    lambda  **kwargs: normalize_signal(**kwargs),\n    lambda **kwargs: filter_signal_with_gaussian(sigma=25, **kwargs),\n    lambda  **kwargs: get_transit_phases2(**kwargs),\n#     lambda **kwargs: polynomial_detrend(**kwargs),\n#     lambda **kwargs: fit_transit_poly(**kwargs),\n#     lambda **kwargs: collect_features(**kwargs),\n]\n\nfull_pipeline = [\n    lambda  **kwargs: bin_wavelengths(**kwargs),\n    lambda  **kwargs: normalize_signal(**kwargs),\n    lambda **kwargs: filter_signal_with_gaussian(sigma=25, **kwargs),\n    lambda **kwargs: polynomial_detrend(**kwargs),\n    lambda **kwargs: fit_transit_poly(**kwargs),\n    lambda **kwargs: collect_features(**kwargs),\n]\n\nwavebins_pipeline = [\n    lambda  **kwargs: bin_wavelengths(**kwargs),\n    lambda  **kwargs: normalize_signal(**kwargs),\n    lambda **kwargs: filter_signal_with_gaussian(sigma=25, **kwargs),\n    lambda **kwargs: polynomial_detrend(**kwargs),\n    lambda **kwargs: fit_transit_poly(**kwargs),\n    lambda **kwargs: collect_wavebin_features(**kwargs),\n]\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\n# initial_kwargs_39 = {'signal': None, 'margin': 1, 'from_wavelength':39, 'to_wavelength':321}\n# initial_kwargs_55 = {'signal': None, 'margin': 1, 'from_wavelength':55, 'to_wavelength':300}\n\n# pipelines = []\n# pipelines.append((full_pipeline, initial_kwargs_39))\n# pipelines.append((full_pipeline, initial_kwargs_55))\n\n# binning = 6\n# for i in range(282):\n#     from_wavelength = 39+i-binning//2\n#     to_wavelength = 39+i+binning//2\n#     initial_kwargs = initial_kwargs_39.copy()\n#     initial_kwargs['from_wavelength'] = from_wavelength\n#     initial_kwargs['to_wavelength'] = to_wavelength\n\n#     pipeline = (wavebins_pipeline, initial_kwargs)\n#     pipelines.append(pipeline)\n        \n# print(len(pipelines))\n\n\ninitial_kwargs_39 = {'signal': None, 'margin': 1, 'from_wavelength':39, 'to_wavelength':321}\ninitial_kwargs_55 = {'signal': None, 'margin': 1, 'from_wavelength':55, 'to_wavelength':300}\n\npipelines = []\npipelines.append((full_pipeline, initial_kwargs_39))\npipelines.append((full_pipeline, initial_kwargs_55))\n\nfor binning in [6, 10, 25]:\n    for i in range(282):\n        from_wavelength = max(0,     39 + i - binning // 2)\n        to_wavelength =   min(356-1, 39 + i + binning // 2)\n        initial_kwargs = initial_kwargs_39.copy()\n        initial_kwargs['from_wavelength'] = from_wavelength\n        initial_kwargs['to_wavelength'] = to_wavelength\n\n        pipeline = (wavebins_pipeline, initial_kwargs)\n        pipelines.append(pipeline)\n        \nprint(len(pipelines))","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"n_planets = len(adc_info)\nn_wavelengths = 282\nall_s = []\nall_td_features = []\nall_other_features = []\n\n# n_planets = 2 # REMOVE\n\nfor planet_id in tqdm(range(n_planets)):\n    signal = preproc_airs(int(planet_ids[planet_id]), adc_info, dataset=mode)\n#     initial_kwargs = initial_kwargs_39.copy()\n#     initial_kwargs['signal'] = signal\n#     result = reduce(lambda kw, func: func(**kw), basic_pipeline_phases, initial_kwargs)\n#     ingress_left, ingress_right, egress_left, egress_right = result['ingress_left'], result['ingress_right'], result['egress_left'], result['egress_right']\n\n#     planet_features = []\n#     for pipeline in pipelines:\n#         initial_kwargs = pipeline[1].copy()\n#         initial_kwargs['signal'] = signal\n#         initial_kwargs['ingress_left'], initial_kwargs['ingress_right'], initial_kwargs['egress_left'], initial_kwargs['egress_right'] = ingress_left, ingress_right, egress_left, egress_right\n#         result = reduce(lambda kw, func: func(**kw), pipeline[0], initial_kwargs)\n\n#         planet_features.extend(result['features'])\n\n#     all_features.append(planet_features)\n\n    initial_kwargs = initial_kwargs_39.copy()\n    initial_kwargs['signal'] = signal\n    result = reduce(lambda kw, func: func(**kw), basic_pipeline_phases, initial_kwargs)\n    ingress_left, ingress_right, egress_left, egress_right = result['ingress_left'], result['ingress_right'], result['egress_left'], result['egress_right']\n\n    td_features = []\n    other_features = []\n    for pipeline in pipelines:\n        initial_kwargs = pipeline[1].copy()\n        initial_kwargs['signal'] = signal\n        initial_kwargs['ingress_left'], initial_kwargs['ingress_right'], initial_kwargs['egress_left'], initial_kwargs['egress_right'] = ingress_left, ingress_right, egress_left, egress_right\n        result = reduce(lambda kw, func: func(**kw), pipeline[0], initial_kwargs)\n\n        td_features.extend(result.get('td_features', []))\n        other_features.extend(result.get('other_features', []))\n\n    all_td_features.append(td_features)\n    all_other_features.append(other_features)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# X = np.array(all_features)\n# scaler = joblib.load('/kaggle/input/ad24-head-10/models/scaler.pkl')\n# X_scaled = scaler.transform(X)\n# X.shape\nX_td = np.array(all_td_features)\nX_other = np.array(all_other_features)\n\nscaler_td = joblib.load('/kaggle/input/ad24-head-11/models/scaler_td.pkl')\nscaler_other = joblib.load('/kaggle/input/ad24-head-11/models/scaler_other.pkl')\nX_td = scaler_td.transform(X_td)\nX_other = scaler_other.transform(X_other)\nX = np.concatenate([X_td, X_other], axis=1)\n\nX_td.shape, X_other.shape, X.shape\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def predict_with_models(k, X, directory='/kaggle/input/ad24-head-11/models/'):\n    all_predictions = []\n    \n    for i in tqdm(range(k)):\n        model = None\n        with open(f'{directory}model_{i}.pkl', 'rb') as file:\n            model = joblib.load(file)\n        preds = model.predict(X)\n        all_predictions.append(preds)\n        del model\n        gc.collect()\n\n    all_predictions = np.array(all_predictions)\n\n    final_predictions = []\n    std_devs = []\n\n    for i in range(all_predictions.shape[1]): \n        sample_preds = all_predictions[:, i, :]\n        sorted_preds = np.sort(sample_preds, axis=0)\n\n        MARGINS = int(k * 0.30)\n        filtered_preds = sorted_preds[MARGINS:-MARGINS] if sorted_preds.shape[0] > MARGINS*2+3 else sorted_preds\n        \n        mean_preds = np.mean(filtered_preds, axis=0)\n        std_preds = np.std(sorted_preds, axis=0)\n\n        final_predictions.append(mean_preds)\n        std_devs.append(std_preds)\n\n    return np.array(final_predictions), np.array(std_devs)\n\ntransit_depths, sigma = predict_with_models(min(1000, 7000), X)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sample_submission = pd.read_csv('/kaggle/input/ariel-data-challenge-2024/sample_submission.csv')\npred = transit_depths.clip(0) \n\n# linear_distribution = np.linspace(0.00016, 0.00013, transit_depths.shape[1])\n# sigma_response = np.tile(linear_distribution, (transit_depths.shape[0], 1))\n# percentage_distribution = np.linspace(0.022, 0.022, pred.shape[1]) # Recommended by calculations\n# sigma = percentage_distribution * pred\n\nCOEFF_0 = 5 # as predicted by the model\nsigma_response = sigma * COEFF_0\n\nsubmission = pd.DataFrame(np.concatenate([pred, sigma_response.clip(0)], axis=1), columns=sample_submission.columns[1:])\nsubmission.index = adc_info.index\n# submission.to_csv('submission.csv')\nsubmission.to_csv('submission.csv', float_format='%.8f')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# !cat submission.csv","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{"trusted":true},"execution_count":null,"outputs":[]}]}