{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.11.11","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"none","dataSources":[{"sourceId":94689,"databundleVersionId":11605086,"sourceType":"competition"},{"sourceId":11800719,"sourceType":"datasetVersion","datasetId":7410737},{"sourceId":11850587,"sourceType":"datasetVersion","datasetId":7446274}],"dockerImageVersionId":31041,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# Foraminifera classification\n\nThe task of this competition is to classify volumetric scans of foraminifera into 14 types, plus a 15th 'unknown' class, which includes scans of foram fragments, sand, or other similar objects. Each image is 128-by-128-by-128 voxels.\n\nKey points:\n- Every 3d input image is converted to six surface views. The six views are rotated and flipped (data augmentation) and then embedded in a 2048-dimensional space by Resnet-50.\n- We apply the same preprocessing and augmentation to the test data (test-time augmentation).\n- The notebook pseudolabels all samples for which four other models agree in their predictions. The pseudolabels are used for training a classifier with 14 classes.\n- The 5 % samples for which the classifier is least confident are labeled as class 14.\n\nReferences:\n- [Forams Classification 2025](https://www.kaggle.com/competitions/forams-classification-2025)\n- [Foraminifera on Wikipedia](https://en.wikipedia.org/wiki/Foraminifera)\n- [Resnet-50 on Tensorflow](https://www.tensorflow.org/api_docs/python/tf/keras/applications/ResNet50)","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19"}},{"cell_type":"code","source":"import pandas as pd\nimport numpy as np\nimport matplotlib\nimport matplotlib.pyplot as plt\nimport matplotlib.image as mpimg\nimport seaborn as sns\nfrom glob import glob\nimport PIL\nfrom tifffile import imread\nfrom tqdm import tqdm\nimport math\nimport pickle\n\nfrom sklearn.pipeline import make_pipeline\nfrom sklearn.preprocessing import StandardScaler, FunctionTransformer, LabelBinarizer\nfrom sklearn.linear_model import LogisticRegression\nfrom sklearn.semi_supervised import LabelPropagation\nfrom sklearn.model_selection import cross_val_score, cross_val_predict, StratifiedGroupKFold\nfrom sklearn.metrics import f1_score, confusion_matrix, log_loss, accuracy_score, classification_report\n\nfrom tensorflow import __version__\nimport keras\nimport keras.applications.resnet50 as resnet50\nfrom keras.utils import image_dataset_from_directory\nfrom keras.layers import Dense, Activation, Flatten, Dropout, Concatenate\nfrom keras.models import Sequential, Model \nfrom keras.optimizers import Adam\nfrom keras.callbacks import ReduceLROnPlateau, EarlyStopping, TerminateOnNaN\n\nprint(f\"Tensorflow version: {__version__}\")\n\nFROM_SCRATCH = False\nPREPROCESSED_DIR = \"/kaggle/input/forams-lighted-from-front\" # lighted from front\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-05-24T08:33:58.212940Z","iopub.execute_input":"2025-05-24T08:33:58.213745Z","iopub.status.idle":"2025-05-24T08:34:15.819645Z","shell.execute_reply.started":"2025-05-24T08:33:58.213716Z","shell.execute_reply":"2025-05-24T08:34:15.818808Z"},"_kg_hide-input":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Preprocessing\n\nWe read the 210 labelled 3d images and process them as follows:\n- Convert the grayscale 3d image to a black-and-white (boolean) 3d image\n- Create six 127*127 grayscale surface views of the 3d image\n- Add two other features (scaling factor and mass of object)\n- Augment the data 8-fold by rotating and mirroring the surfaces\n- Embed the 127*127 surfaces into a 2048-dimensional space with Resnet-50","metadata":{}},{"cell_type":"code","source":"if FROM_SCRATCH:\n    HEIGHT = 127\n    WIDTH = 127\n    \n    resnet50_model = resnet50.ResNet50(include_top=False,\n                                       weights='imagenet',\n                                       pooling='avg',\n                                       input_shape=(HEIGHT, WIDTH, 3))","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-05-24T08:34:15.821192Z","iopub.execute_input":"2025-05-24T08:34:15.821783Z","iopub.status.idle":"2025-05-24T08:34:15.826343Z","shell.execute_reply.started":"2025-05-24T08:34:15.821761Z","shell.execute_reply":"2025-05-24T08:34:15.825573Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"%%time\n# Read the 210 true labels\nlabels = pd.read_csv(f'/kaggle/input/forams-classification-2025/labelled.csv', index_col='id')['label'].values # 210 rows, id and label\n\nif FROM_SCRATCH:\n    # Read the 3d images and convert every 3d image to\n    # 24 2d surface images\n    \n    linear_f = np.arange(128) + 100\n    linear_b = linear_f[::-1]\n    \n    def make_feature_array(vol_list, labels=None, group0=0, plot=False):\n        \"\"\"Create six surface views for every 3d sample in vol_list.\"\"\"\n        X, more_features, groups, y = [], [], [], []\n    \n        def add_surface(surface, i, label):\n            \"\"\"Add this surface to the list.\n            \n            The input surface is a 128*128 array of heights, which will\n            be converted into a 127*127 grayscale image.\"\"\"\n            surface = surface.astype(np.float32)\n            grad = np.sqrt(np.square(surface[1:,1:] - surface[:-1,1:]) + np.square(surface[1:,1:] - surface[1:,:-1]))\n            grad[surface[1:,1:] == 0] = 0\n            grad[surface[:-1,1:] == 0] = 0\n            grad[surface[1:,:-1] == 0] = 0\n            grad = 6 - grad.clip(0, 6)\n            grad *= 42\n            X.append(np.hstack([grad]))\n            more_features.append((scaling_factor, mass))\n            groups.append(group0 + i)\n            y.append(label)\n            \n        for i, filename in enumerate(vol_list):\n            volume = imread(filename)\n    \n            # Make black-white for best contrast and so that we can define surfaces\n            volume_bw = volume > 110\n    \n            # Find total mass of object\n            mass = volume_bw.sum()\n            mass = (mass - 162000) / 56662\n    \n            # Scaling factor taken from filename\n            scaling_factor = float(filename[filename.index('_sc_') + 4 : -4].replace('_', '.'))\n    \n            # Extract the six surface views\n            # A surface is a 128*128 array of heights\n            label = labels[i] if labels is not None else None\n            add_surface((volume_bw * linear_f.reshape(-1, 1, 1)).max(axis=0), i, label)\n            add_surface((volume_bw * linear_b.reshape(-1, 1, 1)).max(axis=0), i, label)\n            add_surface((volume_bw * linear_f.reshape(1, -1, 1)).max(axis=1), i, label)\n            add_surface((volume_bw * linear_b.reshape(1, -1, 1)).max(axis=1), i, label)\n            add_surface((volume_bw * linear_f.reshape(1, 1, -1)).max(axis=2), i, label)\n            add_surface((volume_bw * linear_b.reshape(1, 1, -1)).max(axis=2), i, label)\n    \n            # Plot\n            if plot and i == 0:\n                _, axs = plt.subplots(1, 6, figsize=(7, 1.25))\n                axs = axs.ravel()\n                for j in range(len(axs)):\n                    axs[j].imshow(X[-1-j], cmap='gray')\n                    axs[j].axis('off')\n                plt.show()\n    \n        # Convert grayscale to rgb and center\n        X = np.array(X) # shape (n_samples, 127, 127)\n        X = np.repeat(X, 3).reshape(X.shape + (3, )) # shape (n_samples, 127, 127, 3)\n        X = resnet50.preprocess_input(X)\n    \n        return (X, \n                np.array(more_features).astype(np.float32),\n                np.array(groups),\n                np.array(y) if labels is not None else None)\n    \n    def augment(X, more_features, groups, y):\n        \"\"\"Augment the data 8-fold\"\"\"\n        # Mirror all images\n        X = np.vstack([X, np.flip(X, axis=2)])\n        more_features = np.tile(more_features, (2, 1))\n        groups = np.tile(groups, 2)\n        y = np.tile(y, 2)\n        \n        # Rotate all images by 180 degrees\n        X = np.vstack([X, X[:, ::-1, ::-1, :]])\n        more_features = np.tile(more_features, (2, 1))\n        groups = np.tile(groups, 2)\n        y = np.tile(y, 2)\n        \n        # Rotate all images by 90 degrees\n        X = np.vstack([X, np.rot90(X, axes=(1, 2))])\n        more_features = np.tile(more_features, (2, 1))\n        groups = np.tile(groups, 2)\n        y = np.tile(y, 2)\n    \n        return X, more_features, groups, y\n    \n    # Read and augment the labelled data\n    vol_list_labelled = sorted(glob('/kaggle/input/forams-classification-2025/volumes/volumes/labelled/*.tif'))\n    X_labelled, more_features, groups, y = make_feature_array(vol_list_labelled, labels=labels, plot=True)\n    X_labelled, more_features, groups, y = augment(X_labelled, more_features, groups, y)\n    \n    X_labelled.shape, more_features.shape, groups.shape, y.shape\n    # ((10080, 127, 127, 3), (10080, 2), (10080,), (10080))\n    # 20 s for all 210 labelled forams\n    \n    X_labelled = resnet50_model.predict(X_labelled, batch_size=256, verbose=1)\n    # 6 minutes for 10080 training samples on cpu\n    \n    # Save the features\n    with open(\"X_labelled.pickle\", \"wb\") as f:\n        pickle.dump((X_labelled, more_features, groups, y), f)\n\nelse:\n    with open(PREPROCESSED_DIR + \"/X_labelled.pickle\", \"rb\") as f:\n        X_labelled, more_features, groups, y = pickle.load(f)\n\ngroups = - groups - 1\n\n# Binarize the labels\nlb = LabelBinarizer()\ny_b = lb.fit_transform(y)\n\nX_labelled.shape, more_features.shape, groups.shape, y_b.shape\n# ((10080, 127, 127, 3), (10080, 2), (10080,), (10080, 14))\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-05-24T08:45:46.389731Z","iopub.execute_input":"2025-05-24T08:45:46.390249Z","iopub.status.idle":"2025-05-24T08:45:48.384754Z","shell.execute_reply.started":"2025-05-24T08:45:46.390220Z","shell.execute_reply":"2025-05-24T08:45:48.383947Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Pseudolabeling\n\nWe read the predictions of some other notebooks and determine the trustworthy predictions, i.e., the samples for which all notebooks agree on the classification.","metadata":{}},{"cell_type":"code","source":"predictions = pd.read_csv('/kaggle/input/forams-predictions/forams-predictions.csv', \n                          sep=';',\n                          index_col='id')\npredictions['all_equal'] = predictions.var(axis=1) == 0\npredictions['pseudolabel'] = np.where(predictions.all_equal, predictions['resnet-v4'], -1)\nprint(f\"Agreement: {predictions.all_equal.mean():.0%}\")\npredictions","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-05-24T08:34:17.542927Z","iopub.execute_input":"2025-05-24T08:34:17.543315Z","iopub.status.idle":"2025-05-24T08:34:17.593062Z","shell.execute_reply.started":"2025-05-24T08:34:17.543290Z","shell.execute_reply":"2025-05-24T08:34:17.592212Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"Then we apply to the unlabelled/pseudolabelled data the same preprocessing as for the labelled data:","metadata":{}},{"cell_type":"code","source":"%%time\nif FROM_SCRATCH:\n    # Read and augment the unlabelled data\n    vol_list_ul = sorted(glob('/kaggle/input/forams-classification-2025/volumes/volumes/unlabelled/*.tif'))\n    X_ul_list, more_features_ul_list, groups_ul_list, y_ul_list = [], [], [], [] \n    test_batch_size = 256\n    for batch_start in tqdm(range(0, len(vol_list_ul), test_batch_size)):\n        X_ul, more_features_ul, groups_ul, y_ul = make_feature_array(\n            vol_list_ul[batch_start:batch_start+test_batch_size],\n            labels=predictions.pseudolabel.values[batch_start:batch_start+test_batch_size],\n            group0=batch_start,\n            plot=False)\n        X_ul, more_features_ul, groups_ul, y_ul = augment(X_ul, more_features_ul, groups_ul, y_ul)\n        X_ul = resnet50_model.predict(X_ul, batch_size=256, verbose=0)\n        X_ul_list.append(X_ul)\n        more_features_ul_list.append(more_features_ul)\n        groups_ul_list.append(groups_ul)\n        y_ul_list.append(y_ul)\n    \n    X_ul = np.vstack(X_ul_list)\n    more_features_ul = np.vstack(more_features_ul_list)\n    groups_ul = np.hstack(groups_ul_list)\n    y_ul = np.hstack(y_ul_list)\n    \n    del X_ul_list, more_features_ul_list, groups_ul_list, y_ul_list\n\n    # Save the features (7 GByte)\n    with open(\"X_pseudolabelled.pickle\", \"wb\") as f:\n        pickle.dump((X_ul, more_features_ul, groups_ul, y_ul), f)\n\nelse:\n    with open(PREPROCESSED_DIR + \"/X_pseudolabelled.pickle\", \"rb\") as f:\n        X_ul, more_features_ul, groups_ul, y_ul = pickle.load(f) # 7 GByte\n    \nprint('Both:', X_ul.shape, more_features_ul.shape, groups_ul.shape, y_ul.shape)\n\n# Binarize the labels\ny_ul_b = lb.transform(y_ul)\n\n# Separate pseudolabelled and unlabelled datasets\nX_pl = X_ul[y_ul >= 0]\nX_ul = X_ul[y_ul < 0]\nmore_features_pl = more_features_ul[y_ul >= 0]\nmore_features_ul = more_features_ul[y_ul < 0]\ngroups_pl = groups_ul[y_ul >= 0]\ngroups_ul = groups_ul[y_ul < 0]\ny_pl_b = y_ul_b[y_ul >= 0]\ny_ul_b = y_ul_b[y_ul < 0]\ny_pl = y_ul[y_ul >= 0]\ny_ul = y_ul[y_ul < 0] # always -1, must be assigned last\nprint('Pseudolabelled:', X_pl.shape, more_features_pl.shape, groups_pl.shape, y_pl.shape)\nprint('Unlabelled:', X_ul.shape, more_features_ul.shape, groups_ul.shape, y_ul.shape)\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-05-24T08:34:17.594005Z","iopub.execute_input":"2025-05-24T08:34:17.594703Z","iopub.status.idle":"2025-05-24T08:35:33.586325Z","shell.execute_reply.started":"2025-05-24T08:34:17.594679Z","shell.execute_reply":"2025-05-24T08:35:33.584984Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# The model","metadata":{}},{"cell_type":"code","source":"DROPOUT = 0.5\nFC_LAYERS = [256] # [256]\nN_CLASSES = 14\n\ndef build_model(dropout, fc_layers, n_classes):\n    \"\"\"Some fully connected layers and a final softmax layer\"\"\"\n    embedding = keras.layers.Input(shape=(2048, ), name='embedding')\n    more_features = keras.layers.Input(shape=(2, ), name='more_features')\n    x = Concatenate()([embedding, more_features])\n    for fc in fc_layers:\n        print(fc)\n        x = Dense(fc, activation='relu')(x)\n        x = Dropout(dropout)(x)\n    outputs = Dense(n_classes, activation='softmax')(x)\n    model = Model(inputs={'embedding': embedding, 'more_features': more_features},\n                           outputs=outputs)\n    return model\n\n# build_model(dropout=DROPOUT, fc_layers=FC_LAYERS, n_classes=N_CLASSES).summary()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-05-24T08:35:33.588213Z","iopub.execute_input":"2025-05-24T08:35:33.588562Z","iopub.status.idle":"2025-05-24T08:35:33.596883Z","shell.execute_reply.started":"2025-05-24T08:35:33.588535Z","shell.execute_reply":"2025-05-24T08:35:33.595925Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"# Cross-validation","metadata":{}},{"cell_type":"code","source":"%%time\n# Cross-validate\n\nBATCH_SIZE = 512\nINFERENCE = True\nEARLY_STOPPING = True\nif EARLY_STOPPING:\n    NUM_EPOCHS = 100\nelse:\n    NUM_EPOCHS = 6\nVERBOSE = 0\n\ncv = StratifiedGroupKFold(n_splits=15)\noof_prob = np.full((len(y), N_CLASSES), np.nan)\nepoch_sum = 0\nif INFERENCE:\n    ul_probs = []\n    \nfor fold, (idx_tr, idx_va) in enumerate(cv.split(X_labelled, y, groups=groups)):\n    # Clean up the memory\n    try:\n        del history\n    except NameError:\n        pass\n    try:\n        del my_model\n    except NameError:\n        pass\n\n    # Split into training and validation data\n    X_tr = X_labelled[idx_tr]\n    more_features_tr = more_features[idx_tr]\n    y_tr = y_b[idx_tr]\n    X_va = X_labelled[idx_va]\n    more_features_va = more_features[idx_va]\n    y_va = y_b[idx_va]\n\n    # Add the pseudolabelled data for training\n    X_tr = np.vstack([X_tr, X_pl])\n    more_features_tr = np.vstack([more_features_tr, more_features_pl])\n    y_tr = np.vstack([y_tr, y_pl_b])\n    \n    # Build, compile and fit the model\n    my_model = build_model(dropout=DROPOUT,\n                           fc_layers=FC_LAYERS,\n                           n_classes=N_CLASSES)\n    my_model.compile(optimizer=keras.optimizers.Adam(learning_rate=0.003),\n                     loss=keras.losses.CategoricalCrossentropy(),\n                     metrics=[\"categorical_accuracy\", \"f1_score\"]\n                    )\n    lr = ReduceLROnPlateau(monitor=\"val_loss\", factor=0.7, \n                           patience=4, verbose=VERBOSE)\n    es = EarlyStopping(monitor=\"val_loss\",\n                       patience=10, \n                       verbose=1,\n                       mode=\"min\", \n                       restore_best_weights=True)\n    callbacks = [lr, es, TerminateOnNaN()] if EARLY_STOPPING else [TerminateOnNaN()]\n    history = my_model.fit({'embedding': X_tr, 'more_features': more_features_tr},\n                           y_tr,\n                           batch_size=BATCH_SIZE,\n                           epochs=NUM_EPOCHS,\n                           # steps_per_epoch=3,\n                           validation_data=({'embedding': X_va, 'more_features': more_features_va}, y_va),\n                           validation_batch_size=1024,\n                           verbose=VERBOSE,\n                           callbacks=callbacks)\n    epoch_sum += np.argmin(history.history['val_loss'])\n    \n    # Plot history\n    plt.figure(figsize=(12, 3))\n    plt.plot(np.arange(len(history.history['loss'])) + 1,\n             history.history['loss'],\n             label='train loss')\n    plt.plot(np.arange(len(history.history['val_loss'])) + 1,\n             history.history['val_loss'],\n             label='val loss')\n    plt.title(f'Training history fold {fold}')\n    plt.xlabel('epoch')\n    plt.ylabel('loss')\n    plt.ylim(0, 2)\n    plt.grid(True, axis='y')\n    plt.legend()\n    plt.show()\n\n    # Clean up the memory\n    del lr, es, callbacks, history\n\n    # Compute oof probabilities\n    oof_prob[idx_va] = my_model.predict({'embedding': X_va, 'more_features': more_features_va},\n                                        verbose=0,\n                                        batch_size=1024)\n\n    if INFERENCE:\n        ul_prob = my_model.predict({'embedding': X_ul, 'more_features': more_features_ul}, \n                                  verbose=0,\n                                  batch_size=1024)\n        # ul_prob_grouped = pd.DataFrame(ul_prob).groupby(groups_ul).mean().values # average test-time augmentation\n        ul_probs.append(ul_prob) # add to final ensemble\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-05-24T08:35:33.598185Z","iopub.execute_input":"2025-05-24T08:35:33.598525Z","iopub.status.idle":"2025-05-24T08:44:47.923686Z","shell.execute_reply.started":"2025-05-24T08:35:33.598494Z","shell.execute_reply":"2025-05-24T08:44:47.922766Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"PERCENT_OUTLIERS = 5 # hyperparameter: how many samples do we want to qualify as outliers?\n\nul_prob = np.stack(ul_probs).mean(axis=0) # mean over predictions of five folds (513168, 14)\nul_prob_grouped = pd.DataFrame(ul_prob).groupby(groups_ul).mean().values\nul_pred = np.argmax(ul_prob_grouped, axis=1)\nCLASS_14_THRESHOLD = np.quantile(\n    ul_prob_grouped.max(axis=1),\n    q=PERCENT_OUTLIERS/100/len(ul_prob_grouped) * 18216\n)\nul_pred[ul_prob_grouped.max(axis=1) < CLASS_14_THRESHOLD] = 14 # \"unknown\" class\n\n# Evaluate the results with test-time augmentation\noof_prob_grouped = pd.DataFrame(oof_prob).groupby(groups).mean().values\ny_grouped = pd.Series(y).groupby(groups).mean().values\n\noof = np.argmax(oof_prob_grouped, axis=1)\noof[oof_prob_grouped.max(axis=1) < CLASS_14_THRESHOLD] = 14 # \"unknown\" class\nprint(f\"# F1: {f1_score(y_grouped, oof, labels=np.arange(15), average='macro', zero_division=0):.3f}\", end='   ')\nprint(f\"Acc: {accuracy_score(y_grouped, oof):.3f}\", end='   ')\nprint(f\"Logloss: {log_loss(y_grouped, oof_prob_grouped):.3f} {PERCENT_OUTLIERS}%\")\n\nsns.heatmap(confusion_matrix(y_grouped, oof, labels=np.arange(15)), annot=True, fmt='.0f')\nplt.title(f'Resnet-50 with TTA{\" with early stopping\" if EARLY_STOPPING else \"\"}')\nplt.show()\nprint(classification_report(y_grouped, oof, zero_division=0, labels=np.arange(15)))\n\n# Best epoch\nprint(f\"Best epoch: {epoch_sum / 5:.0f}\")\n\n# Save the oof probabilities\nwith open(\"oof_prob.pickle\", \"wb\") as f:\n    pickle.dump((oof_prob, groups, y), f)\n\n# Save the test probabilities\nwith open(\"ul_prob.pickle\", \"wb\") as f:\n    pickle.dump(ul_prob, f)\n\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-05-24T08:59:38.166713Z","iopub.execute_input":"2025-05-24T08:59:38.167602Z","iopub.status.idle":"2025-05-24T08:59:39.043401Z","shell.execute_reply.started":"2025-05-24T08:59:38.167572Z","shell.execute_reply":"2025-05-24T08:59:39.042512Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"if INFERENCE:\n    # Write the submission file\n    y_pred = predictions.pseudolabel.values.copy()\n    y_pred[y_pred < 0] = ul_pred\n    submission = pd.Series(\n        y_pred,\n        index=pd.RangeIndex(0, len(predictions), name='id'),\n        name='label')\n    submission.to_csv('submission.csv')\n    !head submission.csv\n\n    # Frequency diagram\n    plt.figure(figsize=(12, 4))\n    plt.title('Predicted classes')\n    x, count = np.unique(ul_pred, return_counts=True)\n    _, pl_count = np.unique(predictions.pseudolabel.values, return_counts=True)\n    pl_count = np.array(list(pl_count[1:]) + [0])\n    b = plt.bar(x, count, bottom=pl_count, label='counted')\n    plt.bar(np.arange(14), ul_prob_grouped.sum(axis=0), bottom=pl_count[:-1], alpha=0.6, label='expected')\n    plt.bar_label(b)\n    plt.bar(np.arange(15), pl_count, color='lightgray', label='pseudolabelled')\n    plt.xticks(np.arange(15))\n    plt.xlabel('class')\n    plt.ylabel('count')\n    plt.legend()\n    plt.show()\n    \n    # Embedding into two dimensions diagram\n    tab20 = matplotlib.colormaps['tab20']\n    colors = [tab20.colors[np.round(x).astype(int)] for x in np.linspace(0, 19, 14)] + [(0, 0, 0)]\n    colors = np.array(colors) # shape (15, 3)\n    \n    components = np.vstack([np.cos(np.arange(14) / 14 * 2 * math.pi),\n                            np.sin(np.arange(14) / 14 * 2 * math.pi),\n                           ])\n    tt = oof_prob_grouped @ components.T # train\n    uu = ul_prob_grouped @ components.T # test\n    plt.figure(figsize=(12, 12))\n    plt.title('2d projection of softmax probabilities')\n    plt.scatter(tt[:,0], tt[:,1], s=30, c=colors[labels], marker='x', label='train')\n    plt.scatter(uu[:,0], uu[:,1], s=3, c=colors[ul_pred], label='test')\n    plt.scatter(components[0], components[1], s=100, c='k', alpha=0.3)\n    for cl in range(14):\n        plt.text(components[0, cl] * 1.06, components[1, cl] * 1.06, cl, ha='center', va='center_baseline')\n    plt.xticks([])\n    plt.yticks([])\n    plt.legend()\n    plt.show()\n\nprint((oof==14).mean(), (y_pred==14).mean())","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-05-24T08:59:59.391434Z","iopub.execute_input":"2025-05-24T08:59:59.391829Z","iopub.status.idle":"2025-05-24T09:00:00.321368Z","shell.execute_reply.started":"2025-05-24T08:59:59.391800Z","shell.execute_reply":"2025-05-24T09:00:00.320411Z"},"_kg_hide-input":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"","metadata":{"trusted":true},"outputs":[],"execution_count":null}]}