{"cells":[{"metadata":{},"cell_type":"markdown","source":"# Simple Sound Event Detection (SED) using Max Min filter ","execution_count":null},{"metadata":{},"cell_type":"markdown","source":"- I think the challenage of this competition is neither training the model nor managing the memory. They are just essential for training sound classification models using kaggle kernel. More important is how we handle the Sound Event Detection (SED) to filter out the background noise\n- if I understand it correctly, Hidehisa Arai applies an attention block in NN for SED (https://www.kaggle.com/hidehisaarai1213/introduction-to-sound-event-detection)\n- Alternatively, I use Max Min filter to detect the Sound event on Mel spectrogram. Then I relocate the signals on db spectrogram. \n- This might not be the best way to locate the signal, but just for sharing. Welcome for any comments. \n- In additional, I only apply normalization on freq axis of the spectrogram for each audio for the preprocess. My reasonale is that the birdsong signals are sigifiantly different from the background.\n- **neighborhood_size**, **threshold**, **clip_length**  and **sample** in **getSplitMelMaximaIndex()** are the 4 parameters need to test for different datasets\n- I am using Alex Shonenkov's birdcall dataset (https://www.kaggle.com/shonenkov/birdcall-check) for testing my SED method. The birdcall dataset contains 15 audio files. Most of the audio are less than 1min long, and only one audio is about 4 mins long. \n- Lets see how the length of the audio affect the SED performance.\n- Because the big update of TF from 1.0 to 2.0, I was so tired to learning TF/keras with confusing online documentation. I recently want to pick this up, so this is my first competition using TF/keras 2+ since the update. Not hope for geting a good LB score, but just for learning it.","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"# packages\n## system\nimport os, warnings, sys, pathlib\nfrom pathlib import Path\nwarnings.filterwarnings('ignore')\n\n## data structure\nimport pandas as pd\nimport numpy as np\nimport tensorflow as tf\nfrom tensorflow import keras\nimport tensorflow.keras.backend as K\nfrom tensorflow.python.keras.utils.data_utils import Sequence\n\n## utils\nimport math\nimport scipy\nimport scipy.ndimage as ndimage\nimport scipy.ndimage.filters as filters\nimport random\nimport datetime\nimport time\nfrom tqdm import tqdm\nfrom sklearn.pipeline import Pipeline\nfrom sklearn.utils import class_weight\nimport shutil\n\n## image \nimport cv2\nfrom PIL import Image\nfrom skimage.transform import rotate\nimport imgaug.augmenters as iaa\n\n## audio\nimport soundfile as sf\nimport librosa\nfrom pydub import AudioSegment\n\n## data cleaning\nfrom sklearn.preprocessing import StandardScaler, normalize, LabelEncoder, OneHotEncoder, scale\nfrom sklearn.decomposition import PCA\nfrom sklearn.metrics import confusion_matrix, f1_score, ConfusionMatrixDisplay\n\n## graphing\nimport librosa.display\nimport seaborn as sns\nimport matplotlib.pyplot as plt\nfrom matplotlib.collections import PatchCollection\nfrom matplotlib.patches import Rectangle\n%matplotlib inline","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# print(tf.test.is_gpu_available(), tf.config.list_physical_devices('GPU'))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# traning encoder","execution_count":null},{"metadata":{},"cell_type":"markdown","source":"## training data path","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"# get image path\n# define audio dataset path\n\nDATASET_AUDIO = \"train_audio/\"\nDATASET_SPEC = \"/kaggle/input/\" + \"birdsong-recognition/\" + DATASET_AUDIO\n\nprint(DATASET_SPEC)\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## lable encoder ","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"# get index of label and freq \nlabel_names = sorted(set(os.listdir(DATASET_SPEC)))\n\n# onehot encode label and freq\nlabel_to_onehot = OneHotEncoder(sparse = False)\nlabel_to_onehot.fit(np.array(label_names).reshape(-1, 1))\nlabel_dict = dict(zip(label_names, label_to_onehot.fit_transform(np.array(label_names).reshape(-1, 1))))\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"pd.DataFrame(label_dict)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## freq encoder ","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"freq_bin = np.arange(0,5000)\ndef split(word): \n    return [char for char in word] \nfreq_bin_encode = [split(str(freq).zfill(4)) for freq in freq_bin]\nfreq_bin_encode = [[int(digit) for digit in freq]  for freq in freq_bin_encode]\nfreq_dict = dict(zip(list(freq_bin), freq_bin_encode))\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"pd.DataFrame(freq_dict)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Model","execution_count":null},{"metadata":{"scrolled":true,"trusted":true},"cell_type":"code","source":"initializer = tf.keras.initializers.VarianceScaling()\nimages_shape = [128,128,1]\nfreqs_shape = len(freq_dict.get(0))\noutput_size = len(label_dict)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"model = tf.keras.models.load_model(\"/kaggle/input/model-v2-e45/model_v2_e45.h5\")","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# for layer in model.layers:\n#     print(layer.output_shape)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# audio Loader","execution_count":null},{"metadata":{},"cell_type":"markdown","source":"## readWave","execution_count":null},{"metadata":{"scrolled":true,"trusted":true},"cell_type":"code","source":"# load the audio as a wave, Store the sampling rate as `sr`\ndef readWave(path, displayWave=False, windowSize = 200):\n#     def decode_mp3(path, normalized=False):\n#         \"\"\"MP3 to numpy array\"\"\"\n#         AS = AudioSegment.from_mp3(path)\n#         wave = np.array(AS.get_array_of_samples()).astype(float)\n#         sr = AS.frame_rate\n#         if AS.channels ==2:\n#             sr = sr * 2\n#         return wave, sr\n\n    def decode_mp3(path):\n        \"\"\"MP3 to numpy array\"\"\"\n        wave, sr = librosa.core.load(path, sr = 44100)\n        return wave, sr\n    \n    \n    wave, sr = decode_mp3(path)\n        \n#     add_time = int(windowSize * hop_length /sr /2) # about 6 second for windowSize = 2**9\n\n    if displayWave:\n        plt.figure(figsize=(10, 4))\n        librosa.display.waveplot(y = wave, sr=sr)\n        plt.title(path +\", sr=\" + str(sr) + \" (length: {})\".format(len(wave)))\n        plt.show()\n    return wave, sr\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## wave2spectrogram","execution_count":null},{"metadata":{"scrolled":false,"trusted":true},"cell_type":"code","source":"# waveform -> spectrogram\ndef wave2spectrogram(wave, display = False, specType = \"db\"):\n    def denoise(spectrogram):\n        # prepare for extrapolation\n        spectrogram[spectrogram == 0] = 'nan'\n        # remove background white noise\n        spectrogram = scale(spectrogram, axis=1)\n        # repalce nan with global mean\n        spectrogram[np.isnan(spectrogram)] = np.nanmean(spectrogram).mean()\n        # amplify the signal # not so good\n#         spectrogram = spectrogram**(2)\n        return spectrogram\n    \n    def remap(x, out_min, out_max):\n        return (x - x.min()) * (out_max - out_min) / (x.max() - x.min()) + out_min\n\n    hop_length = 512 # number of samples between each successive FFT window\n    n_fft = 2048 # number of samples in a single window\n    \n    # db\n    stft_wave = librosa.stft(wave, n_fft=n_fft, hop_length=hop_length)\n    spectrogram = librosa.amplitude_to_db(abs(stft_wave))\n    spectrogram = remap(spectrogram, 0, 255)\n    spectrogram = denoise(spectrogram) \n    spectrogram = remap(spectrogram, 0, 255)\n    \n    if display and specType == \"db\":\n        plt.figure(figsize=(10, 4))\n        librosa.display.specshow(spectrogram, sr=sr, y_axis='hz', x_axis='time')\n        plt.title(\"db spectrogram: {}\".format(spectrogram.shape))\n        plt.colorbar()\n        plt.show()\n    # melspectrogram  \n    if specType == \"mel\":\n        spectrogram = librosa.feature.melspectrogram(y=wave, n_fft=n_fft, hop_length=hop_length)\n        spectrogram = remap(spectrogram, 0, 255)\n    if display and specType == \"mel\":\n        plt.figure(figsize=(10, 4))\n        librosa.display.specshow(spectrogram, sr=sr, y_axis='mel', x_axis='time')\n        plt.title(\"Mel spectrogram: {})\".format(spectrogram.shape))\n        plt.colorbar()\n        plt.show()\n        \n    return spectrogram\n\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def meanboundary(spectrogram, sr = 44100, hop_length = 512, time = 0.2, display = False):\n    height = spectrogram.shape[0]\n    width = spectrogram.shape[1]\n    width_remove = round(time * sr / hop_length)\n    # left and right\n    replacewith = spectrogram.mean()\n    spectrogram[:, 0: width_remove] = replacewith\n    spectrogram[:, width - width_remove: width] = replacewith\n    spectrogram[0: 3, :] = replacewith\n    spectrogram[height - 3: height, :] = replacewith\n    \n    if display:\n        plt.figure(figsize=(10, 4))\n        sns.heatmap(spectrogram)\n        plt.title(\"spectrogram: {})\".format(spectrogram.shape))\n        plt.gca().invert_yaxis()\n        plt.show()\n    return spectrogram\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## spectrogram2image (Sound Event Detection by Max, Min filters)","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"def remap(x, out_min, out_max):\n    return (x - x.min()) * (out_max - out_min) / (x.max() - x.min()) + out_min\n\ndef findLocalMax(spectrogram, neighborhood_size, threshold):\n    spectrogram_min = filters.minimum_filter(spectrogram, neighborhood_size)\n    spectrogram_max = filters.maximum_filter(spectrogram, neighborhood_size)\n    diff_spectrogram = (spectrogram_max - spectrogram_min)\n    diff_spectrogram = remap(diff_spectrogram, 0, 1)\n    diff = ((diff_spectrogram) > threshold)\n    maxima = np.logical_and(spectrogram == spectrogram_max, diff)\n    \n    \n    plt.figure(figsize=(10, 4))\n    sns.heatmap(diff_spectrogram)\n    plt.gca().invert_yaxis()\n    plt.show()\n    plt.figure(figsize=(10, 4))\n    sns.heatmap(diff)\n    plt.gca().invert_yaxis()\n    plt.show()\n    return maxima\n","execution_count":null,"outputs":[]},{"metadata":{"scrolled":false,"trusted":true},"cell_type":"code","source":"def getMelMaximaIndex(mel_spectrogram, sr, clip_length = 60, displaySpectralPoint = False):\n    def getSplitMelMaximaIndex(spectrogram, neighborhood_size = 32, threshold = 0.4, sample = 30):\n        \"\"\"\n        modified from \n        https://stackoverflow.com/questions/9111711/get-coordinates-of-local-maxima-in-2d-array-above-certain-value\n        \"\"\"\n        def denoise(spectrogram):\n            # prepare for extrapolation\n            spectrogram[spectrogram < 0] = 'nan'\n            # remove background white noise\n            spectrogram = scale(spectrogram, axis=1)\n            # repalce nan with global mean\n            spectrogram[np.isnan(spectrogram)] = np.nanmean(spectrogram).mean()\n            # amplify the signal # not so good\n            #spectrogram = spectrogram**2\n            return spectrogram\n\n        def remap(x, out_min, out_max):\n            return (x - x.min()) * (out_max - out_min) / (x.max() - x.min()) + out_min\n        \n        def findLocalMax(spectrogram, neighborhood_size, threshold):\n            spectrogram_min = filters.minimum_filter(spectrogram, neighborhood_size)\n            spectrogram_max = filters.maximum_filter(spectrogram, neighborhood_size)\n            diff_spectrogram = (spectrogram_max - spectrogram_min)\n            diff_spectrogram = remap(diff_spectrogram, 0, 1)\n            diff = ((diff_spectrogram) > threshold)\n            maxima = np.logical_and(spectrogram == spectrogram_max, diff)\n            return maxima\n        \n        maxima_reg = findLocalMax(spectrogram = spectrogram, \n                                  neighborhood_size = neighborhood_size, \n                                  threshold = threshold)\n        \n        spectrogram = remap(spectrogram, -1, 1)\n        spectrogram = denoise(spectrogram) \n        spectrogram = remap(spectrogram, 0, 1)\n        maxima_denoise = findLocalMax(spectrogram = spectrogram, \n                                  neighborhood_size = neighborhood_size, \n                                  threshold = threshold)\n\n        maxima = np.logical_or(maxima_reg, maxima_denoise)\n        \n\n        labeled, num_objects = ndimage.label(maxima)\n        xy = np.array(ndimage.center_of_mass(maxima, labeled, range(1, num_objects+1)))\n        \n#         print(\"num of valid maxima_reg: \", maxima_reg.sum())\n#         print(\"num of valid maxima_denoise: \", maxima_denoise.sum())\n#         print(\"num of valid maxima: \", maxima.sum())\n#         print(np.argwhere(labeled != 0))\n        x = []\n        y = []\n        if len(xy.shape)>1:\n            xy = xy[~np.isnan(xy).any(axis=1)]\n            y = list(np.around(xy[:,0]).astype(int))\n            x = list(np.around(xy[:,1]).astype(int))\n            # limite max num of samples\n            if len(x) > sample:\n                sampled_indice = random.sample(list(range(len(x))), sample)\n                x = [ x[i] for i in sampled_indice ]\n                y = [ y[i] for i in sampled_indice ]\n        return x, y\n    \"\"\"\n    modified from \n    https://stackoverflow.com/questions/37999150/how-to-split-a-wav-file-into-multiple-wav-files\n    \"\"\"\n    # define croping area\n    hop_length = 512\n    totalIndex = mel_spectrogram.shape[1]\n    # scan for every 1 mins as defaul \n    min_fft_per_split = int(clip_length * sr / hop_length)\n    max_num_split = totalIndex // min_fft_per_split\n    if max_num_split ==0: \n        # in case max_num_split ==0\n        xs, ys = getSplitMelMaximaIndex(spectrogram = mel_spectrogram)\n    else:\n        min_fft_per_split = math.ceil(totalIndex / max_num_split)\n        # overLap with 1/4 * mins_per_split\n        split_overLap = math.ceil(min_fft_per_split*3/4)\n        # split spectrogram in to piece and get local mel xy\n        xs = []\n        ys = []\n        for i in range(0, totalIndex, split_overLap):\n            # output for each clip\n            if i + min_fft_per_split <= totalIndex:\n                index_clip = np.arange(totalIndex)[i: i + min_fft_per_split]\n            if i + min_fft_per_split > totalIndex:\n                index_clip = np.arange(totalIndex)[i: totalIndex]\n            # calculate x and y at mel_spectrogram \n            x, y = getSplitMelMaximaIndex(spectrogram = mel_spectrogram[:,index_clip])\n            \n            x = list(np.array(x, dtype = type(x) ) + i)\n            xs += x\n            ys += y\n        \n    if displaySpectralPoint and len(ys)>1:\n        plt.figure(figsize=(10, 4))\n        sns.heatmap(mel_spectrogram)\n        plt.plot(xs, ys, 'ro', mfc='none', color='white')\n        plt.title(\"mel spectrogram {} (split: {}), (points: {})\".format(mel_spectrogram.shape, max_num_split, len(ys)))\n        plt.gca().invert_yaxis()\n        plt.show()\n    return xs, ys\n\n\n\n# there is some issue, but so clear. find a way to combine reg and denoised maximum ","execution_count":null,"outputs":[]},{"metadata":{"scrolled":true,"trusted":true},"cell_type":"code","source":"def yMelIndex2yhz(y_mel_index, sr, n_fft=2048):\n    hz_scale = librosa.core.fft_frequencies(sr=sr, n_fft=n_fft)\n    mel_scale = librosa.core.mel_frequencies(n_mels=128, fmin=hz_scale.min(), \n                                             fmax=hz_scale.max(), htk=False)\n    y_hz = mel_scale[int(y_mel_index)] \n    return y_hz/10\n","execution_count":null,"outputs":[]},{"metadata":{"scrolled":false,"trusted":true},"cell_type":"code","source":"def ymel2hzindex(y_mel_index, sr, n_fft = 2048, printSummary = False):\n    def find_nearest(a, a0):\n        \"\"\"\n        https://stackoverflow.com/questions/2566412/find-nearest-value-in-numpy-array\n        \"\"\"\n        \"Element in nd array `a` closest to the scalar value `a0`\"\n        idx = np.abs(a - a0).argmin()\n        return a.flat[idx]\n    \"\"\"\n    http://man.hubwiz.com/docset/LibROSA.docset/Contents/Resources/Documents/_modules/librosa/display.html#specshow\n    \"\"\"\n    hz_scale = librosa.core.fft_frequencies(sr=sr, n_fft=n_fft)\n    mel_scale = librosa.core.mel_frequencies(n_mels=128, \n                                               fmin=hz_scale.min(), fmax=hz_scale.max(), htk=False)\n\n    y_hz = mel_scale[int(y_mel_index)] \n    y_hz_nearest = find_nearest(hz_scale, y_hz)\n    y_hz_index = list(hz_scale).index(y_hz_nearest)\n    if printSummary:\n        print(\"y_mel_index: {}, \".format(y_mel_index), \n              \"y_hz: {}\".format(y_hz), \n              \"y_hz_index: {}, \".format(y_hz_index), \n              \"(y_hz_nearest: {}) \".format(y_hz_nearest))\n    return y_hz_index\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def cropResize(spectrogram, xdata, ydata, \n               window_height_haft = 200, window_width_haft = 200,  \n               resizeShape = (128,128),\n               displayCropedImage = False):\n    \n    # define cropping window\n    h_spec = spectrogram.shape[0]\n    w_spec = spectrogram.shape[1]\n\n    w_l = xdata - window_width_haft\n    w_h = xdata + window_width_haft\n    h_l = ydata - window_height_haft\n    h_h = ydata + window_height_haft\n    w_l_offset = - w_l if w_l <= 0 else 0\n    w_h_offset = window_height_haft*2 - (w_h - w_spec) if w_h >= w_spec else window_height_haft*2\n    h_l_offset = - h_l if h_l <= 0 else 0\n    h_h_offset = window_height_haft*2 - (h_h - h_spec) if h_h >= h_spec else window_height_haft*2\n    \n    w_l = w_l if w_l >= 0 else 0\n    w_h = w_h if w_h <= w_spec else w_spec\n    h_l = h_l if h_l >= 0 else 0\n    h_h = h_h if h_h <= h_spec else h_spec\n\n    # crop image\n    image = spectrogram[ h_l: h_h, w_l: w_h ]   \n    img = np.full((window_height_haft*2, window_width_haft*2), image.mean() ) # create a single channel  image  \n    img[h_l_offset: h_h_offset, w_l_offset: w_h_offset ] = image[:,:]\n\n    # resize\n    image = cv2.resize(img, resizeShape, cv2.INTER_AREA)  \n    \n    def remap(x, out_min, out_max):\n        return (x - x.min()) * (out_max - out_min) / (x.max() - x.min()) + out_min\n    image = remap(image, 0, 255)\n    \n    if displayCropedImage:\n        plt.figure(figsize=(4, 4))\n        sns.heatmap(img, cbar = False)\n        plt.title(\"croped image {}  at ({}, {})\".format(image.shape, xdata, ydata))\n        plt.gca().invert_yaxis()\n        plt.show()\n    \n    return np.expand_dims(image, -1)\n","execution_count":null,"outputs":[]},{"metadata":{"scrolled":true,"trusted":true},"cell_type":"code","source":"def plotSpecWin(x, y, spectrogram, window_height_haft = 200, window_width_haft = 200):\n    def makeBoxes(ax, x, y, window_width_haft, window_height_haft,\n                  linewidth=1, edgecolor='white', facecolor='none'):\n        \"\"\"\n        modify from\n        https://matplotlib.org/api/_as_gen/matplotlib.patches.Rectangle.html#matplotlib.patches.Rectangle\n        \"\"\"\n        # define cropping window\n        h_spec = spectrogram.shape[0]\n        w_spec = spectrogram.shape[1]\n        \n        w_l = [i - window_width_haft if i - window_width_haft >= 0 else 0 for i in x]\n        w_h = [i + window_width_haft if i + window_width_haft <= w_spec else w_spec for i in x]\n        h_l = [i - window_height_haft if i - window_height_haft >= 0 else 0 for i in y]\n        h_h = [i + window_height_haft if i + window_height_haft <= h_spec else h_spec for i in y]\n\n        boxes = zip(w_l,\n                    h_l,\n                    list(np.array(w_h)- np.array(w_l)), \n                    list(np.array(h_h)- np.array(h_l)))\n        # Loop over data points; create box from errors at each point\n        errorboxes = [Rectangle((x, y), xe, ye) for x, y, xe, ye in boxes]\n        # Create patch collection with specified colour/alpha\n        pc = PatchCollection(errorboxes,\n                             linewidth=linewidth, \n                             edgecolor = edgecolor, \n                             facecolor = facecolor)\n        return pc\n    fig = plt.figure(figsize=(10, 4))\n    ax = fig.subplots()\n    sns.heatmap(spectrogram)\n    pc = makeBoxes(ax, x, y, window_width_haft, window_height_haft)\n    ax.add_collection(pc) # # Add collection to axes\n    plt.title(\"spectrogram {} with cropping windows\".format(spectrogram.shape))\n    fig.gca().invert_yaxis()\n    fig.show()\n","execution_count":null,"outputs":[]},{"metadata":{"scrolled":false,"trusted":true},"cell_type":"code","source":"def xIndex2xTime(x_index, sr = 44100):\n    hop_length = 512\n    x_time = x_index / sr * hop_length\n    x_time = round(x_time, 1)\n    \n    return x_time\n","execution_count":null,"outputs":[]},{"metadata":{"scrolled":false,"trusted":true},"cell_type":"code","source":"# spectrogram -> images\ndef spectrogram2image(db_spectrogram, mel_spectrogram, sr, clip_length= 60,\n                      window_height_haft = 200, window_width_haft = 200, \n                      displaySpectralPoint = False, displayPointCoord = False, \n                      displayCropedImage = False, displayCropedSpec = False):\n    # get the coordinatate center of proposed windows\n    x,y = getMelMaximaIndex(mel_spectrogram = mel_spectrogram, sr = sr,clip_length=clip_length, displaySpectralPoint = displaySpectralPoint)  \n    # find corresponding y_hzs\n    y_hz = map(lambda y: yMelIndex2yhz(y_mel_index = y, sr = sr), y)\n    # find the time in the original audio\n    x_time = map(lambda x: xIndex2xTime(x_index = x, sr = sr), x)\n    # cast the mel scale coordination back to hz scale\n    y = list(map(lambda y: ymel2hzindex(y_mel_index = y, sr = sr), y))\n    # crop and resized images \n    image_db = map(lambda x, y: cropResize(xdata = x, ydata = y,\n                                         spectrogram = db_spectrogram, \n                                         window_height_haft = window_height_haft, \n                                         window_width_haft = window_width_haft,\n                                         displayCropedImage = displayCropedImage), x, y)  \n    if displayCropedSpec:\n        plotSpecWin(x = x, y = y, spectrogram = db_spectrogram, \n                    window_width_haft = window_width_haft, window_height_haft = window_height_haft)      \n    return list(image_db), len(x), list(y_hz), list(x_time)\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## composited loader","execution_count":null},{"metadata":{"scrolled":true,"trusted":true},"cell_type":"code","source":"def dataLoader(filename, clip_length = 60,\n               window_height = 200, window_width = 200,\n               displayWave = False, displaySpectrogram = False, \n               displaySpectralPoint = False, displayPointCoord = False,\n               displayCropedImage = False, displayCropedSpec = False):\n    window_height_haft = int(window_height/2)\n    window_width_haft = int(window_width/2)\n    # read audio file\n    wave, sr = readWave(filename, displayWave = displayWave, windowSize = window_height)\n    # generate spectrogram\n    db_spectrogram = wave2spectrogram(wave, display = displaySpectrogram, specType = \"db\")\n    mel_spectrogram = wave2spectrogram(wave, display = displaySpectrogram, specType = \"mel\")\n    db_spectrogram = meanboundary(db_spectrogram, sr = sr)\n    mel_spectrogram = meanboundary(mel_spectrogram, sr = sr)  \n    # generate images\n    images, num_images, y_hz, x_time = spectrogram2image(db_spectrogram = db_spectrogram,\n                                           mel_spectrogram = mel_spectrogram,\n                                           sr = sr, \n                                           clip_length = clip_length,\n                                           window_height_haft = window_height_haft, \n                                           window_width_haft = window_width_haft,\n                                           displaySpectralPoint = displaySpectralPoint,\n                                           displayPointCoord = displayPointCoord,\n                                           displayCropedImage = displayCropedImage,\n                                           displayCropedSpec = displayCropedSpec)  \n    return images, num_images, y_hz, x_time","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Data generator ","execution_count":null},{"metadata":{},"cell_type":"markdown","source":"## load audio and produce images","execution_count":null},{"metadata":{"scrolled":false,"trusted":true},"cell_type":"code","source":"def paths_to_image(audio_paths, i, clip_length = 60, verbose = 0):\n    \"\"\"\n    modified from\n    https://keras.io/examples/audio/speaker_recognition_using_cnn/\n    Constructs a dataset of audios and labels.\n    \"\"\"        \n    image_ds = []\n    freqs_ds = []\n    times_ds =[]\n    image_audio_paths = []\n    len_audio_path = len(audio_paths)\n    for num_audio_path in np.arange(len_audio_path):\n        images, num_images, freqs, times = dataLoader(audio_paths[num_audio_path], clip_length = clip_length) # crop\n        image_ds += images\n        freqs_ds += freqs\n        times_ds += times\n        image_audio_paths += [audio_paths[num_audio_path]] * num_images\n        if verbose == 1:\n            print(\"processing {}\".format(audio_paths[num_audio_path]))\n            print(\"batch_{} ({}/{}):  +{} , total sample: {}\".format(\n                i, num_audio_path+1, len_audio_path, len(images), len(times_ds)))\n        if verbose == 0:\n            print(\"batch_{} ({}/{}): processing {}\".format(\n                i, num_audio_path+1, len_audio_path, audio_paths[num_audio_path]))\n    print(\"Total sample: {}\".format(len(times_ds)))\n    return image_audio_paths, image_ds, freqs_ds, times_ds\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## save crop spec image","execution_count":null},{"metadata":{"scrolled":true,"trusted":true},"cell_type":"code","source":"# loop over the label and save files \ndef saveImage(image_audio_paths, image_ds, freqs_ds, times_ds, save_dir):\n    \"\"\"\n    modified from\n    https://stackoverflow.com/questions/2659312/how-do-i-convert-a-numpy-array-to-and-display-an-image\n    \"\"\"\n    # create folder for input_data\n    if not os.path.exists(save_dir):\n        os.mkdir(save_dir)\n    len_sample = len(freqs_ds)\n    times_ds = [str(times).replace('.', '') for times in times_ds]\n\n    for num_smaple in np.arange(len_sample):\n        # create label folder if not exist\n        if not os.path.exists(save_dir):\n            os.mkdir(save_dir)\n        # write png file\n        img = image_ds[num_smaple][:,:,0]\n        img = img.astype(np.uint8)\n        img = Image.fromarray(img, 'L')\n        freqs_ds = list(map(int, freqs_ds))\n        audio_filename = image_audio_paths[num_smaple].split('/')[-1].split('.')[0]\n        img.save(save_dir + \"/{}_{}_{}.jpg\".format(audio_filename, freqs_ds[num_smaple], times_ds[num_smaple]))\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Locate the test audio path","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"test_root = \"/kaggle/input/birdsong-recognition/\"\ntest_dir = test_root + 'test_audio/'\nprint(test_dir, test_dir)\n\nif os.path.isdir(test_dir):\n    print('Running actual submission!')\n    test_root = test_root\n    test_dir = test_root + \"test_audio/\"\n    test_df = test_root + 'test.csv' \n    test_df = pd.read_csv(test_df)\n    print(\"test_root: \", test_root)\n    print(\"test_dir: \", test_dir)\nelse:\n#     test_folder = 'faketest/'\n    test_folder = 'birdcall-check/'\n    print('Running debug submission')\n    test_root = '/kaggle/input/'  + test_folder\n    test_dir = test_root + \"test_audio\"\n    test_df = test_root + 'test.csv'\n    test_df = pd.read_csv(test_df)\n    print('Naturally no test dir found') \n    print(\"test_root: \", test_root)\n    print(\"test_dir: \", test_dir)\nsub = pd.read_csv(\"/kaggle/input/birdsong-recognition/sample_submission.csv\")\nsub.to_csv(\"submission.csv\", index=False)  # this will be overwritten if everything goes well","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"DATASET_AUDIO_PATH = test_dir\ndir_path = Path(DATASET_AUDIO_PATH)\naudio_paths = [\n    os.path.join(dir_path, filepath)\n    for filepath in os.listdir(dir_path)\n    if filepath.endswith(\".mp3\")\n]\nlen_audio_paths = len(audio_paths)\nprint(\"count of audio: \", len_audio_paths)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# Visualize Spectrogram\n\nlen_audio_paths = len(audio_paths)\ni=0\nfor num_audio_path in np.arange(len_audio_paths):\n    i+=1\n    print(\"{}/{} processing:{} \".format(i,len_audio_paths, audio_paths[num_audio_path]))\n    images, num_images, freqs, times = dataLoader(audio_paths[num_audio_path], \n                                                  clip_length= 30, \n                                                  displaySpectralPoint = True)\n    \n    ","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"The 4min audio is much longer then the other audio in the birdcall dataset, so the current set of SED parameters is not fit for this long audio, but it seem good for the the short audios. For the private test datatset, the audio are all about 10 minutes long, so there is probably one set of SED parameters for the private test.","execution_count":null},{"metadata":{},"cell_type":"markdown","source":"# Composited data generator","execution_count":null},{"metadata":{"scrolled":false,"trusted":true},"cell_type":"code","source":"num_process = 2\nbatch_size = 5\nnum_batch = math.ceil(len(audio_paths) / batch_size)\nprint(\"num_process: \", num_process)\nprint(\"num_batch\", num_batch)\nprint(\"batch_size: \", batch_size)\naudio_paths_batches = [audio_paths[x: x + batch_size] for x in range(0, len(audio_paths), batch_size)]\nprint(\"total number of audio file: \", len(audio_paths))\nsave_dir = \"/kaggle/working/img_test/\"\nprint(\"save to: \", save_dir)","execution_count":null,"outputs":[]},{"metadata":{"scrolled":true,"trusted":true},"cell_type":"code","source":"def loadSave(audio_paths_batches, i, save_dir):\n    batch_size = len(audio_paths_batches)\n    print(\"processing batch: {} for audio #{} to #{}\".format(i, 0, batch_size ))\n    image_audio_paths, image_ds, freqs_ds, times_ds = paths_to_image(audio_paths_batches,\n                                                                     clip_length = 30,\n                                                                     i=i,\n                                                                     verbose = 0)\n    print(\"batch_{} is processed\".format(i))\n    time.sleep(0.2)\n    \n    print(\"saving to batch_{} to [{}]\".format(i,save_dir))\n    saveImage(image_audio_paths = image_audio_paths, \n                     image_ds = image_ds, \n                     freqs_ds = freqs_ds,\n                     times_ds = times_ds, \n                     save_dir = save_dir)\n    print(\"saving batch_{} successfully\".format(i))\n    time.sleep(0.2)\n    del image_ds, freqs_ds, times_ds\n    \n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"for i in range(num_batch):\n    loadSave(audio_paths_batches[i], i, save_dir)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Prediction\n- the audio will be evaluate in a a clip with  3/4 overlap for high density detection. \n- make prediction on all sampling using local maximum window with a fixed threhold to target on the high power sound,\n- then assign the prediction to each signal ","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"# load dataset and preprocess images（the preprocess function in map  can preprocess all previous pictures)\ndef load_and_preprocess_image(img_path):\n    # read pictures\n    img_raw = tf.io.read_file(img_path)\n    # decode pictures\n    img_tensor = tf.image.decode_jpeg(img_raw, channels=1) # change channels to 3 result in 3-channel image\n    # resize\n#     img_tensor = tf.image.resize(img_tensor, [128, 128])\n    #tf.cast() function is a type conversion function that converts the data format of x into dtype\n    img_tensor = tf.cast(img_tensor, tf.float32)\n    # normalization\n    img_tensor = img_tensor / 255.0\n    # flip left or right\n#     img_tensor = tf.image.random_flip_left_right(img_tensor)\n    return img_tensor\n","execution_count":null,"outputs":[]},{"metadata":{"scrolled":false,"trusted":true},"cell_type":"code","source":"# define valid path\ntest_image_root = pathlib.Path(save_dir)\ntest_image_path = test_image_root.rglob('*.jpg')\ntest_image_path = [str(pathlib.Path(path)) for path in test_image_path]\nprint(\"total test images: \", len(test_image_path))\n\n# get index of label and freq \ntest_freq = [int(path.split('/')[-1].split('_')[1]) for path in test_image_path]\nencode_test_freq = [freq_dict.get(freq) for freq in test_freq]\n\n# convert to tensor and zip\ntest_image_ds = tf.data.Dataset.from_tensor_slices(test_image_path).map(load_and_preprocess_image)\ntest_freq_ds = tf.data.Dataset.from_tensor_slices(encode_test_freq)\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"batch_size = 128\ntest_ds = tf.data.Dataset.zip(((test_image_ds, test_freq_ds), )).batch(batch_size)\niteration_test = iter(test_ds)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"scrolled":true},"cell_type":"code","source":"audio_id = [path.split('/')[-1].split('_')[0] for path in test_image_path]\nraw_seconds = [int(path.split('/')[-1].split('_')[2].split('.')[0])/10 for path in test_image_path]\nseconds = [5 * math.ceil(s/5) for s in raw_seconds]\n\nresult_pred_class = []\nresult_pred_class_prob = []\nfor i in range(0, len(test_image_path), batch_size):\n    pred_prob = model.predict(iteration_test.get_next())\n    pred_class_prob = np.amax(pred_prob, 1)\n    pred_class_integer = np.argmax(pred_prob, 1)\n    pred_class = [list(label_dict)[i] for i in pred_class_integer]\n    result_pred_class += pred_class\n    result_pred_class_prob += list(pred_class_prob)\n    \nresult_raw = pd.DataFrame(data= {'audio_id': audio_id,\n                             'raw_seconds': raw_seconds, \n                             'seconds': seconds, \n                             'freq': test_freq, \n                             'birds': result_pred_class,\n                             'pred_class_prob': result_pred_class_prob}).sort_values(by = ['audio_id', 'raw_seconds'])","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"threshold = 0.8\nresult = result_raw[result_raw['pred_class_prob'] > threshold]\nsns.distplot(result['pred_class_prob'])\nresult","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"##  Site 1 & 2\n- each signal will be aggregate to the 5 seconds time frame,\n- the final prediction are determined by the weighted probablity of each bird prediction within each 5s clip. hopefully this will combine the effect of count of signal and the confident of each signal  \n","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"def appandLabel(x):\n    x = list(x)\n    separator = ' '\n    x = separator.join(x)\n    return x","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sub_test_12 = test_df[test_df.site.isin(['site_1', 'site_2'])]\nsub_test_12 = sub_test_12.merge(result, on = ['audio_id','seconds'], how ='left')\nsub_test_12 = sub_test_12[['row_id','birds','pred_class_prob']].groupby(['row_id','birds'])['pred_class_prob'].sum().reset_index()\nsub_test_12 = sub_test_12[sub_test_12['pred_class_prob'] > 1]\nsub_test_12 = sub_test_12[['row_id','birds']].drop_duplicates().groupby(['row_id'])['birds'].apply(appandLabel).reset_index()\nsub_test_12.head(25)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## site 3\n- sum pred_class_prob by group of 'row_id','birds'\n- remove the sum prob less than 15","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"sub_test_3 = test_df[test_df.site.isin(['site_3'])]\nsub_test_3 = sub_test_3[['row_id', 'audio_id']].merge(result[['audio_id','birds','pred_class_prob']], on=['audio_id'], how='left')\nsub_test_3 = sub_test_3[['row_id','birds','pred_class_prob']].groupby(['row_id','birds'])['pred_class_prob'].sum().reset_index()\nsub_test_3 = sub_test_3[sub_test_3['pred_class_prob'] > 15]\nsub_test_3 = sub_test_3[['row_id','birds']].drop_duplicates().groupby(['row_id'])['birds'].apply(appandLabel).reset_index()\nsub_test_3.head(25)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## submission","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"try:\n    shutil.rmtree(save_dir)\nexcept:\n    pass","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sub = sub_test_12.append(sub_test_3)\nsub = test_df.merge(sub, how='left')[['row_id','birds']]\nsub['birds'] = sub['birds'].fillna('nocall')\nsub.to_csv('submission.csv', index=False)\nsub.iloc[random.sample(range(test_df.shape[0]),50),:]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"plt.figure(figsize=(15, 5))\nsns.countplot(x=\"birds\", data = sub)\nplt.title(\"Distribution of the birdsong prediction\")\nplt.xticks(rotation=45)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"","execution_count":null,"outputs":[]}],"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":4,"nbformat_minor":4}