{"cells":[{"metadata":{},"cell_type":"markdown","source":"![](https://lh3.googleusercontent.com/-ChT8E8LooOg/X7V8Mi8Uv_I/AAAAAAAAHZY/c_kzqGN4Z686UBICVsec5SkCBfNycxE_gCLcBGAsYHQ/header.png)"},{"metadata":{},"cell_type":"markdown","source":"### We will mainly use two libraries for audio acquisition and playback:\n# Librosa\nIt is a Python module to analyze audio signals in general but geared more towards music. It includes the nuts and bolts to build a MIR(Music information retrieval) system. It has been very well [documented](https://librosa.org/librosa/) along with a lot of examples and tutorials.\n# IPython.display.Audio\nIPython.display.Audio lets you play audio directly in a jupyter notebook."},{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"cell_type":"code","source":"!pip install librosa","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","trusted":true},"cell_type":"code","source":"import numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\nimport os","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Loading an audio file:"},{"metadata":{"trusted":true},"cell_type":"code","source":"import librosa\naudio_data = '../input/rfcx-species-audio-detection/train/00204008d.flac'\nx , sr = librosa.load(audio_data)\nprint(type(x), type(sr))\nprint(x.shape, sr)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"#### This returns an audio time series as a numpy array with a default sampling rate(sr) of 22KHZ mono. We can change this behavior by resampling at 44.1KHz."},{"metadata":{"trusted":true},"cell_type":"code","source":"librosa.load(audio_data, sr=44100)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Playing Audio:\n\n#### Using,IPython.display.Audio you can play the audio in your jupyter notebook."},{"metadata":{"trusted":true},"cell_type":"code","source":"import IPython.display as ipd\nipd.Audio(audio_data)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"%matplotlib inline\nimport matplotlib.pyplot as plt\nimport librosa.display\nplt.figure(figsize=(14, 5))\nlibrosa.display.waveplot(x, sr=sr)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### Spectrogram\n#### A spectrogram is a visual way of representing the signal strength, or “loudness”, of a signal over time at various frequencies present in a particular waveform. Not only can one see whether there is more or less energy at, for example, 2 Hz vs 10 Hz, but one can also see how energy levels vary over time.\n\n#### A spectrogram is usually depicted as a [heat map](https://en.wikipedia.org/wiki/Heat_map), i.e., as an image with the intensity shown by varying the color or brightness.\n\n### We can display a spectrogram using. librosa.display.specshow."},{"metadata":{"trusted":true},"cell_type":"code","source":"X = librosa.stft(x)\nXdb = librosa.amplitude_to_db(abs(X))\nplt.figure(figsize=(14, 5))\nlibrosa.display.specshow(Xdb, sr=sr, x_axis='time', y_axis='hz')\nplt.colorbar()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"#### .stft() converts data into short term Fourier transform. [STFT](https://www.youtube.com/watch?v=g1_wcbGUcDY) converts signals such that we can know the amplitude of the given frequency at a given time. Using STFT we can determine the amplitude of various frequencies playing at a given time of an audio signal. .specshow is used to display a spectrogram.\n\n#### The vertical axis shows frequencies (from 0 to 10kHz), and the horizontal axis shows the time of the clip. Since we see that all action is taking place at the bottom of the spectrum, we can convert the frequency axis to a logarithmic one."},{"metadata":{"trusted":true},"cell_type":"code","source":"librosa.display.specshow(Xdb, sr=sr, x_axis='time', y_axis='log')\nplt.colorbar()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Create an Audio Signal:"},{"metadata":{"trusted":true},"cell_type":"code","source":"# In new Librosa version, you need to use soundfile\nimport soundfile as sf","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sr = 22050 # sample rate\nT = 5.0    # seconds\nt = np.linspace(0, T, int(T*sr), endpoint=False) # time variable\nx = 0.5*np.sin(2*np.pi*220*t)# pure sine wave at 220 Hz\n#Playing the audio\nipd.Audio(x, rate=sr) # load a NumPy array\n#Saving the audio\nsf.write('tone_220.wav', x, sr)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"import sklearn\nspectral_centroids = librosa.feature.spectral_centroid(x, sr=sr)[0]\nspectral_centroids.shape\n# Computing the time variable for visualization\nplt.figure(figsize=(12, 4))\nframes = range(len(spectral_centroids))\nt = librosa.frames_to_time(frames)\n# Normalising the spectral centroid for visualisation\ndef normalize(x, axis=0):\n    return sklearn.preprocessing.minmax_scale(x, axis=axis)\n#Plotting the Spectral Centroid along the waveform\nlibrosa.display.waveplot(x, sr=sr, alpha=0.4)\nplt.plot(t, normalize(spectral_centroids), color='b')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"### .spectral_centroid will return an array with columns equal to a number of frames present in your sample."},{"metadata":{},"cell_type":"markdown","source":"# Spectral Rolloff\n### It is a measure of the shape of the signal. It represents the frequency at which high frequencies decline to 0. To obtain it, we have to calculate the fraction of bins in the power spectrum where 85% of its power is at lower frequencies.\n\n### librosa.feature.spectral_rolloff computes the rolloff frequency for each frame in a signal:"},{"metadata":{"trusted":true},"cell_type":"code","source":"spectral_rolloff = librosa.feature.spectral_rolloff(x+0.01, sr=sr)[0]\nplt.figure(figsize=(12, 4))\nlibrosa.display.waveplot(x, sr=sr, alpha=0.4)\nplt.plot(t, normalize(spectral_rolloff), color='r')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Spectral Bandwidth\n\n#### The spectral bandwidth is defined as the width of the band of light at one-half the peak maximum (or full width at half maximum [FWHM]) and is represented by the two vertical red lines and λSB on the wavelength axis."},{"metadata":{"trusted":true},"cell_type":"code","source":"spectral_bandwidth_2 = librosa.feature.spectral_bandwidth(x+0.01, sr=sr)[0]\nspectral_bandwidth_3 = librosa.feature.spectral_bandwidth(x+0.01, sr=sr, p=3)[0]\nspectral_bandwidth_4 = librosa.feature.spectral_bandwidth(x+0.01, sr=sr, p=4)[0]\nplt.figure(figsize=(15, 9))\nlibrosa.display.waveplot(x, sr=sr, alpha=0.4)\nplt.plot(t, normalize(spectral_bandwidth_2), color='r')\nplt.plot(t, normalize(spectral_bandwidth_3), color='g')\nplt.plot(t, normalize(spectral_bandwidth_4), color='y')\nplt.legend(('p = 2', 'p = 3', 'p = 4'))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"x, sr = librosa.load(audio_data)\n#Plot the signal:\nplt.figure(figsize=(14, 5))\nlibrosa.display.waveplot(x, sr=sr)\n# Zooming in\nn0 = 9000\nn1 = 9100\nplt.figure(figsize=(14, 5))\nplt.plot(x[n0:n1])\nplt.grid()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## How many zero crossings?"},{"metadata":{"trusted":true},"cell_type":"code","source":"zero_crossings = librosa.zero_crossings(x[n0:n1], pad=False)\nprint(sum(zero_crossings))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Mel-Frequency Cepstral Coefficients(MFCCs)"},{"metadata":{"trusted":true},"cell_type":"code","source":"fs=10\nmfccs = librosa.feature.mfcc(x, sr=fs)\nprint(mfccs.shape)\n(20, 97)\n#Displaying  the MFCCs:\nplt.figure(figsize=(15, 7))\nlibrosa.display.specshow(mfccs, sr=sr, x_axis='time')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Chroma feature"},{"metadata":{"trusted":true},"cell_type":"code","source":"hop_length=12\nchromagram = librosa.feature.chroma_stft(x, sr=sr, hop_length=hop_length)\nplt.figure(figsize=(15, 5))\nlibrosa.display.specshow(chromagram, x_axis='time', y_axis='chroma', hop_length=hop_length, cmap='coolwarm')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Modeling\n### I will use some existing kernels. [[Inference] ResNest RFCX Audio Detection](https://www.kaggle.com/kneroma/inference-resnest-rfcx-audio-detection) by [kkiller](https://www.kaggle.com/kneroma)"},{"metadata":{"trusted":true},"cell_type":"code","source":"!pip install  resnest > /dev/null","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"from pathlib import Path\nimport librosa as lb\n\nimport torch\nfrom  torch.utils.data import Dataset, DataLoader\n\nfrom tqdm.notebook import tqdm\n\nfrom resnest.torch import resnest50","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"NUM_CLASSES = 24\nSR = 16_000\nDURATION =  60\nDATA_ROOT = Path(\"../input/rfcx-species-audio-detection\")\nTRAIN_AUDIO_ROOT = Path(\"../input/rfcx-species-audio-detection/train\")\nTEST_AUDIO_ROOT = Path(\"../input/rfcx-species-audio-detection/test\")","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"class MelSpecComputer:\n    def __init__(self, sr, n_mels, fmin, fmax):\n        self.sr = sr\n        self.n_mels = n_mels\n        self.fmin = fmin\n        self.fmax = fmax\n\n    def __call__(self, y):\n\n        melspec = lb.feature.melspectrogram(\n            y, sr=self.sr, n_mels=self.n_mels, fmin=self.fmin, fmax=self.fmax,\n        )\n\n        melspec = lb.power_to_db(melspec).astype(np.float32)\n        return melspec","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def mono_to_color(X, eps=1e-6, mean=None, std=None):\n    X = np.stack([X, X, X], axis=-1)\n\n    # Standardize\n    mean = mean or X.mean()\n    std = std or X.std()\n    X = (X - mean) / (std + eps)\n\n    # Normalize to [0, 255]\n    _min, _max = X.min(), X.max()\n\n    if (_max - _min) > eps:\n        V = np.clip(X, _min, _max)\n        V = 255 * (V - _min) / (_max - _min)\n        V = V.astype(np.uint8)\n    else:\n        V = np.zeros_like(X, dtype=np.uint8)\n\n    return V\n\n\ndef normalize(image, mean=None, std=None):\n    image = image / 255.0\n    if mean is not None and std is not None:\n        image = (image - mean) / std\n    return np.moveaxis(image, 2, 0).astype(np.float32)\n\n\ndef crop_or_pad(y, length, sr, is_train=True):\n    if len(y) < length:\n        y = np.concatenate([y, np.zeros(length - len(y))])\n    elif len(y) > length:\n        if not is_train:\n            start = 0\n        else:\n            start = np.random.randint(len(y) - length)\n\n        y = y[start:start + length]\n\n    y = y.astype(np.float32, copy=False)\n\n    return y","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"class RFCXDataset(Dataset):\n\n    def __init__(self, data, sr, n_mels=128, fmin=0, fmax=None,  is_train=False,\n                 num_classes=NUM_CLASSES, root=None, duration=DURATION):\n\n        self.data = data\n        \n        self.sr = sr\n        self.n_mels = n_mels\n        self.fmin = fmin\n        self.fmax = fmax or self.sr//2\n\n        self.is_train = is_train\n\n        self.num_classes = num_classes\n        self.duration = duration\n        self.audio_length = self.duration*self.sr\n        \n        self.root =  root or (TRAIN_AUDIO_ROOT if self.is_train else TEST_AUDIO_ROOT)\n\n        self.wav_transfos = get_wav_transforms() if self.is_train else None\n\n        self.mel_spec_computer = MelSpecComputer(sr=self.sr, n_mels=self.n_mels, fmin=self.fmin, fmax=self.fmax)\n\n\n    def __len__(self):\n        return len(self.data)\n    \n    def read_index(self, idx, fill_val=1.0, offset=None, use_offset=True):\n        d = self.data.iloc[idx]\n        record, species = d[\"recording_id\"], d[\"species_id\"]\n        try:\n            if use_offset and (self.duration < d[\"duration\"]+1):\n                offset = offset or np.random.uniform(1, int(d[\"duration\"]-self.duration))\n\n            y, _ = lb.load(self.root.joinpath(record).with_suffix(\".flac\").as_posix(),\n                           sr=self.sr, duration=self.duration, offset=offset)\n            \n            if self.wav_transfos is not None:\n                y = self.wav_transfos(y, self.sr)\n            y = crop_or_pad(y, self.audio_length, sr=self.sr)\n            t = np.zeros(self.num_classes)\n            t[species] = fill_val\n        except Exception as e:\n#             print(e)\n            raise ValueError()  from  e\n            y = np.zeros(self.audio_length)\n            t = np.zeros(self.num_classes)\n        \n        return y,t\n            \n        \n\n    def __getitem__(self, idx):\n\n        y, t = self.read_index(idx)\n        \n        \n        melspec = self.mel_spec_computer(y) \n        image = mono_to_color(melspec)\n        image = normalize(image, mean=None, std=None)\n\n        return image, t","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def get_duration(audio_name, root=TEST_AUDIO_ROOT):\n    return lb.get_duration(filename=root.joinpath(audio_name).with_suffix(\".flac\"))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"data = pd.DataFrame({\n    \"recording_id\": [path.stem for path in Path(TEST_AUDIO_ROOT).glob(\"*.flac\")],\n})\ndata[\"species_id\"] = [[] for _ in range(len(data))]\n\nprint(data.shape)\ndata[\"duration\"] = data[\"recording_id\"].apply(get_duration)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Inference"},{"metadata":{"trusted":true},"cell_type":"code","source":"TEST_BATCH_SIZE = 40\nTEST_NUM_WORKERS = 2","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"test_data = RFCXDataset(data=data, sr=SR)\ntest_loader = DataLoader(test_data, batch_size=TEST_BATCH_SIZE, num_workers=TEST_NUM_WORKERS)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"device = torch.device(\"cuda\" if torch.cuda.is_available() else \"cpu\")\n\nnet = resnest50(pretrained=True).to(device)\nn_features = net.fc.in_features\nnet.fc = torch.nn.Linear(n_features, NUM_CLASSES)\nnet = net.to(device)\nnet.load_state_dict(torch.load(\"../input/kkiller-rfcx-species-detection-public-checkpoints/rfcx_resnest50.pth\", map_location=device))\nnet = net.eval()\nnet","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"preds = []\nnet.eval()\nwith torch.no_grad():\n    for (xb, yb) in  tqdm(test_loader):\n        xb, yb = xb.to(device), yb.to(device)\n        o = net(xb)\n        o = torch.sigmoid(o) \n        preds.append(o.detach().cpu().numpy())\npreds = np.vstack(preds)\npreds.shape","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sub = pd.DataFrame(preds, columns=[f\"s{i}\" for i in range(24)])\nsub[\"recording_id\"] = data[\"recording_id\"].values[:len(sub)]\nsub = sub[[\"recording_id\"] + [f\"s{i}\" for i in range(24)]]\nprint(sub.shape)\nsub.head()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"## Trying to do some ensembling with the result from [Giba](https://www.kaggle.com/titericz) is notebook [Tabular XGboost GPU + FFT GPU + Cuml = FAST](https://www.kaggle.com/titericz/0-525-tabular-xgboost-gpu-fft-gpu-cuml-fast)"},{"metadata":{"trusted":true},"cell_type":"code","source":"sub1 = pd.read_csv('../input/tabular-xgboost-gpu-fft-gpu-cuml-fast/submission - 2020-11-20T161440.346.csv')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def get_ensebles(data1,data2,weight):\n    df1=data1.copy()\n    df2=data1.copy()\n    df3=data1.copy()\n    df4=data1.copy()\n    for i in data1.columns[1:]:\n        res1 = []\n        res2 = []\n        res3 = []\n        res4 = []\n        l1=[]\n        l2=[]\n        l1=data1[i].tolist()\n        l2=data2[i].tolist()\n        # get min value\n        for j in range(len(data1)):\n            res1.append(max(l1[j],l2[j]))\n            res2.append(min(l1[j],l2[j]))\n            res3.append((l1[j]+l2[j])/2)\n            res4.append(l1[j]*weight+l2[j]*0.1)\n        df1[i]=res1\n        df2[i]=res2\n        df3[i]=res3\n        df4[i]=res4\n    return df1,df2,df3,df4\n\ndf1,df2,df3,df4= get_ensebles(sub,sub1,0.9)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sub.to_csv(\"submission.csv\", index=False)\ndf1.to_csv(\"submission1.csv\", index=False)\ndf2.to_csv(\"submission2.csv\", index=False)\ndf3.to_csv(\"submission3.csv\", index=False)\ndf4.to_csv(\"submission4.csv\", index=False)","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}