update

.gitignore

@@ -13,3 +13,4 @@ annot_distrib.pdf
 annot_distrib.tex
 humpback/annot
 humpback_CARIMAM/
+dolphin/zips

good_species.txt

@@ -4,4 +4,5 @@ california_thrashers
 cassin_vireo
 black-headed_grosbeaks
 humpback
-zebra_finch
+dolphin
+otter

models.py

@@ -41,12 +41,12 @@ meta = {
     'sampleDur': 2
   },
   'humpback':{
-    'nfft': 512,
+    'nfft': 1024,
     'sr': 11025,
     'sampleDur': 2
   },
   'humpback2':{
-    'nfft': 512,
+    'nfft': 2048,
     'sr': 11025,
     'sampleDur': 2
   },
@@ -57,12 +57,17 @@ meta = {
   },
   'otter':{
     'nfft':2048,
-    'sr':48000,
-    'sampleDur':0.5
+    'sr':96000,
+    'sampleDur':1
+  },
+  'dolphin':{
+    'nfft':512,
+    'sr':96000,
+    'sampleDur':2
   }
 }
 
-vgg16 = torchmodels.vgg16(pretrained=True) # weights=torchmodels.VGG16_Weights.DEFAULT)
+vgg16 = torchmodels.vgg16(weights=torchmodels.VGG16_Weights.DEFAULT)
 vgg16 = vgg16.features[:13]
 for nm, mod in vgg16.named_modules():
     if isinstance(mod, nn.MaxPool2d):
@@ -82,6 +87,13 @@ frontend = {
     nn.InstanceNorm2d(1),
     u.Croper2D(n_mel, 128)
   ),
+  'logMel_vggish': lambda sr, nfft, sampleDur, n_mel : nn.Sequential(
+    STFT(nfft, int((sampleDur*sr - nfft)/96)),
+    MelFilter(sr, nfft, n_mel, 0, sr//2),
+    Log1p(7, trainable=False),
+    nn.InstanceNorm2d(1),
+    u.Croper2D(n_mel, 128)
+  ),
   'logSTFT': lambda sr, nfft, sampleDur, n_mel : nn.Sequential(
     STFT(nfft, int((sampleDur*sr - nfft)/128)),
     Log1p(7, trainable=False),
@@ -117,8 +129,28 @@ sparrow_VQ_encoder = lambda nfeat, shape : nn.Sequential(
   u.Reshape(nfeat * shape[0] * shape[1])
 )
 
+sparrow_encoder_maxPool = lambda nfeat, shape : nn.Sequential(
+  nn.Conv2d(1, 32, 3, bias=False, padding=1),
+  nn.BatchNorm2d(32),
+  nn.ReLU(True),
+  nn.Conv2d(32, 32, 3, bias=False, padding=1),
+  nn.BatchNorm2d(32),
+  nn.MaxPool2d((4, 4)),
+  nn.ReLU(True),
+  nn.Conv2d(32, 32, 3, bias=False, padding=1),
+  nn.BatchNorm2d(32),
+  nn.ReLU(True),
+  nn.Conv2d(32, 32, 3, bias=False, padding=1),
+  nn.BatchNorm2d(32),
+  nn.ReLU(True),
+  nn.Conv2d(32, nfeat, 3, bias=False, padding=1),
+  nn.MaxPool2d((4, 4)),
+  u.Reshape(nfeat * shape[0] * shape[1])
+)
+
+
 sparrow_encoder = lambda nfeat, shape : nn.Sequential(
-  nn.Conv2d(1, 32, 3, stride=2, bias=False, padding=(1)),
+  nn.Conv2d(1, 32, 3, stride=2, bias=False, padding=1),
   nn.BatchNorm2d(32),
   nn.ReLU(True),
   nn.Conv2d(32, 64, 3, stride=2, bias=False, padding=1),
@@ -177,7 +209,7 @@ sparrow_decoder = lambda nfeat, shape : nn.Sequential(
   nn.Conv2d(32, 1, (3, 3), bias=False, padding=1),
   nn.BatchNorm2d(1),
   nn.ReLU(True),
-  nn.Conv2d(1, 1, (3, 3), bias=False, padding=1),
+  nn.Conv2d(1, 1, (3, 3), bias=False, padding=1)
 )
 
 sparrow_decoder_old = lambda nfeat, shape : nn.Sequential(

plot_annot_distrib.py

@@ -2,6 +2,15 @@ import matplotlib.pyplot as plt
 import pandas as pd, numpy as np
 
 species = np.loadtxt('good_species.txt', dtype=str)
+info = {
+    'bengalese_finch1': ['bengalese finch', 'nicholson2017bengalese', 'bird'],
+    'bengalese_finch2': ['bengalese finch', 'koumura2016birdsongrecognition', 'bird'],
+    'california_thrashers': ['california trashers', 'arriaga2015bird', 'bird'],
+    'cassin_vireo': ['cassin vireo', 'arriaga2015bird', 'bird'],
+    'black-headed_grosbeaks': ['black-headed grosbeaks', 'arriaga2015bird', 'bird'],
+    'zebra_finch': ['zebra finch', 'elie2018zebra', 'bird'],
+    'humpback': ['humpback whale', 'malige2021use', 'cetacean'],
+}
 
 fig, ax = plt.subplots(nrows=4, ncols=3, figsize=(10, 10))
 for i, specie in enumerate(species):
@@ -12,10 +21,10 @@ plt.tight_layout()
 plt.savefig('annot_distrib.pdf')
 
 
-a = "Specie & \# Classes & \# Annotated samples & \# Samples & Proportion of annotations\\\\ \hline \n"
+a = "Specie & \# Classes & \# Samples & Annotations \% \\\\ \hline \n"
 for specie in species:
     df = pd.read_csv(f'{specie}/{specie}.csv')
-    a += f"{specie.replace('_',' ')} & {df.label.nunique()} & {(~df.label.isna()).sum()} & {len(df)} & {int(100*(~df.label.isna()).sum()/len(df))} \\\\ \hline \n"
+    a += f"{info[specie][0]} \cite{{{info[specie][1]}}} & {df.label.nunique()} & {len(df)} & {int(100*(~df.label.isna()).sum()/len(df))} \\\\ \hline \n"
 f = open('annot_distrib.tex', 'w')
 f.write(a)
 f.close()

plot_results_hdbcsan.py

@@ -5,7 +5,7 @@ import numpy as np
 from sklearn import metrics
 
 species = np.loadtxt('good_species.txt', dtype=str)
-frontends = ['16_pcenMel128', '16_logMel128', '16_logSTFT', '16_Mel128', '8_pcenMel64', '32_pcenMel128']
+frontends = ['16_logMel128', '16_logSTFT', '16_Mel128', '16_pcenMel128', '8_pcenMel64', '32_pcenMel128', '64_pcenMel128']
 plt.figure()
 for specie in species:
     df = pd.read_csv(f'{specie}/{specie}.csv')
@@ -30,3 +30,4 @@ plt.grid()
 plt.tight_layout()
 plt.legend()
 plt.savefig('NMIs_hdbscan.pdf')
+plt.close()

print_annot.py

@@ -8,7 +8,7 @@ import torch
 parser = argparse.ArgumentParser()
 parser.add_argument("specie", type=str)
 parser.add_argument("-frontend", type=str, default='logMel')
-parser.add_argument("-nMel", type=int, default=64)
+parser.add_argument("-nMel", type=int, default=128)
 args = parser.parse_args()
 
 meta = models.meta[args.specie]

run_baseline.py

@@ -12,8 +12,8 @@ args = parser.parse_args()
 
 df = pd.read_csv(f'{args.specie}/{args.specie}.csv')
 
-def norm(arr):
-    return (arr - np.mean(arr) ) / np.std(arr)
+norm = lambda arr: (arr - np.mean(arr) ) / np.std(arr)
+
 meta = models.meta[args.specie]
 
 feats = ['fund', 'cvfund', 'maxfund', 'minfund', 'meansal', 'meanspect', 'stdspect', 'skewspect',\
@@ -32,7 +32,7 @@ def process(idx):
     if fs != meta['sr']:
         sig = resample(sig, int(len(sig)/fs*meta['sr']))
     sound = BioSound(soundWave=norm(sig), fs=fs)
-    sound.spectroCalc(max_freq=meta['sr']//2)
+    sound.spectroCalc(max_freq=meta['sr']//2, spec_sample_rate=128//meta['sampleDur'])
     sound.rms = sound.sound.std() 
     sound.ampenv(cutoff_freq = 20, amp_sample_rate = 1000)
     sound.spectrum(f_high=meta['sr']//2 - 1)
@@ -43,10 +43,10 @@ def process(idx):
     
     return [sound.__dict__[f] for f in feats]
 
-res = p_tqdm.p_map(process, df.index[:10])
+res = p_tqdm.p_map(process, df.index[:100], num_cpus=16)
 
-for i, mr in zip(df.index[:10], res):
+for i, mr in zip(df.index[:100], res):
     for f, r in zip(feats, mr):
         df.loc[i, f] = r
 
-df.to_csv(f'{args.specie}/{args.specie}_biosound.csv', index=False)
\ No newline at end of file
+df.to_csv(f'{args.specie}/{args.specie}_biosound.csv', index=False)

run_vggish.py

+from sklearn import metrics
+import matplotlib.pyplot as plt
+import umap, hdbscan
+from tqdm import tqdm
+import argparse, os
+import models, utils as u
+import pandas as pd, numpy as np, torch
+torch.multiprocessing.set_sharing_strategy('file_system')
+
+
+parser = argparse.ArgumentParser()
+parser.add_argument("specie", type=str)
+args = parser.parse_args()
+
+df = pd.read_csv(f'{args.specie}/{args.specie}.csv')
+
+meta = models.meta[args.specie]
+
+if not os.path.isfile(f'{args.specie}/encodings_vggish.npy'):
+    gpu = torch.device('cuda')
+    frontend = models.frontend['logMel_vggish'](meta['sr'], meta['nfft'], meta['sampleDur'], 64)
+    vggish = torch.hub.load('harritaylor/torchvggish', 'vggish')
+    # vggish.preprocess = False
+    vggish.postprocess = False
+    model = torch.nn.Sequential(frontend, vggish).to(gpu)
+    model.eval()
+    loader = torch.utils.data.DataLoader(u.Dataset(df, f'{args.specie}/audio/', 16000, 1), batch_size=1, shuffle=True, num_workers=8, collate_fn=u.collate_fn)
+    with torch.no_grad():
+        encodings, idxs = [], []
+        for x, idx in tqdm(loader, desc='test '+args.specie, leave=False):
+            # encoding = model(x.to(gpu))
+            encoding = vggish(x.numpy().squeeze(0), fs=16000)
+            idxs.extend(idx)
+            encodings.extend(encoding.cpu().detach())
+
+    idxs, encodings = np.array(idxs), np.stack(encodings)
+    X = umap.UMAP(n_jobs=-1).fit_transform(encodings)
+    np.save(f'{args.specie}/encodings_vggish.npy', {'idxs':idxs, 'encodings':encodings, 'umap':X})
+else:
+    dic = np.load(f'{args.specie}/encodings_vggish.npy', allow_pickle=True).item()
+    idxs, encodings, X = dic['idxs'], dic['encodings'], dic['umap']
+
+clusters = hdbscan.HDBSCAN(min_cluster_size=50, min_samples=5, cluster_selection_epsilon=0.05, core_dist_n_jobs=-1, cluster_selection_method='leaf').fit_predict(X)
+df.loc[idxs, 'cluster'] = clusters.astype(int)
+mask = ~df.loc[idxs].label.isna()
+
+#print('Found clusters : \n', pd.Series(clusters).value_counts())
+
+plt.figure(figsize=(20, 10))
+plt.scatter(X[clusters==-1,0], X[clusters==-1,1], s=2, alpha=.2, color='Grey')
+plt.scatter(X[clusters!=-1,0], X[clusters!=-1,1], s=2, c=clusters[clusters!=-1], cmap='tab20')
+plt.tight_layout()
+plt.savefig(f'{args.specie}/vggish_projection_clusters.png')
+
+plt.figure(figsize=(20, 10))
+plt.scatter(X[~mask,0], X[~mask,1], s=2, alpha=.2, color='Grey')
+for l, grp in df.groupby('label'):
+    plt.scatter(X[df.loc[idxs].label==l, 0], X[df.loc[idxs].label==l, 1], s=4, label=l)
+plt.legend()
+plt.tight_layout()
+plt.savefig(f'{args.specie}/vggish_projection_labels.png')
+
+
+clusters, labels = clusters[mask], df.loc[idxs[mask]].label
+print('Silhouette', metrics.silhouette_score(encodings[mask], clusters))
+print('NMI', metrics.normalized_mutual_info_score(labels, clusters))
+print('Homogeneity', metrics.homogeneity_score(labels, clusters))
+print('Completeness', metrics.completeness_score(labels, clusters))
+print('V-Measure', metrics.v_measure_score(labels, clusters))
+
+labelled = df[~df.label.isna()]
+for l, grp in labelled.groupby('label'):
+    best = (grp.groupby('cluster').fn.count() / labelled.groupby('cluster').fn.count()).idxmax()
+    print(f'Best precision for {l} is for cluster {best} with {(df.cluster==best).sum()} points, \
+with precision {((labelled.cluster==best)&(labelled.label==l)).sum()/(labelled.cluster==best).sum():.2f} and recall {((labelled.cluster==best)&(labelled.label==l)).sum()/(labelled.label==l).sum():.2f}')

test_AE.py

@@ -7,6 +7,7 @@ from sklearn import metrics
 import umap, hdbscan
 torch.multiprocessing.set_sharing_strategy('file_system')
 
+
 parser = argparse.ArgumentParser()
 parser.add_argument("specie", type=str)
 parser.add_argument("-bottleneck", type=int, default=16)
@@ -21,8 +22,6 @@ df = pd.read_csv(f'{args.specie}/{args.specie}.csv')
 print(f'Tests for model {modelname}')
 print(f'{len(df)} available vocs')
 
-loader = torch.utils.data.DataLoader(u.Dataset(df, f'{args.specie}/audio/', meta['sr'], meta['sampleDur']), batch_size=64, shuffle=True, num_workers=8, collate_fn=u.collate_fn)
-
 if os.path.isfile(f'{args.specie}/encodings_{modelname[:-4]}npy'):
     dic = np.load(f'{args.specie}/encodings_{modelname[:-4]}npy', allow_pickle=True).item()
     idxs, encodings, X = dic['idxs'], dic['encodings'], dic['umap']
@@ -34,6 +33,7 @@ else:
     model = torch.nn.Sequential(frontend, encoder, decoder).to(gpu)
     model.load_state_dict(torch.load(f'{args.specie}/{modelname}'))
     model.eval()
+    loader = torch.utils.data.DataLoader(u.Dataset(df, f'{args.specie}/audio/', meta['sr'], meta['sampleDur']), batch_size=64, shuffle=True, num_workers=8, collate_fn=u.collate_fn)
     with torch.no_grad():
         encodings, idxs = [], []
         for x, idx in tqdm(loader, desc='test '+args.specie, leave=False):
@@ -45,7 +45,7 @@ else:
     X = umap.UMAP(n_jobs=-1).fit_transform(encodings)
     np.save(f'{args.specie}/encodings_{modelname[:-4]}npy', {'idxs':idxs, 'encodings':encodings, 'umap':X})
 
-clusters = hdbscan.HDBSCAN(min_cluster_size=10, min_samples=5, cluster_selection_epsilon=0.0, core_dist_n_jobs=-1, cluster_selection_method='eom').fit_predict(X)
+clusters = hdbscan.HDBSCAN(min_cluster_size=50, min_samples=5, cluster_selection_epsilon=0.05, core_dist_n_jobs=-1, cluster_selection_method='leaf').fit_predict(X)
 df.loc[idxs, 'cluster'] = clusters.astype(int)
 mask = ~df.loc[idxs].label.isna()
 
@@ -74,7 +74,15 @@ print('Completeness', metrics.completeness_score(labels, clusters))
 print('V-Measure', metrics.v_measure_score(labels, clusters))
 
 labelled = df[~df.label.isna()]
+goodClusters = []
 for l, grp in labelled.groupby('label'):
-    best = (grp.groupby('cluster').fn.count() / labelled.groupby('cluster').fn.count()).idxmax()
+    precisions = grp.groupby('cluster').fn.count() / labelled.groupby('cluster').fn.count()
+    best = precisions.idxmax()
+    goodClusters.extend(precisions[precisions > 0.9].index)
     print(f'Best precision for {l} is for cluster {best} with {(df.cluster==best).sum()} points, \
-with precision {((labelled.cluster==best)&(labelled.label==l)).sum()/(labelled.cluster==best).sum():.2f} and recall {((labelled.cluster==best)&(labelled.label==l)).sum()/(labelled.label==l).sum():.2f}')
+    with precision {((labelled.cluster==best)&(labelled.label==l)).sum()/(labelled.cluster==best).sum():.2f}\
+     and recall {((labelled.cluster==best)&(labelled.label==l)).sum()/(labelled.label==l).sum():.2f}')
+
+
+print(f'{len(goodClusters)} clusters would sort {df.cluster.isin(goodClusters).sum()/len(df)*100:.0f}% of samples')
+print(f'{len(goodClusters)/df.label.nunique():.1f} cluster per label in avg)')
\ No newline at end of file

train_AE.py

@@ -67,8 +67,10 @@ for epoch in range(100_000//len(loader)):
 
         if len(loss) > 2000 and np.median(loss[-2000:-1000]) < np.median(loss[-1000:]):
             print('Early stop')
+            torch.save(model.state_dict(), f'{args.specie}/{modelname}')
             exit()
-
+        step += 1
+        continue
         # TEST ROUTINE
         if step % 500 == 0:
             # Plot reconstructions