{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# Finding duplicate recordings using embeddings\n\nThis notebook analyzes the embedding distance between each audio recording in the BirdCLEF 2023 dataset in order to find duplicate recordings. \"embedding distance\" is the euclidean distance between the embeddings of the first and last 5 seconds of the audio (concatenated). The embeddings were precomputed using the Google Bird Vocalization Classifier [here](https://www.kaggle.com/code/robbynevels/google-bird-model-embeddings-predictions-score?scriptVersionId=125329121). The logic below finds 78 duplicates using distance thresholds + noise filtering + audio magnitude threshold + some manual filtering, making a total of 79 known duplicates as of this writing. Check out some [lhanhsin](https://www.kaggle.com/competitions/birdclef-2023/discussion/398229)'s post listing all duplicates, and [matt op](https://www.kaggle.com/code/mattop/birdclef-2023-eda)'s work on this as well.\n\nI've saved the 16941x16941 distance matrix in a `distances.npy` file which you can find in the output of this notebook. It takes a while to compute without a GPU, so feel free to use it if you're running low on GPU hours :)\n\nAlso, if you're interested in exploring the embeddings themselves, you might find [my exploration of projecting them into 3D space](https://www.kaggle.com/code/robbynevels/birdclef2023-eda-with-3d-embeddings#Musing) inspiring. I used the noise analysis in that notebook to filter out many recordings in this one.","metadata":{}},{"cell_type":"code","source":"import pandas as pd\nfrom pathlib import Path\nimport torch\nimport matplotlib.pyplot as plt\nimport numpy as np\nfrom tqdm import tqdm\nimport plotly.express as px","metadata":{"execution":{"iopub.status.busy":"2023-04-11T03:24:13.082450Z","iopub.execute_input":"2023-04-11T03:24:13.083528Z","iopub.status.idle":"2023-04-11T03:24:13.088783Z","shell.execute_reply.started":"2023-04-11T03:24:13.083457Z","shell.execute_reply":"2023-04-11T03:24:13.087439Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df = pd.read_csv('/kaggle/input/birdclef-2023/train_metadata.csv')\n\nTRAIN_PATH = Path('/kaggle/input/birdclef-2023/train_audio')\nDEVICE = 'cuda' if torch.cuda.is_available() else 'cpu'","metadata":{"execution":{"iopub.status.busy":"2023-04-11T03:24:13.249815Z","iopub.execute_input":"2023-04-11T03:24:13.250928Z","iopub.status.idle":"2023-04-11T03:24:13.312525Z","shell.execute_reply.started":"2023-04-11T03:24:13.250882Z","shell.execute_reply":"2023-04-11T03:24:13.311499Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Fetch precomputed embeddings from Google Bird Vocalization Classifier\n# combine first and last embedding of all recordings (corresponds to first and last 5 seconds of audio)\nembeddings = torch.load('/kaggle/input/gbvc-embeddings/embeddings.pt')\nfirst5sec_embeddings = torch.stack([torch.tensor(embeddings[filename][0]) for filename in df.filename])\nlast5sec_embeddings = torch.stack([torch.tensor(embeddings[filename][-1]) for filename in df.filename])\ncombined_embeddings = torch.concat([first5sec_embeddings, last5sec_embeddings], -1)\ncombined_embeddings.shape","metadata":{"execution":{"iopub.status.busy":"2023-04-11T03:24:13.385267Z","iopub.execute_input":"2023-04-11T03:24:13.385624Z","iopub.status.idle":"2023-04-11T03:24:29.543298Z","shell.execute_reply.started":"2023-04-11T03:24:13.385592Z","shell.execute_reply":"2023-04-11T03:24:29.542073Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# compute euclidean distances between each cell of a 2d matrix\ndef get_distances(x):\n    result = []\n    # loop over each row, since doing it all at once would cause a memory error\n    for row in tqdm(x):\n        result.append(((row[:, None] - x.T) ** 2).T.sum(-1))\n    return torch.stack(result)\n\ndistances = get_distances(combined_embeddings.to(DEVICE))","metadata":{"execution":{"iopub.status.busy":"2023-04-11T03:24:29.787289Z","iopub.execute_input":"2023-04-11T03:24:29.787642Z","iopub.status.idle":"2023-04-11T03:24:59.182548Z","shell.execute_reply.started":"2023-04-11T03:24:29.787612Z","shell.execute_reply":"2023-04-11T03:24:59.181099Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# save distances so it can be analyzed by others\nnp.save('distances.npy', distances.to('cpu').numpy())","metadata":{"execution":{"iopub.status.busy":"2023-04-11T03:25:19.313771Z","iopub.execute_input":"2023-04-11T03:25:19.314173Z","iopub.status.idle":"2023-04-11T03:25:21.525163Z","shell.execute_reply.started":"2023-04-11T03:25:19.314139Z","shell.execute_reply":"2023-04-11T03:25:21.516550Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize=(9, 9))\nplt.title('All embedding distances (notice 0 along the diagonal, as expected)')\nplt.imshow(distances.to('cpu'))\nplt.colorbar();","metadata":{"execution":{"iopub.status.busy":"2023-04-11T03:26:10.141969Z","iopub.execute_input":"2023-04-11T03:26:10.142544Z","iopub.status.idle":"2023-04-11T03:26:20.605353Z","shell.execute_reply.started":"2023-04-11T03:26:10.142445Z","shell.execute_reply":"2023-04-11T03:26:20.604512Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Fetch outliers from a previous analysis on 3D projections of embeddings\noutliers = pd.read_csv('/kaggle/input/birdclef2023-eda-with-3d-embeddings/outliers.csv')\n\nfor outlier_type in [1, 2]:\n    outlier_distances = get_distances(combined_embeddings[df.filename.isin(set(outliers[outliers.outlier_type == outlier_type].filename))].to(DEVICE))\n    plt.title('outliers on line ' + str(outlier_type))\n    plt.imshow(outlier_distances.to('cpu'))\n    plt.colorbar()\n    plt.show()","metadata":{"execution":{"iopub.status.busy":"2023-04-11T03:26:35.231341Z","iopub.execute_input":"2023-04-11T03:26:35.232219Z","iopub.status.idle":"2023-04-11T03:26:35.893195Z","shell.execute_reply.started":"2023-04-11T03:26:35.232182Z","shell.execute_reply":"2023-04-11T03:26:35.891279Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Known duplicates, copied from https://www.kaggle.com/competitions/birdclef-2023/discussion/398229\nknown_dupes = \\\n    [['carcha1/XC324665.ogg', 'carcha1/XC324666.ogg'],\n     ['wbrcha2/XC613380.ogg', 'wbrcha2/XC613384.ogg'],\n     ['comsan/XC613127.ogg', 'comsan/XC613128.ogg'],\n     ['tafpri1/XC443724.ogg', 'tafpri1/XC443725.ogg'],\n     ['colsun2/XC755891.ogg', 'colsun2/XC755892.ogg'],\n     ['grewoo2/XC527938.ogg', 'grewoo2/XC527939.ogg'],\n     ['sccsun2/XC609477.ogg', 'sccsun2/XC609478.ogg'],\n     ['subbus1/XC603421.ogg', 'subbus1/XC603426.ogg'],\n     ['combul2/XC748220.ogg', 'combul2/XC748221.ogg'],\n     ['wtbeat1/XC234928.ogg', 'wtbeat1/XC234929.ogg'],\n     ['afrthr1/XC652880.ogg', 'afrthr1/XC652884.ogg'],\n     ['laudov1/XC405374.ogg', 'laudov1/XC405375.ogg'],\n     ['litegr/XC411319.ogg', 'litegr/XC411320.ogg'],\n     ['cohmar1/XC564020.ogg', 'cohmar1/XC564021.ogg'],\n     ['egygoo/XC358927.ogg', 'egygoo/XC528135.ogg'],\n     ['combuz1/XC647786.ogg', 'combuz1/XC647787.ogg'],\n     ['combuz1/XC144257.ogg', 'combuz1/XC144258.ogg'],\n     ['wlwwar/XC478705.ogg', 'wlwwar/XC478767.ogg'],\n     ['woosan/XC740798.ogg',  'woosan/XC742927.ogg'],\n     ['gnbcam2/XC530150.ogg', 'gnbcam2/XC530151.ogg'],\n     ['litswi1/XC443712.ogg', 'litswi1/XC443713.ogg'],\n     ['combul2/XC650878.ogg', 'combul2/XC447669.ogg'],\n     ['gobbun1/XC394478.ogg', 'gobbun1/XC395111.ogg'],\n     ['fislov1/XC503794.ogg', 'fislov1/XC526237.ogg'],\n     ['cibwar1/XC395511.ogg', 'cibwar1/XC432840.ogg'],\n     ['combul2/XC447668.ogg', 'combul2/XC650877.ogg']]","metadata":{"execution":{"iopub.status.busy":"2023-04-11T03:26:36.306995Z","iopub.execute_input":"2023-04-11T03:26:36.307665Z","iopub.status.idle":"2023-04-11T03:26:36.316292Z","shell.execute_reply.started":"2023-04-11T03:26:36.307628Z","shell.execute_reply":"2023-04-11T03:26:36.314685Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# make sure none of the known dupes are in our list of outliers\noutliers.filename.isin(set([x for li in known_dupes for x in li])).sum()","metadata":{"execution":{"iopub.status.busy":"2023-04-11T03:26:41.144879Z","iopub.execute_input":"2023-04-11T03:26:41.146264Z","iopub.status.idle":"2023-04-11T03:26:41.155610Z","shell.execute_reply.started":"2023-04-11T03:26:41.146196Z","shell.execute_reply":"2023-04-11T03:26:41.154535Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# find the euclid distance between the embeddings of all known duplicates\nfilename_to_index = {f: i for i, f in enumerate(df.filename)}\nknown_dupes = [(a, b, distances[filename_to_index[a], filename_to_index[b]].item()) for a, b in known_dupes]\nsorted(known_dupes, key=lambda a: a[2])","metadata":{"execution":{"iopub.status.busy":"2023-04-11T03:26:44.363888Z","iopub.execute_input":"2023-04-11T03:26:44.364505Z","iopub.status.idle":"2023-04-11T03:26:44.386414Z","shell.execute_reply.started":"2023-04-11T03:26:44.364444Z","shell.execute_reply":"2023-04-11T03:26:44.385501Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"n = len(combined_embeddings)\nnot_eye = (1 - torch.triu(torch.ones(n, n)).to(DEVICE)).bool()\n\ndef get_counts_for_threshold_range(start, end, n):\n    counts = []\n    thresholds = np.linspace(start, end, n)\n    for threshold in tqdm(thresholds):\n        counts.append(((distances < threshold) & not_eye).sum())\n    return thresholds, torch.stack(counts)\n\n# plot distribution of pairs under the max known dupe distance\nthresholds, counts = get_counts_for_threshold_range(0, max([x[2] for x in known_dupes]), 1000)\npx.line(\n    x=thresholds,\n    y=counts.to('cpu'),\n    title=\"Pairs closer than the MAX distance between known dupe pairs\",\n    labels={\"y\": \"number of recording pairs beneath threshold\", \"x\": \"embedding distnace threshold\"}\n)","metadata":{"execution":{"iopub.status.busy":"2023-04-11T03:26:48.366939Z","iopub.execute_input":"2023-04-11T03:26:48.367319Z","iopub.status.idle":"2023-04-11T03:27:05.440376Z","shell.execute_reply.started":"2023-04-11T03:26:48.367288Z","shell.execute_reply":"2023-04-11T03:27:05.439419Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# plot distribution of pairs under the min known dupe distance -- interesting that there are some!\nthresholds, counts = get_counts_for_threshold_range(0, min([x[2] for x in known_dupes]), 100)\npx.line(\n    x=thresholds,\n    y=counts.to('cpu'),\n    title=\"Pairs closer than the MIN distance between known dupe pairs\",\n    labels={\"y\": \"number of recording pairs beneath threshold\", \"x\": \"embedding distnace threshold\"}\n)","metadata":{"execution":{"iopub.status.busy":"2023-04-11T03:27:05.442295Z","iopub.execute_input":"2023-04-11T03:27:05.442771Z","iopub.status.idle":"2023-04-11T03:27:06.817532Z","shell.execute_reply.started":"2023-04-11T03:27:05.442734Z","shell.execute_reply":"2023-04-11T03:27:06.815410Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# helpers for displaying bird recordings for manual inspection\n\nfrom IPython.display import Audio\nimport torchaudio\nimport matplotlib.pyplot as plt\n\nSAMPLE_RATE = 32_000\n\ncompute_melspec = torchaudio.transforms.MelSpectrogram(\n    sample_rate=SAMPLE_RATE,\n    n_mels=128,\n    n_fft=2048, \n    hop_length=512,\n    f_min=0,\n    f_max=SAMPLE_RATE // 2,\n)\n\npower_to_db = torchaudio.transforms.AmplitudeToDB(\n    stype=\"power\",\n    top_db=80.0,\n)\n\n\ndef show_bird(index, start=0, secs=5):\n    audio = torchaudio.load(TRAIN_PATH / df.filename[index], start, start+32_000*secs)[0][0]\n    display(Audio(audio, rate=SAMPLE_RATE))\n    plt.figure(figsize=(12, 2.5))\n    plt.subplot(121)\n    plt.plot(audio)\n    plt.gca().get_xaxis().set_visible(False)\n    plt.subplot(122)\n    plt.imshow(power_to_db(compute_melspec(audio)))\n    plt.show()\n    return df.iloc[index]","metadata":{"execution":{"iopub.status.busy":"2023-04-11T03:27:18.001927Z","iopub.execute_input":"2023-04-11T03:27:18.002302Z","iopub.status.idle":"2023-04-11T03:27:18.537883Z","shell.execute_reply.started":"2023-04-11T03:27:18.002269Z","shell.execute_reply":"2023-04-11T03:27:18.536770Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# these are have low embedding distance to one another but are actually distinct\n# I found them by manually inspecting the spectrograms that were output from the loop below\nexclude = [\n    'comsan/XC554039.ogg',\n    'comsan/XC559524.ogg',\n    'comsan/XC582935.ogg',\n    'comsan/XC587730.ogg',\n    'comsan/XC646310.ogg',\n    'comsan/XC672721.ogg',\n    'comsan/XC580678.ogg',\n    'comsan/XC699920.ogg',\n    'wlwwar/XC638950.ogg',\n    'wlwwar/XC511683.ogg',\n    'wlwwar/XC568662.ogg',\n    'wlwwar/XC545476.ogg',\n    'barswa/XC651329.ogg',\n    'barswa/XC289357.ogg',\n    'comsan/XC659279.ogg',\n    'comsan/XC580678.ogg',\n    'quailf1/XC200244.ogg',\n    'eaywag1/XC653296.ogg',\n    'gargan/XC710595.ogg',\n    'eaywag1/XC636799.ogg',\n    'wlwwar/XC635230.ogg',\n    'wlwwar/XC372877.ogg',\n    'wlwwar/XC365486.ogg',\n    'wlwwar/XC213445.ogg',\n    'wlwwar/XC715016.ogg',\n    'wlwwar/XC511760.ogg',\n    'wlwwar/XC635230.ogg',\n    'wlwwar/XC581923.ogg',\n    'sltnig1/XC436062.ogg',\n    'sltnig1/XC436061.ogg',\n    'wlwwar/XC715016.ogg',\n    'wlwwar/XC635230.ogg',\n    'litegr/XC577784.ogg',\n    'litegr/XC576988.ogg',\n    'wlwwar/XC298820.ogg',\n    'wlwwar/XC113737.ogg',\n    'fotdro5/XC195989.ogg',\n    'gnbcam2/XC195528.ogg',\n]","metadata":{"execution":{"iopub.status.busy":"2023-04-11T03:27:19.981480Z","iopub.execute_input":"2023-04-11T03:27:19.982194Z","iopub.status.idle":"2023-04-11T03:27:19.988591Z","shell.execute_reply.started":"2023-04-11T03:27:19.982155Z","shell.execute_reply":"2023-04-11T03:27:19.987442Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"possible_dupes = ((distances < 5) & not_eye).nonzero().to('cpu').numpy()\npossible_dupes = [(a, b, distances[a][b].item()) for a, b in possible_dupes]\npossible_dupes = sorted(possible_dupes, key=lambda a: a[2])\n\noutlier_names = set(outliers.filename)\n\nembedding_dupes = []\n\nfor a, b, dist in possible_dupes:\n    a_filename = df.filename[a]\n    b_filename = df.filename[b]\n\n    if a_filename in outlier_names or b_filename in outlier_names:\n        # exclude noisy outliers\n        continue\n    if a_filename in exclude or b_filename in exclude:\n        # exclude any manually identified as distinct\n        continue\n\n    a_audio = torchaudio.load(TRAIN_PATH / a_filename)[0]\n    a_audio = torch.concat([a_audio[:, :5*32_000], a_audio[:, -5*32_000:]], -1)\n    b_audio = torchaudio.load(TRAIN_PATH / b_filename)[0]\n    b_audio = torch.concat([b_audio[:, :5*32_000], b_audio[:, -5*32_000:]], -1)\n\n    if dist > 1 and (a_audio.max() < 0.2 or b_audio.max() < 0.2):\n        # audio magnitude threshold for distances > 1\n        continue\n    if dist > 4 and (a_audio.max() < 0.5 or b_audio.max() < 0.5):\n        # audio magnitude threshold for distances > 4\n        continue\n\n    print(a, a_filename, b, b_filename, dist)\n    plt.figure(figsize=(12, 5))\n    plt.subplot(121)\n    plt.axis('off')\n    plt.imshow(power_to_db(compute_melspec(a_audio[0])))\n    plt.subplot(122)\n    plt.imshow(power_to_db(compute_melspec(b_audio[0])))\n    plt.axis('off')\n    plt.show()\n\n    embedding_dupes.append((b_filename, a_filename, dist))","metadata":{"execution":{"iopub.status.busy":"2023-04-11T03:27:20.791943Z","iopub.execute_input":"2023-04-11T03:27:20.792707Z","iopub.status.idle":"2023-04-11T03:31:08.422854Z","shell.execute_reply.started":"2023-04-11T03:27:20.792673Z","shell.execute_reply":"2023-04-11T03:31:08.420711Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"to_set = lambda dupes: set([','.join(sorted([str(s) for s in li])) for li in dupes])\nto_list = lambda dupes: [x.split(',')[1:] for x in dupes]\n\nknown_dupes_set = to_set(known_dupes)\nembedding_dupes_set = to_set(embedding_dupes)\n\nprint('known dupes not found from analyzing embeddings:')\ns = to_list(known_dupes_set - embedding_dupes_set)\nprint(len(s))\ndisplay(s)\n\nprint('\\nadditional dupes found from analyzing embeddings, not previously known:')\ns = to_list(embedding_dupes_set - known_dupes_set)\nprint(len(s))\ndisplay(s)\n\nprint('\\nall dupes, including embedding analysis:')\ns = to_list(embedding_dupes_set | known_dupes_set)\nprint(len(s))\ndisplay(s)","metadata":{"execution":{"iopub.status.busy":"2023-04-11T03:31:08.428678Z","iopub.execute_input":"2023-04-11T03:31:08.431107Z","iopub.status.idle":"2023-04-11T03:31:08.461908Z","shell.execute_reply.started":"2023-04-11T03:31:08.431066Z","shell.execute_reply":"2023-04-11T03:31:08.461021Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}