Commit 20d84de8 authored by lucas_miranda's avatar lucas_miranda
Browse files

Added CustomStopper class to train_utils.py, to start early stopping only after annealing is over

parent ab55fbf2
Pipeline #89197 passed with stage
in 29 minutes and 38 seconds
......@@ -658,7 +658,9 @@ class SEQ_2_SEQ_GMVAE:
generator = Model_B3(generator)
generator = Model_D5(generator)
generator = Model_B4(generator)
generator = Dense(tfpl.IndependentNormal.params_size(input_shape[2:]))(generator)
generator = Dense(tfpl.IndependentNormal.params_size(input_shape[2:]))(
generator
)
x_decoded_mean = tfpl.IndependentNormal(
event_shape=input_shape[2:],
convert_to_tensor_fn=tfp.distributions.Distribution.mean,
......@@ -675,15 +677,6 @@ class SEQ_2_SEQ_GMVAE:
model_metrics = {"vae_reconstruction": ["mae", "mse"]}
loss_weights = [1.0]
# x_decoded_mean = TimeDistributed(Dense(input_shape[2]), name="vae_prediction")(
# generator
# )
#
# model_outs = [x_decoded_mean]
# model_losses = [Huber(delta=self.delta, reduction="sum_over_batch_size")]
# model_metrics = {"vae_reconstruction": ["mae", "mse"]}
# loss_weights = [1.0]
if self.predictor > 0:
# Define and instantiate predictor
predictor = Dense(
......
......@@ -104,8 +104,8 @@ def get_callbacks(
("P" if predictor > 0 and variational else ""),
("_Pheno" if phenotype_class > 0 else ""),
("_loss={}".format(loss) if variational else ""),
("_encoding={}".format(logparam["encoding"])),
("_k={}".format(logparam["k"])),
("_encoding={}".format(logparam["encoding"]) if logparam is not None else ""),
("_k={}".format(logparam["k"]) if logparam is not None else ""),
(datetime.now().strftime("%Y%m%d-%H%M%S")),
)
......
%% Cell type:code id: tags:
``` python
%load_ext autoreload
%autoreload 2
```
%% Cell type:code id: tags:
``` python
import os
os.chdir(os.path.dirname("../"))
```
%% Cell type:code id: tags:
``` python
import deepof.data
import deepof.models
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
import seaborn as sns
from sklearn.preprocessing import StandardScaler, MinMaxScaler
import tensorflow as tf
import tqdm.notebook as tqdm
from ipywidgets import interact
```
%% Cell type:markdown id: tags:
# Retrieve phenotypes
%% Cell type:code id: tags:
``` python
flatten = lambda t: [item for sublist in t for item in sublist]
```
%% Cell type:code id: tags:
``` python
# Load first batch
dset11 = pd.ExcelFile(
"../../Desktop/deepof-data/tagged_videos/Individual_datasets/DLC_batch_1/DLC_single_CDR1_1/1.Openfield_data-part1/JB05.1-OF-SI-part1.xlsx"
)
dset12 = pd.ExcelFile(
"../../Desktop/deepof-data/tagged_videos/Individual_datasets/DLC_batch_1/DLC_single_CDR1_1/2.Openfielddata-part2/AnimalID's-JB05.1-part2.xlsx"
)
dset11 = pd.read_excel(dset11, "Tabelle2")
dset12 = pd.read_excel(dset12, "Tabelle2")
dset11.Test = dset11.Test.apply(lambda x: "Test {}_s1.1".format(x))
dset12.Test = dset12.Test.apply(lambda x: "Test {}_s1.2".format(x))
dset1 = {"CSDS":list(dset11.loc[dset11.Treatment.isin(["CTR+CSDS","NatCre+CSDS"]), "Test"]) +
list(dset12.loc[dset12.Treatment.isin(["CTR+CSDS","NatCre+CSDS"]), "Test"]),
"NS": list(dset11.loc[dset11.Treatment.isin(["CTR+nonstressed","NatCre+nonstressed"]), "Test"]) +
list(dset12.loc[dset12.Treatment.isin(["CTR+nonstressed","NatCre+nonstressed"]), "Test"]),}
dset1inv = {}
for i in flatten(list(dset1.values())):
if i in dset1["CSDS"]:
dset1inv[i] = "CSDS"
else:
dset1inv[i] = "NS"
assert len(dset1inv) == dset11.shape[0] + dset12.shape[0], "You missed some labels!"
```
%% Cell type:code id: tags:
``` python
# Load second batch
dset21 = pd.read_excel(
"../../Desktop/deepof-data/tagged_videos/Individual_datasets/DLC_batch_2/Part1/2_Single/stressproject22.04.2020genotypes-openfieldday1.xlsx"
)
dset22 = pd.read_excel(
"../../Desktop/deepof-data/tagged_videos/Individual_datasets/DLC_batch_2/Part2/2_Single/OpenFieldvideos-part2.xlsx"
)
dset21.Test = dset21.Test.apply(lambda x: "Test {}_s2.1".format(x))
dset22.Test = dset22.Test.apply(lambda x: "Test {}_s2.2".format(x))
dset2 = {"CSDS":list(dset21.loc[dset21.Treatment == "Stress", "Test"]) +
list(dset22.loc[dset22.Treatment == "Stressed", "Test"]),
"NS": list(dset21.loc[dset21.Treatment == "Nonstressed", "Test"]) +
list(dset22.loc[dset22.Treatment == "Nonstressed", "Test"])}
dset2inv = {}
for i in flatten(list(dset2.values())):
if i in dset2["CSDS"]:
dset2inv[i] = "CSDS"
else:
dset2inv[i] = "NS"
assert len(dset2inv) == dset21.shape[0] + dset22.shape[0], "You missed some labels!"
```
%% Cell type:code id: tags:
``` python
# Load third batch
dset31 = pd.read_excel(
"../../Desktop/deepof-data/tagged_videos/Individual_datasets/DLC_batch_3/1.Day2OF-SIpart1/JB05 2Female-ELS-OF-SIpart1.xlsx"
)
dset32 = pd.read_excel(
"../../Desktop/deepof-data/tagged_videos/Individual_datasets/DLC_batch_3/2.Day3OF-SIpart2/JB05 2FEMALE-ELS-OF-SIpart2.xlsx"
)
dset31.Test = dset31.Test.apply(lambda x: "Test {}_s3.1".format(x))
dset32.Test = dset32.Test.apply(lambda x: "Test {}_s3.2".format(x))
dset3 = {"CSDS":[],
"NS": list(dset31.loc[:, "Test"]) +
list(dset32.loc[:, "Test"])}
dset3inv = {}
for i in flatten(list(dset3.values())):
if i in dset3["CSDS"]:
dset3inv[i] = "CSDS"
else:
dset3inv[i] = "NS"
assert len(dset3inv) == dset31.shape[0] + dset32.shape[0], "You missed some labels!"
```
%% Cell type:code id: tags:
``` python
# Load fourth batch
dset41 = os.listdir("../../Desktop/deepof-data/tagged_videos/Individual_datasets/DLC_batch_4/JB05.4-OpenFieldvideos/")
# Remove empty video!
dset41 = [vid for vid in dset41 if "52" not in vid]
dset4 = {"CSDS":[],
"NS": [i[:-4]+"_s4" for i in dset41]}
dset4inv = {}
for i in flatten(list(dset4.values())):
if i in dset4["CSDS"]:
dset4inv[i] = "CSDS"
else:
dset4inv[i] = "NS"
assert len(dset4inv) == len(dset41), "You missed some labels!"
```
%% Cell type:code id: tags:
``` python
# Merge phenotype dicts and serialise!
aggregated_dset = {**dset1inv, **dset2inv, **dset3inv, **dset4inv}
```
%% Cell type:code id: tags:
``` python
from collections import Counter
print(Counter(aggregated_dset.values()))
print(115+52)
```
%%%% Output: stream
Counter({'NS': 115, 'CSDS': 52})
167
%% Cell type:markdown id: tags:
# Define and run project
%% Cell type:code id: tags:
``` python
%%time
deepof_main = deepof.data.project(path=os.path.join("..","..","Desktop","deepof_single_topview"),
deepof_main = deepof.data.project(path=os.path.join("..","..","Desktop","deepof-data","deepof_single_topview"),
smooth_alpha=0.99,
arena_dims=[380],
exclude_bodyparts=["Tail_1", "Tail_2", "Tail_tip", "Tail_base", "Spine_2"]
#exclude_bodyparts=["Tail_1", "Tail_2", "Tail_tip", "Tail_base", "Spine_2"]
#exp_conditions=dset2inv
)
```
%%%% Output: stream
CPU times: user 28.4 s, sys: 5 s, total: 33.4 s
Wall time: 10.8 s
CPU times: user 28.1 s, sys: 4.99 s, total: 33.1 s
Wall time: 7.5 s
%% Cell type:code id: tags:
``` python
%%time
deepof_main = deepof_main.run(verbose=True)
print(deepof_main)
```
%%%% Output: stream
Loading trajectories...
Smoothing trajectories...
Computing distances...
Computing angles...
Done!
deepof analysis of 167 videos
CPU times: user 23.5 s, sys: 1.76 s, total: 25.3 s
Wall time: 29.7 s
CPU times: user 46.6 s, sys: 5.09 s, total: 51.7 s
Wall time: 53 s
%% Cell type:code id: tags:
``` python
all_quality = pd.concat([tab for tab in deepof_main.get_quality().values()]).droplevel("scorer", axis=1)
```
%% Cell type:code id: tags:
``` python
all_quality.boxplot(rot=45)
plt.ylim(0.99985, 1.00001)
plt.show()
```
%%%% Output: display_data
![](data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAZ8AAAEnCAYAAAB2e06MAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO2debgdRZn/P28WtkAgEb2EgEQlCCG4cMMmKnsIixNUxgGXBDVkFFRc+A1xUCPgaBgZF0YEFCKJOqLDqES2GGKuyww4iStqyCQ4KnFARoIwwXEb3t8fbx1v53DOPUsv956b7+d5+jlddarfqu6urrfqrberzd0RQgghqmTMcBdACCHE9oeUjxBCiMqR8hFCCFE5Uj5CCCEqR8pHCCFE5Uj5CCGEqBwpnxaY2V+a2Y/N7AkzmzXc5RFCiNGAlE8GMzvWzG6oi/4R8DLgG9WXSAghRifjhrsAIx13Xw9gZsNdFCGEGDVo5COEEKJyNPIBzOzbwI7ArsBkM/t++usid185fCUTQojRiZQP4O5HQMz5AOe4+znDWiAhhBjlyOwmhBCicqR8WmBmLzWzzcBRwK1mJjOcEELkxPRJBSGEEFWjkY8QQojKkfIRQghROdu9t9uee+7p06ZNazv9448/zoQJE0opS5myJV/yJV/yi5L9ne9859fu/tRcmbr7dr319/d7J6xZs6aj9CNFtuRLvuRLflGygXWes+2V2U0IIUTlSPkIIYSoHCkfIYQQlSPlI4QQonKkfIQQQlSOlI8QQojKkfIRQghROVI+QgghKqcQ5WNmc8xsg5ltMrNFDf7f0cw+n/7/tplNy/z3zhS/wcxObiXTzJ6RZGxKMndolYcQQoiRRW7lY2ZjgauAU4AZwNlmNqMu2euBR9x9f+DDwOXp2BnAWcDBwBzg42Y2toXMy4EPJ1mPJNlN8xBCCDHyKGLkcziwyd1/6u5/AG4E5talmQssS/s3ASeYmaX4G9399+7+n8CmJK+hzHTM8UkGSeYZLfIQQggxwihC+UwF7s+EN6e4hmnc/U/Ao8BThji2WfxTgN8kGfV5NctDCCHECGO7XNXazBYCCwH6+voYGBh4UprzVz/O43+M/Z9ffnpTWftddAsAE8bDVSd0vuLscccd1/S/NWvWdCxP8ltT9r2V/PbpxfozEuT3QtlbUYTy+SWwbya8T4prlGazmY0DdgcebnFso/iHgT3MbFwa3WTTN8vjSbj7J4BPAMyaNcuPPfbYJ6V5/I5b+dmS0yKwZPBrrwMDAzRKP23RrQ3jW+GZL8lOW5TJsyDKkP/cS77Ko/8brVOtAWrEOXc8DsDuO4/nB4tnd5VXGeUv+96O+fkh7Jb2Z94wc4iUg745xx57T9vyyy5/2fJ7vf6MFvlll70VRSiftcB0M3sGoQDOAl5Zl2YFMB+4CzgT+Jq7u5mtAP7JzD4E7A1MB/4dsEYy0zFrkowbk8ybh8qjgPMTdTz6v39sWFGHapw6Ids4tSMrT+NUBv+zfkmp12e3gxZxyLInOZUGy54ctdtBANU2LENRdv0RvUFu5ePufzKzNwErgbHAUnf/sZldSnzzYQVwPfBpM9sEbCGUCSndF4CfAH8Cznf3/wNoJDNleRFwo5m9D/hekk2zPLql1x/wXqbsxqnX762UmxgNFDLn4+63AbfVxb0ns/874C+bHPt3wN+1IzPF/5TwhquPb5pHN5T9gJdJr48cyqaX720V9LpyK7v+6/kqhu3S4WC4Kbvyyqwx/DS9pnc0vr/bE2Urt7Lrv56vYpDyGQZ6vfLKLDM0zSZuh2NSV4iRipSP6JheN8sIIYYfKR8x4tCcjBCjHykfIXoQzSmJXkfKR4iSyS4xaGm52zyvoGlOSXTLSPLUk/IRokSarW1rZrkUkBi9dKIgetkTVspHCCFGEJ0oiF6e79SXTIUQQlSOlI8QQojKkdlNbJfIW2z0UvZ7YnoPrRikfEYhejiGRt5io5uy3xPTe2jFIOUzDJStHPRwCCFGOlI+w4CUgxjpyCwpykbKRwixDTJLiiqQ8hGixyl6BQUhqkCu1kL0MEOtoCDESEbKRwghROXkUj5mNtnMVpnZxvQ7qUm6+SnNRjObn4nvN7N7zGyTmV1pqbvWTK4FV6b0PzSzQ1P888zsLjP7cYr/qzznJYQQolzyjnwWAavdfTqwOoW3wcwmA4uBI4DDgcUZJXU1cC4wPW1zWsg9JZN2YToe4LfAPHc/OMn4iJntkfPchBBClERe5TOXwTdTlgFnNEhzMrDK3be4+yPAKmCOmU0BJrr73R6zo8szxzeTOxdY7sHdwB5mNsXd/8PdNwK4+38BDwFPzXluQgghSiKvt1ufuz+Q9h8E+hqkmQrcnwlvTnFT0359/FBym8mqpcXMDgd2AO7r9GTEyEHvmYxudH9FS+VjZncCezX46+JswN3dzAr37+xEbhpNfRqY7+5PDJFuIWG2o6+vj4GBgYbpGsVv3bq1o/TNKFN2FfI7aTwmjO9M/g1zJjSMP+eOx5v+12n5m1GUnJGQVxnyi5BZxf3t9eerE/kjrext4+5db8AGYEranwJsaJDmbODaTPjaFDcFuLdRumZya8c2yX8i8F3gzE7Oob+/3xux30W3NIxfs2ZNR+mrll2F/E7zLYpelA803XpBfpZeuf69/nx1In+4yg6s8xy6w91zz/msAGrea/OBmxukWQnMNrNJydFgNrDSw6z2mJkdmbzc5mWObyZ3BTAveb0dCTzq7g+Y2Q7Al4j5oJtynpMQQoiSyat8lgAnmdlG4MQUxsxmmdl1AO6+BbgMWJu2S1McwHnAdcAmYo7m9qHkArcBP03pP5mOB3gF8GLgHDP7ftqel/PchBBClEQuhwN3fxg4oUH8OmBBJrwUWNok3cwO5DpwfoP4zwCf6bD4QgghhgmtcCCEEKJytLCo2O7Rwpxie2EkfWhSykds1wy1MKcUkBhtjKRvicnsJoQQonKkfIQQQlSOlI8QQojKkfIRQghROVI+QgghKkfebkIIMYLoxB26TFfospHyEUIMid6DqpZO3KHLdIUuG5ndhBBNGeo9KCHyoJHPEOiDV0IIUQ5SPk1oNOyFUEjN/hNCCNEeMrsJIYSoHCkfIYQQlSPlI4QQonKkfIQQQlSOlI8QQojKya18zGyyma0ys43pd1KTdPNTmo1mNj8T329m95jZJjO70tILBM3kWnBlSv9DMzu0Lp+JZrbZzD6W99yEEEKUQxEjn0XAanefDqxO4W0ws8nAYuAI4HBgcUZJXQ2cC0xP25wWck/JpF2Yjs9yGfCNAs5LCCFESRShfOYyuOLQMuCMBmlOBla5+xZ3fwRYBcwxsynARHe/22O9juWZ45vJnQss9+BuYI8kBzPrB/qArxZwXkIIIUqiCOXT5+4PpP0Hica/nqnA/Znw5hQ3Ne3Xxw8lt6EsMxsD/ANwYZfnIYQQoiLaWuHAzO4E9mrw18XZgLu7mRW+4mCbcs8DbnP3za3WnTKzhYTJjr6+PgYGBjoqT6fp25WxdevWprKLKGOR8jvJt0jKll9lXpJfnsxef746kT/Syt427p5rAzYAU9L+FGBDgzRnA9dmwtemuCnAvY3SNZNbO7Y+f+CzwC+AnwG/Bh4DlrQqf39/v3fCfhfd0lH6TmSsWbOmkDzLlt9pvkVRhnyg6Sb55cvPUnY97JXnqxP5w1V2YJ3n1B1FrO22ApgPLEm/NzdIsxJ4f8bJYDbwTnffYmaPmdmRwLeBecA/tpC7AniTmd1IODA86mGee1UtMzM7B5jl7k0+iiGEGM2UvSiwFh3OTxHKZwnwBTN7PfBz4BUAZjYLeIO7L0hK5jJgbTrmUnffkvbPA24AdgZuT1tTucBtwKnAJuC3wGsLOAchxCih7EWBtehwMeRWPu7+MHBCg/h1wIJMeCmwtEm6mR3IdeD8FmW6gVBoQgghRiBa4UAIIUTl6Hs+YsSjzzgLMfqQ8hEjmqE+4ywFNDpQ52L7RMpnmBhN3jhqPES3qHOx/SLlMwyMJm8cNR5CiG6Qw4EQQojKkfIRQghROVI+QgghKkfKRwghROVI+QghhKgcKR8hhBCVI+UjhBCicqR8hBBCVI6UjxBCiMqR8hFCCFE5Uj5CCCEqR8pHCCFE5Uj5CCGEqJxcysfMJpvZKjPbmH4nNUk3P6XZaGbzM/H9ZnaPmW0ysystLZHcTK4FV6b0PzSzQzOynm5mXzWz9Wb2EzOblufchBBClEfekc8iYLW7TwdWp/A2mNlkYDFwBHA4sDijpK4GzgWmp21OC7mnZNIuTMfXWA580N0PSvk8lPPchBBClERe5TMXWJb2lwFnNEhzMrDK3be4+yPAKmCOmU0BJrr73R4fflmeOb6Z3LnAcg/uBvYwsylmNgMY5+6rANx9q7v/Nue5CSGEKIm8yqfP3R9I+w8CfQ3STAXuz4Q3p7ipab8+fii5zWQdAPzGzL5oZt8zsw+a2dguz0kIIUTJtPySqZndCezV4K+LswF3dzMr/NOVbcodB7wIeD7wC+DzwDnA9Y0Sm9lCwmxHX18fAwMDHZWp0/QjRXYV8qvMS/IlfyTILEN+Izlbt25tGF9E+9VMdjfy28bdu96ADcCUtD8F2NAgzdnAtZnwtSluCnBvo3TN5NaOrc8fOBL4eib+NcBV7ZxDf3+/d8J+F93SUfqRIrss+UDTTfIlf7jlZ+mV56uZnDVr1uTOsxPZQ6UH1nkO3eHuuc1uK4Ca99p84OYGaVYCs81sUnI0mA2s9DCrPWZmRyYvt3mZ45vJXQHMS15vRwKPJjlrifmfp6Z0xwM/yXluQgghSiKv8lkCnGRmG4ETUxgzm2Vm1wG4+xbgMkJBrAUuTXEA5wHXAZuA+4Dbh5IL3Ab8NKX/ZDoed/8/4EJgtZndA1j6XwghxAik5ZzPULj7w8AJDeLXAQsy4aXA0ibpZnYg14Hzm5RlFfCcDoovhBBimNAKB0IIISpHykcIIUTlSPkIIYSoHCkfIYQQlSPlI4QQonKkfIQQQlSOlI8QQojKkfIRQghROVI+QgghKkfKRwghROVI+QghhKgcKR8hhBCVI+UjhBCicqR8hBBCVI6UjxBCiMqR8hFCCFE5Uj5CCCEqJ9eXTLcnzGxw//L4jQ+rCiGE6BSNfNogq3jaiRdCCDE0uZWPmU02s1VmtjH9TmqSbn5Ks9HM5mfi+83sHjPbZGZXWmrRm8m14MqU/odmdmhG1t+b2Y/NbH1WlhBCiJFFESOfRcBqd58OrE7hbTCzycBi4AjgcGBxRkldDZwLTE/bnBZyT8mkXZiOx8xeABwNPAeYCRwGHFPA+QkhhCiYIpTPXGBZ2l8GnNEgzcnAKnff4u6PAKuAOWY2BZjo7nd7TKAszxzfTO5cYLkHdwN7JDkO7ATsAOwIjAd+VcD5CSGEKJgiHA763P2BtP8g0NcgzVTg/kx4c4qbmvbr44eS21CWu99lZmuABwADPubu6xsV2MwWEqMm+vr6GBgYaHWOTclzbBXyqpZfZV6SL/kjQWYZ8hvJ2bp1a8P4TvPsRHY38tvG3VtuwJ3Ajxpsc4Hf1KV9pMHxFwLvyoTfneJmAXdm4l8E3JL2G8oFbgFemIlfneTsD9wK7Jq2u4AXtTq3/v5+bwUxqmq4Fcl+F91SqLwq5Jd9bSRf8ouiV56vZnLWrFmTO89OZA+VHljnbeiOoba2Rj7ufmKz/8zsV2Y2xd0fSOavhxok+yVwbCa8DzCQ4vepi/9l2m8m95fAvg2OeTVwt7tvTeW6HTgK+GY75yiEEKI6ipjzWQHUvNfmAzc3SLMSmG1mk5KjwWxgpYdZ7TEzOzJ5ps3LHN9M7gpgXvJ6OxJ4NMn5BXCMmY0zs/GEs0FDs5sQQojhpQjlswQ4ycw2AiemMGY2y8yuA3D3LcBlwNq0XZriAM4DrgM2AfcBtw8lF7gN+GlK/8l0PMBN6fh7gB8AP3D3rxRwfkIIIQomt8OBuz8MnNAgfh2wIBNeCixtkm5mB3IdOL9B/P8Bf91h8YUQQgwDWuFACCFE5Uj5CCGEqBwpHyGEEJUj5SOEEKJypHyEEEJUjr7nI4QQ2xHTFt3a+I87nhy/+87jSyuHlI8QQmwn/GzJaQ3jpy26tel/ZSGzmxBCiMqR8hFCCFE5Uj5CCCEqR8pHCCFE5Uj5CCGEqBwpHyGEEJUj5SOEEKJypHyEEEJUjpSPEEKIytEKB0KIUY2ZDe5fHr/xTUoxnEj5CCFGLVnFUx9flAKScuuOXGY3M5tsZqvMbGP6ndQk3fyUZqOZzc/E95vZPWa2ycyutHQXm8k1swPN7C4z+72ZXViXxxwz25BkLcpzXkII0Q5DKTcxNHnnfBYBq919OrA6hbfBzCYDi4EjgMOBxRkldTVwLjA9bXNayN0CvAW4oi6PscBVwCnADOBsM5uR89yEEEKURF7lMxdYlvaXAWc0SHMysMrdt7j7I8AqYI6ZTQEmuvvdHmPU5ZnjG8p194fcfS3wx7o8Dgc2uftP3f0PwI1JhhBCiBFIXuXT5+4PpP0Hgb4GaaYC92fCm1Pc1LRfH9+u3HbyEEIIMQJp6XBgZncCezX46+JswN3dzAqfZStDrpktBBYC9PX1MTAw0LWsPMdWIa9q+VXmJfmSP1rlNzp269atDeOLOo8q2wYgvDK63YANwJS0PwXY0CDN2cC1mfC1KW4KcG+jdK3kAu8FLsyEjwJWZsLvBN7Zzjn09/d7K4CmW5Hsd9EthcqrQn7Z10byJX97k9/sOV2zZk3baYvKsxnAOs+hO9w9t9ltBVDzXpsP3NwgzUpgtplNSo4GswlF8QDwmJkdmbzc5mWOb0dulrXAdDN7hpntAJyVZAghhBiB5H3PZwnwBTN7PfBz4BUAZjYLeIO7L3D3LWZ2GaEgAC519y1p/zzgBmBn4Pa0DSV3L2AdMBF4wszeCsxw98fM7E2EohsLLHX3H+c8NyGEECWRS/m4+8PACQ3i1wELMuGlwNIm6WZ2IPdBYJ8mZbkNuK2D4gshhBgmtLabEEKIypHyEUIIUTlSPkIIISpHC4sKIcQIY9qiWxv/cce28bvvPL6C0pSDlM8IQKviCiFq/GzJaQ3jpy26tel/vYiUzzCjJd+FENsjmvMZ5WjJdyHESETKRwghROVI+QghhKgcKR8hhBCVI+UjhBCicqR8hBBCVI6UjxBCiMqR8hFCCFE5Uj5CCCEqR8pHCCFE5Uj5CCGEqBwpHyGEEJWTS/mY2WQzW2VmG9PvpCbp5qc0G81sfia+38zuMbNNZnalpQXHmsk1swPN7C4z+72ZXZiRs6+ZrTGzn5jZj83sgjznJYQQolzyjnwWAavdfTqwOoW3wcwmA4uBI4DDgcUZJXU1cC4wPW1zWsjdArwFuKIumz8B73D3GcCRwPlmNiPnuQkhhCiJvMpnLrAs7S8DzmiQ5mRglbtvcfdHgFXAHDObAkx097s91vdfnjm+oVx3f8jd1wJ/zGbg7g+4+3fT/v8A64GpOc9NCCFESeRVPn3u/kDafxDoa5BmKnB/Jrw5xU1N+/Xx7cptiJlNA54PfLvdY4QQQlRLy4/JmdmdwF4N/ro4G3B3N7PCv1DWiVwz2xX4F+Ct7v7YEOkWAgsB+vr6GBgY6Lp8eY4dTtmSL/mS33vye7m9eRLu3vUGbACmpP0pwIYGac4Grs2Er01xU4B7G6VrJRd4L3BhXdx4YCXw9k7Oob+/31sBNN3yUqZsyZd8ye9t+Vn2u+iWwmV2KxtY5zl0h7vnNrutAGrea/OBmxukWQnMNrNJydFgNrDSw6z2mJkdmbzc5mWOb0fun0nHXw+sd/cP5TkhIYQQ5ZNX+SwBTjKzjcCJKYyZzTKz6wDcfQtwGbA2bZemOIDzgOuATcB9wO0t5O5lZpuBtwPvMrPNZjYROBp4DXC8mX0/bafmPDchhBAl0XLOZyjc/WHghAbx64AFmfBSYGmTdDM7kPsgsE+DonwLsE7KLoQQYvjQCgdCCCEqR8pHCCFE5Uj5CCGEqBwpHyGEEJUj5SOEEKJypHyEEEJUjpSPEEKIypHyEUIIUTlSPkIIISpHykcIIUTlSPkIIYSoHCkfIYQQlSPlI4QQonKkfIQQQlSOlI8QQojKkfIRQghROVI+QgghKkfKRwghROXkUj5mNtnMVpnZxvQ7qUm6+SnNRjObn4nvN7N7zGyTmV1pZjaUXDM70MzuMrPfm9mFDfIZa2bfM7Nb8pyXEEKIcsk78lkErHb36cDqFN4GM5sMLAaOAA4HFmeU1NXAucD0tM1pIXcL8BbgiibluQBYn/OchBBClExe5TMXWJb2lwFnNEhzMrDK3be4+yPAKmCOmU0BJrr73e7uwPLM8Q3luvtD7r4W+GN9Jma2D3AacF3OcxJCCFEyeZVPn7s/kPYfBPoapJkK3J8Jb05xU9N+fXy7cuv5CPA3wBPtFV0IIcRwMa5VAjO7E9irwV8XZwPu7mbmRRWsE7lmdjrwkLt/x8yObSXTzBYCCwH6+voYGBjounx5jh1O2ZIv+ZLfe/J7ub15Eu7e9QZsAKak/SnAhgZpzgauzYSvTXFTgHsbpWslF3gvcGEm/AFi5PQzYqT0W+Az7ZxDf3+/twJouuWlTNmSL/mS39vys+x30S2Fy+xWNrDOc+gOd89tdlsB1LzX5gM3N0izEphtZpOSo8FsYKWHWe0xMzsyebnNyxzfjtw/4+7vdPd93H0acBbwNXd/dY7zEkIIUSJ5lc8S4CQz2wicmMKY2Swzuw7A3bcAlwFr03ZpigM4j3AQ2ATcB9zeQu5eZrYZeDvwLjPbbGYTc56DEEKIimk55zMU7v4wcEKD+HXAgkx4KbC0SbqZHch9ENinRZkGgIGWhRdCCDFsaIUDIYQQlSPlI4QQonKkfIQQQlSOlI8QQojKkfIRQghROVI+QgghKkfKRwghROVI+QghhKgcKR8hhBCVY7FG3PaLmf038PMWyfqH+O87OYtQpmzJl3zJ7235WfYEfl2wzG5l7+fuT82T4XavfDrFzNa5+6xeky35ki/5kj8cspshs5sQQojKkfIRQghROVI+nfOJHpUt+ZIv+ZI/HLIbojkfIYQQlaORjxBCiMqR8hFCCFE5Uj6icszMypJnZrsUKbudPEd6fmZWynNuZuPMbHwZssXoR8qnS6pufBrlbWbTzWy3CvLpLyIfM9sZwN29qOtnZuZp4tLMFgDzzGxsEbIb5LWnmR1kZjt5iZOlmes+2cz2gLhmHcp4ppnNT8c+UbQCMrMZwHXAP5vZydlyV8Foy6vMDllVdJqnlE8X1DV4Z5jZuWY2y8x2qipvM3sR8DFgUll5ZfK5Ctg3j6zUWH3JzF6ekZ37Acnch6OAucA/ufv/5ZVbj5kdBHwN+ABwt5k9O8UX/pCna/MS4GbgKjP7ZIdlPQD4JvBGM/vbJLMwBZTO/TPA14GvAFeY2ZFlKeSMMj7UzOaa2f4VKf/x0Lni7yYvYPe6cF55exQhr1U+ZtZXszakett2HZPy6YJMg3ce8DfAeOAO4Lgq8jazWcCZwPXu/ouy8kqNzBuBT7j7T7ptvMysj2hIfwO80szOhGIUkAUzCVfR3wF/yiOvSR77Ap8HPujuZwArU36lNExmdjhwMfAq4NvAMWa2awciDgO+CJwHHFykAkqjylcRSv5T7n498Cng5en/onvwY1I9OYk4pzOAL5rZa1O9KpyU3+mE4r8hjXYnlJjXqcQI8lPA681sz25kZTqmpwIrzOzjwDtrFoeiyORzOrAGuNbMvgSd1TEpny5IDd7TgaOBE4HfAt8HvlpyvrX7dSJwKrCPlWBzzzQgM4CnACeZWZ+7P9GlyP8GLgcuBL4AnJNVQDnKhwc/At4L7A0cWYLZbTJwjbt/OoX/Fnik6Ic6gwFLgCOAVwInu/tWM3t+Owe7+2eB9wD3EEryYDO7OP33hJnt0G3B0qjys8DS9ByMAX4FTEv/F6KMzWz3JO+J1Ak6F5jn7q8lrv8xwKEpbWEm3PR7GPA+Qtn9HrgAOL7gvMak3yOAtxEj6rXAM4B3dNjZyCqEw4AFwPuBAWLNto/muef1ZU75zCAsDW9K5f+dmd2V/m+vnXB3bW1spHei6uLeDdxE9ITHpLg3Ac8qI2/gKZm4NwK3AbNqeReYz9RM3NHAx9N57dmFzDF14QnAXwG3AK+onRewUxeyFxKN9GJCQbyG6AAcC4wt8PrvBOxVOx9gR+AHwLNT3G6N6kcX1/2gdB7HAj8B7gJ2Tf8dB3wZ6Ouknqayvhj4J+ANwOHA64DxOcpbf09nAjek/aOA03Nej12BK4jOxBjgImA98M5Mmr8mzH47FnB/9wUOSfvTgU8DH878f356xovI69nA89N+H3AncGOmbr0g5T+jTXnPAg5I+08F7gU+m8I7APsDy4AX5Cz33sApqT7tCWwGbgd2yaT5EvDX7crUyKcNar2KtH++mS1II44dganABR49tL8iHopuRwgNcXc3s9OA5Wb2ATM71d2vJhra9wBHFNEjS/mcQszNvN/M3gfcDfwzcAAxYuloJVuv6wW5++NEpf0s8DIz+0CS39HclZm9GXgFocReArzZY2RyE/D3xENcCO7+O3d/MAXHuvvvidHuA6nnej3Jxt6lfDezvyCU/H7uPkCY+fYGDjCzVwH/CFzn7r9qJasu/HuiR30JYS67C3jY3f+Yo7z19XsM8AczOwZYDvyxvhydZkF0KsYRHZTLgWuBvc1sbkqzFng0pema1Js/ChiTnunHgMeB55nZCwDc/SriHIdawbpd9gcmmtnO6V5+CTjazM509yfc/d+AnYFD2pR3FDDZwgnmv4F/AE43sznu/gd33wSMBfbJWe6ZxOr/O7v7r4nR1QHA7Eya76a82iOvJt+eNuAdwL8Cz03hpwE3ED2LL6WLP7OEfI8FfggcCKwghtML038XEb2y3QvI52jCVHMQYQa4l2hYx6VKdjXROBZ1Xn8H/C9wZhtp63vbf5/K9TZiBLgjsEP679XA00uuC1cR5phvAy/NKes56brXRlK10c7rgGsIRX1yims4oqDFSI8w1T4EnDaUnCGObyifaJQPSbJ/WCtnl9dhp1o+wBRiBLWKmLNP1EEAABS8SURBVE8aS8yvDgCfIzpFf1HQvRxHjL5vA55HdIQuJxTgSwnz8ybaHI20kd8kYCtwTArPJ0ZAbyEa+fXAUR3IeyrwH8BhGXmbiI7w89Nz/KICyr1HuvYLUvglwH8ClxFzcZs6uf+5L+T2shHmkC8TJoE+wnT0dkL7H0zYhKcWnOcYYuj8tvRQnEx8I+TNhINDTQFNy5mPpYf7lFT5TwLWEb2q2xhUQJMLPLfpwC+BubUyDJH2KQyavY5Jv58hevE3ZRqsN5BTEbRzrdLvt4j5gBNblb8NmacRXmPPI+YzvgrcX6tPNDGREfMsH8yEmyqgdG3OyNzvluVtVz5hSv0mcHrOun488FZC6XyImP/4i/TcvTyV+21EZ29BAfdyTGZ/EtGR+3J6nvcCPgz8mHCWmV1/TLd5pfB5xPdzXpDCC4m50TXA4a3yqr8XxBTAPQya9BYQHbuVpA5xN2VvkM/ZhFlwXgrPSXV1GXBQJ89CaQ9pr2/1F5Cw63+D6IXeRNik/x14d1l5A+MycbsQE6BTUngNMerquoefyWfHTNx4YCmDjfxVqQIfUvA5jgFm1coxVIVNjdLHiZ7ohhR3GPAjko0ZOIfoMe5fcr0Yk8nvhJzXfUIm7gZC4b+WsKm/H3h5o7qYOWYi8CDwj5m4ViOgthRPJ/KJjsnencpvIGcvYDXRKB+X4vYgJrb/hTAbjiFGQFcDJ+W5h2l/JjGy2S2F30aYcg8m3J8vJ3r2B3aZV3ZOpB94GTAxU4ceZVABvTKd5ylDyNs1s38YMcqfkMIXEKOc56XwmalOHTNUPWqSzx6Z/RMIxdyfwqcRI6BXp/BsYuQ15Oj8SXl0c0FH+5a9eEQvbA7hWfOUVGFqE3xnEm6mO3b7wDXLmxh93JAq1Nz00H0/NU7PJhThcwrI7zSit3cN8KYU9wVikvXYlE9bD16m7H0MMTHOk3tT49qQ/TnCVHFq7RjClHQvMZn+7+Q0i7Rb/pRmp8x+y/I3yOM0QslfRnJQqclJjcq9wBFNZOwP/FXa351Quh9v9/q2qqtdyLdGv11c/52IxncFoXx3ysS/LNXT3Qkl9XbgaV3ksQ+hYGqm5AeInvy/Mtixu4B4p+uQ9Kx9lHB978gphrCWvBd4UXqWfpbO4T7g4JTmNcQc8QsJK8drgRtp4MSS4j5AWCiOAjYQo+Sv1+oK4Ri0GTg0hV+Xzm1Cu/cllWMlMVrenzCnfoYY3bwlXbtT0rnMz5zH98gox5b5dFNJRvuWeYjOI+Zx/i41Boszac4nhrkHl5D/iYQ31cnEaGd5ij8+lWct8LIC8jk2yZpJmDnWpvhjCFPDauAvO5R5RnpwVwCXknrEdWlqZrJdhpBT/+C9MN2HL5DMEil+Yto69sQrs/xt3N8fEo3bdwmT25zUQPQDGxnChEX0yo8mNb7E6KChgsj8TiKcAVp6bJUtf6j7TXTk9iM6Q1emuL3S9XlqJn235q9nEPNFlxBmtdqo4zKic1dTQO9gcGQ+I5t3B3ntTSiLJUQHqSbvEsJkXDOHvZZBs96upJFRnawZhDPBWwhrxG2Z4xcTDio1BfQ20min9ox0UOZDGOx4fouYY6uV+2XARxhUQKeTlFz6f7eOrk83N3C0bsDTGRzCPi1d+JodcxIxofZmwgT2KQqagGxQjjcCzyUa3LXAvil+F6L382cTR858Xkp4hZ2WHoZp2UrE4MR3uz2m5xC2/93SA3Z3fYVksLHag5i/OqiBnOzIcx7Rm6sN+S8iGusDiPmBRQVe90LKP4T8MUQv/h+IieCT0/39EHAr0cDuOVS9YtDstzNha39Hpjzr2dZENi7z3yqSKatVGcuU3yC/sZn9mgIaQ8x/fZLogf+QNJme8/7Wzm3/dL2/UatXKf5S4hl/Uoeji7xq57Iv8C5iZP7Xmf9r72G1tF6k67uOaOynE6OZe4DzMmneRZgLX1Bfhg7KvCvhQLOAGP3MJMyu78+kmUt4Hr6jwXXtLL+8F3m0bISp5aPEi5C1RvcmMu/sEN4dS9L+DiWW5TxipPVtBnufpxA9jiLfXzmX6GV/g+RMQDSIH6XDXkw69oh0/eYB/wY8M8XXFHit4d6d8O4Z0gOH6Gl9jzBPfQp4XYq/iJh7+wHJvl3Q9Si0/Bm5tUY627vfk2iwa4r+XmIeo2EPu9GDTYxONhJu5rVy3Q9cm0kziTDNtLrWpcqvkzuNFo4MxNzR04heeFdza3Xyap6QO6ff2qob7wEmZdJ9oJNzaZFnrR0ZR5jtriCNcFL8pbTx/g0xgrqeUAqfIUYnbybW1pubSXcJyeGgy/KOTfX/I6k+TmRwlP7GTLqXUUDHu5CHdjRsRG/r1UQv9M0MvmW+NtN4vDk1emMbPawFluVpxLsvH0nho4le55wS8vocYTMeR6ya8BPSvEobx9abxp5F9L7WMjiPcSoxmsiacL4FvLCF7JcTZr+aUjybWMvunHRvJpPT+67M8mfu4/i0fyLwQWJSeSbRs/w2oewPIHr4h7Yh8wjC07L2UuQswt31/Ez5jkn74wjzYdsjkrLlp+NKdZSoO24ayRRI9Nq/Qrjpn0Eo01uIUcNTOpXdIK+9GTSF/UWqN3cQCnRSyueD7T5fdbKvJrwra/dhCjEncw0dmsZb5PMO4j2n7L15MTFye2tR+bhL+UAMY2vvVxgxuvk4g27MVxO94GuIXnihprZGDxSh3F5MTIR+kzCJvaTgfLMeP58jlN2qTh8MYkR2aa1iEqawfyKUxCsIV9XTM+nPIWOPbnYdCJPgH4E3pPA44CwyI6CCrkMh5W8gd4fU0FxDzKF9NzU+y1PcSUQj/y1C4Td8Z4Ww9Z+VKesG4P8RjfcrU/yhhIfYBQ3KMOSkfNny69J35ShBDkcG4uXc+wkz9leId2BeTZhUzyU6MV8nXJXbdhxpkte7iE7EyQyu9Xh4ur8XpOt1GdHBbWuOksHR9uuJOap/AY5OcX2ENWRp2i9ilY0XEmbu9xFzR7VXHI4nTH1Pp6gVVYoQ0qsb4b32BPGC3PlET6I2Anovg268R6QG5BkF5Zu1cw9Z4YnJ1j2zFaSDfHagbv6mUYVL+zsy6ALayhuqVlGfmx6stxOOEZ9P8WcRJoZrqHvprNH51pXjIJKnGWHjvof0EirhBn4mLTzR2rguhZa/WR6E88CHCG+jWqM7jcG1tyYQ3lfPzpYrI+MAosPzeqKnexcxOjse+AWhuF6f0vbT+QikVPkN8qvckSHJuIZY1PY9mbgDCdPp1HQdDu9Wfl1elxLKbHkm7unper6YGB11/ToAsIhwNjgyhfeihBeqiZHhRwlrT80JY1KheRRd6F7b0oP2RLrInyBGAJ8iJtX+hVBIudd0yuS3GzGncBCh0JYRvfr6hievM8FYogc2h/Cm+TwNvLMYnCxs2ZshlHXtxcfD0/Wan8LjiDmyz2fSZ98fauelxv9H9Bi/kB7iPQlTyXdJ7xTkvCallr/RtSUcC+4g5tVqPfmDUgPV9KVkwsV3I/CBTNwBqdzrkuz5xIuEr+y0nGXLH6KeVeLIUHsGMvtXA1vY9t25ZaRGvIC6lc3rnYQin1V75lJ9ntuBvPr2ICv/QsIZKtd6bUPdp7T/EqIdfBslzHEXKqxXN8IEsp4YKeybHrrbgYeJFxlzL12T8pmRKuYriTWk/pPkxtgkfa33N44uHA0Ic8ld6WFvahfO5LMjTUYVhKfd20nL6xCjwe8TpoDaxOpYwoa+Kit3iHz7U2M3iXgX4qsp/gvE3FpthHIm0QvvegHPMsrfQd4HE5PDVxOjt/1T3s8cop6sTQ3vIrZ1m/1LBheinEXMi3U0yVy2/Lq8KnNkGKpup/3PpnxnpWfjFxQ04knysw33+4nR9PnEvOFmWphraWERqZO/iC49AOnMInIGZXn1liG0FzfC3fg/GJzgnkRMGE8rSP5uZNajIibPf8ngekzZHtkYtjU7fJn2bcRWV3k+SJgXzia5bNelz7oOf3mo801pphC27YnEy5ADxAtmu2TK3rKxIkZka4lR4L7EW9R/S4x+bmdwkrj2tnbbL69VUf4h8mjmufVcwr33/vTb8D0eYmSwhnibf3eix7yE5OBAzFHeTHg9/ZjOTW2lyh8i39IdGbLXu9k9IbzGniBMSscWcW51eWUVxMWE63bLvGjfIpJrvoX2LSKlOVT9OY+yM+iljZhw/Q8K8HxpIv/dRM/6R8Qo4+XATxl8a38mmRfC0gP51W4ewlSJn5Yq9SHEJPp5qUI/mzCxjMnks6rZA1L3QB1H9OIXEXMWRxPKbUGjStxE3jFED/SwTNyBRC/7Gwx6iL2FmCTu6mXOssrfQP402lsD7RDgSpqsXJBJt1dm/9mEgvgAg3b+5xCOEd0u71Oq/CSjMkeG+mtO41FD9v+P0eXSPJnyDTVyyNa3S2gxuqJzi0jNLDmu0bm2Uf52LSK1fJpaRPJshQobDRsxx/A9CvLoSDJr5qNXAX8Absn890pCAV0E/BdpNVs6NDsQb25/KO0fQywbciPR65pCmLg+m8KP1OQSvd+BZvlkyj6dUGY7pQbrY8RIZQIxkfot2lxYlTB/XZD2axV8AqGcr0jX4jzi/YJcq4SXUf4GebS7BtoY2nTqqKXPlP1SwpTzgro0eTycSpFP9Y4M7Y4aGi4L1GFeHc2ltln2nraIdH3fihY4GjYKMPFkb3pm/1nEJN4VRO+7Zuo5JTW4x2cq2ec6eSiJOY1fEwrmUsKkNJMY+i8lecUQLswvyuRzKcl1cwjZc1Kj8RliMvoZqcH6CNGzm0BmIcJW14Jwf31fLS7TCE4ieoDXEUvpFLWEfSHlbyC3U9fhbd4E76TxS3ldTryjUtjq4kXLp3pHhkpHDenYbkYOew2RrictIrnrXBlCtf35RtcamdnpBtfeJZlOrM90DYNvXNe/09D2IoaZY3dJD/j6zH8HEz38z5FZIHSofIje6jPS/nOJ5WBqCuutRM9sKmHH/0c6/HIr0QO+k8Elc8ZkzuHNhAmu65Fn2eWvu7aVuQ6nejO9xPqaSz4VOjIkOZWMGtIxhY8c6GGLSCH1rSzB2v5882cTnk1HEz2zd6eKfCAxIrm+QcXu6k3u9LsLMcl5Xea/Q1K+7awjdSDx7stZxItr3ycmnp+VeaCWkN6ZoIv5MWKk8V6il51dW+ssYr20Jz3MHcguvfzpuMpdh0fyxvA5MlQ2akjHFzJyYBRYRHLfu+GutKN1Y/ADbZ9MlfQkYqj+9Eya/ckxp8Fgz+kI4l2YD6fwBGJtpqzraktTIjFx/iOSPT7FPY3ozS7KxJ1P5hv3XZZ9KrGu1teJkcn7UpnzXI/Sy0+DjgEVug6P5I0KHBnq7wMljhrSMYWPHOhRi0jh9WW4K+xo23iynf8S4i33NQy+yf4a0ofC8uZD9I5WpYf8vxjs0e9CvFuwtAOZrwU+mvbHEOaRuYRN/YH0UCwkJpNzf8KY6C2/kLSaBOk7STnkVVZ+KnQd7rWN8h0lKhk1pONKGTnQYxaRUurJcFfU0bIRPaTd037WrfN1wO8Y7P0dSvTOj+kyn30y+7sRb8+/LIVnEmvBXZHCE+hgspDo2f0r4c2zND2s64nJ/68T5qyvkuYFKNAjsKB7UFr5qdh1eLRsFOwoQUWjhrrjCxs50IMWkdLqxnBXztGyEasWP0LymiK9q5L2350qwfXpt6ted6q4l5D5pDUxqjqHwa8+Hkf0pLJLoLe75MouhFnk+8RSMy8iTBWHEvb7ZxNzMu8d7utdZfmp2HV4tG0U7ChBRaOGdNy0TN3KM5fasxaR0urFcBdgNG2ES+99pAX42PZTy8cDzyS5D+d4GIz40uNXUngesRZd7bO8BxErBLRczmOIPCbXhY8llv+31BD/GwV9ObSk+1BY+anYdVjbkPei9FFD9t5RwMiBUWARKe1+DncBRtuWehz3ZRtA4iXGq0hfSe1C5t6E2af2RdOdCffN2ue130UM/z9PmIP2T3En5jyX8cS6VD8ATsvE51p6vsJ7kav8VOw6rK3pfahk1JDNi4JGDowCi0hp93W4K9Zo3FLF/WnaPxj4FfDSLmUdSJiK7iTcVGufF9iNmAT9dAo/k/gEwbOIddLW0+X7K0neeAaXnnlJijPqzBkjdctbfobJdVjbNvegklFDklHmXOqosIgUfn+HuwCjdUsK6H+JCekzahWkQxkziLmG44jRz/nAJZn/JyQFdCuDpoL9iWVpDingHMYz+DGpEa9wii4/FboOa2t4/UsfNdTqBiWPHBhFFpHC7u9wF2A0b6lXU+s9ddP4vRB4IhPen/iQ1PPYdiJ0KZklRSj4o0/b+0bJrsPahrz2pY8aasdS/lxqz1tEitxqvWVRImZm3uWFNrM5xHItzzSzs4mX2X5BfM/950TP75vu/gczG+PuTxRWcPEkzGx/4vPLBixx9y3DXKRRj5mdQiydc1jtepvZi4n3rP7G3R/vUu7exPtY/+Pu95vZzsRI4TF3n2dm7yKUxBiiw3casXLG3e5+Z45z+SLwKPGJ+C932j6Y2QxCwbydGNG8lHDlX5z+n0B8KHEP4tMdnurtF4FXufs93ZS9aKR8egAzO5X4wNp6dz/MzCYDuxIv0V3n7t8b1gJuZ5jZdAB33zjcZdleSI32VakTdjDwNaLx/lKX8g4kGvBHCNP4l939JjPbjXhX6Al3f42ZPZMwba0nVtD4GNGg35fjXI4nTIlf7KZjamYvBL7h7mNSeH/iUx1/C/zG3X9mZruksn7c3deldJPc/ZFuy100Uj49Qqqwy919n+EuixDDQRGjhiRnRIwctneLiJRPD5Eq3HLCvXTE9GCEqIq8o4YkY3SMHHrcIiLl02OkCvdbdx8Y7rIIMVzkGTWk43t+5AC9bREZN9wFEJ3h7rdB/odPiF4mb9139zvM7E1mtpUYOTytbuSwxd3/kNKOSMUD4O5fM7MFZvYQPWYR0chHCLHd0ssjhyy9aBGR8hFCbNeMprnUXrKISPkIIbZ7enHk0OtI+QghRKKXRg69jpSPEEKIyhkz3AUQQgix/SHlI4QQonKkfIQQQlSOlI8QQojKkfIRQghROVI+QgghKuf/A6DKGclUAvvWAAAAAElFTkSuQmCC)
%% Cell type:code id: tags:
``` python
@interact(quality_top=(0., 1., 0.01))
def low_quality_tags(quality_top):
pd.DataFrame(pd.melt(all_quality).groupby("bodyparts").value.apply(
lambda y: sum(y<quality_top) / len(y) * 100)
).sort_values(by="value", ascending=False).plot.bar(rot=45)
plt.xlabel("body part")
plt.ylabel("Tags with quality under {} (%)".format(quality_top))
plt.tight_layout()
plt.legend([])
plt.show()
```
%%%% Output: display_data
%% Cell type:markdown id: tags:
# Generate coords
%% Cell type:code id: tags:
``` python
%%time
deepof_coords = deepof_main.get_coords(center="Center", polar=False, speed=0, align="Spine_1", align_inplace=True, propagate_labels=False)
#deepof_dists = deepof_main.get_distances(propagate_labels=False)
#deepof_angles = deepof_main.get_angles(propagate_labels=False)
```
%%%% Output: stream
CPU times: user 49.3 s, sys: 510 ms, total: 49.8 s
Wall time: 49.5 s
CPU times: user 54.3 s, sys: 475 ms, total: 54.8 s
Wall time: 54.8 s
%% Cell type:markdown id: tags:
# Visualization
%% Cell type:code id: tags:
``` python
%%time
tf.keras.backend.clear_session()
print("Preprocessing training set...")
deepof_train = deepof_coords.preprocess(
window_size=11,
window_step=11,
conv_filter=None,
scale="minmax",
shuffle=True,
test_videos=0,
)[0]
print("Loading pre-trained model...")
encoder, _, grouper, gmvaep, = deepof.models.SEQ_2_SEQ_GMVAE(
encoder, decoder, grouper, gmvaep, = deepof.models.SEQ_2_SEQ_GMVAE(
loss="ELBO",
number_of_components=10,
compile_model=True,
kl_warmup_epochs=20,
montecarlo_kl=333,
neuron_control=True,
montecarlo_kl=10,
encoding=6,
mmd_warmup_epochs=0,
mmd_warmup_epochs=20,
predictor=0,
phenotype_prediction=0,
).build(deepof_train.shape)[:4]
gmvaep.load_weights("../../Desktop/deepof-data/trained_weights/GMVAE_loss=ELBO_encoding=6_run_1_final_weights.h5")
gmvaep.load_weights("../../Desktop/GMVAE_loss=ELBO_encoding=6_run_0_final_weights.h5")
```
%%%% Output: stream
Preprocessing training set...
Loading pre-trained model...
(227686, 11, 16)
%%%% Output: error
---------------------------------------------------------------------------
ValueError Traceback (most recent call last)
<timed exec> in <module>
~/opt/anaconda3/envs/Machine_Learning/lib/python3.6/site-packages/tensorflow/python/keras/engine/training.py in load_weights(self, filepath, by_name, skip_mismatch, options)
2209 f, self.layers, skip_mismatch=skip_mismatch)
2210 else:
-> 2211 hdf5_format.load_weights_from_hdf5_group(f, self.layers)
2212
2213 def _updated_config(self):
~/opt/anaconda3/envs/Machine_Learning/lib/python3.6/site-packages/tensorflow/python/keras/saving/hdf5_format.py in load_weights_from_hdf5_group(f, layers)
704 str(len(symbolic_weights)) +
705 ' weights, but the saved weights have ' +
--> 706 str(len(weight_values)) + ' elements.')
707 weight_value_tuples += zip(symbolic_weights, weight_values)
708 K.batch_set_value(weight_value_tuples)
ValueError: Layer #12 (named "kl_divergence_layer" in the current model) was found to correspond to layer kl_divergence_layer in the save file. However the new layer kl_divergence_layer expects 3 weights, but the saved weights have 1 elements.
%% Cell type:code id: tags:
``` python
samples = 10000
all_clusters = grouper.predict(deepof_train[:samples])
all_encodings = encoder.predict(deepof_train[:samples])
```
%% Cell type:code id: tags:
``` python
video_key = np.random.choice(list(deepof_coords.keys()), 1)[0]
video_key
```
%%%% Output: execute_result
'Test 42_s12'
(227686, 11, 26)
CPU times: user 3.35 s, sys: 656 ms, total: 4.01 s
Wall time: 4.07 s
%% Cell type:code id: tags:
``` python
video_key = np.random.choice(list(deepof_coords.keys()), 1)[0]
video_input = deepof.data.table_dict({video_key:deepof_coords[video_key]}, typ="coords").preprocess(
window_size=11,
window_step=1,
conv_filter=None,
scale="minmax",
scale="standard",
shuffle=False,
test_videos=0,
)[0]
scaler = MinMaxScaler()
scaler.fit(deepof_coords[video_key])
scaler = StandardScaler()
scaler.fit(np.array(pd.concat(list(deepof_coords.values()))))
# Get reconstruction
video_pred = gmvaep.predict(video_input)[:, 6, :]
# Get encodings
video_clusters = grouper.predict(video_input)
video_encodings = encoder.predict(video_input)
# video_clusters = grouper.predict(video_input)
# video_encodings = encoder.predict(video_input)
scaled_video_pred = scaler.inverse_transform(video_pred)
scaled_video_input = scaler.inverse_transform(video_input[:, 6, :])
scaled_video_input = pd.DataFrame(scaled_video_input, columns=deepof_coords[video_key].columns)
scaled_video_pred = pd.DataFrame(scaled_video_pred, columns=deepof_coords[video_key].columns)
```
%% Cell type:code id: tags:
``` python
tf.reduce_mean( tf.keras.losses.mean_absolute_error(video_pred, video_input[:, 6, :]))
```
%%%% Output: execute_result
<tf.Tensor: shape=(), dtype=float64, numpy=0.6589403550173625>
<tf.Tensor: shape=(), dtype=float64, numpy=0.8289839462658529>
%% Cell type:code id: tags:
``` python
# Draft: function to produce a video with the animal in motion using cv2
import cv2
w = 400
h = 400
factor = 2.5
# Instantiate video
writer = cv2.VideoWriter()
writer.open(
"test_video.avi",
cv2.VideoWriter_fourcc(*"MJPG"),
24,
(int(w * factor), int(h * factor)),
True,
)
for frame in tqdm.tqdm(range(500)):
for frame in tqdm.tqdm(range(100)):
image = np.zeros((h, w, 3), np.uint8) + 30
for bpart in scaled_video_input.columns.levels[0]:
try:
pos = (
(-int(scaled_video_input[bpart].loc[frame, "x"]) + w // 2),
(-int(scaled_video_input[bpart].loc[frame, "y"]) + h // 2),
)
pos_pred = (
(-int(scaled_video_pred[bpart].loc[frame, "x"]) + w // 2),
(-int(scaled_video_pred[bpart].loc[frame, "y"]) + h // 2),
)
cv2.circle(image, pos, 2, (0, 0, 255), -1)
cv2.circle(image, pos_pred, 2, (0, 255, 0), -1)
except KeyError:
continue
# draw skeleton
def draw_line(start, end, df, col):
for bpart in end:
cv2.line(
image,
tuple(-df[start].loc[frame, :].astype(int) + w // 2),
tuple(-df[bpart].loc[frame, :].astype(int) + h // 2),
col,
1,
)
for df, col in zip([scaled_video_input, scaled_video_pred], [(0,0,255),(0,255,0)]):
draw_line("Nose", ["Left_ear", "Right_ear"], df, col)
draw_line("Spine_1", ["Left_ear", "Right_ear", "Left_fhip", "Right_fhip"], df, col)
#draw_line("Spine_2", ["Spine_1", "Tail_base", "Left_bhip", "Right_bhip"], df, col)
#draw_line("Tail_1", ["Tail_base", "Tail_2"], df, col)
#draw_line("Tail_tip", ["Tail_2"], df, col)
draw_line("Spine_2", ["Spine_1", "Tail_base", "Left_bhip", "Right_bhip"], df, col)
draw_line("Tail_1", ["Tail_base", "Tail_2"], df, col)
draw_line("Tail_tip", ["Tail_2"], df, col)
image = cv2.resize(image, (0, 0), fx=factor, fy=factor)
writer.write(image)
writer.release()
cv2.destroyAllWindows()
```
%%%% Output: display_data
%%%% Output: stream
%% Cell type:code id: tags:
``` python
samples = 100000
all_clusters = grouper.predict(deepof_train[:samples])
all_encodings = encoder.predict(deepof_train[:samples])
```
%% Cell type:code id: tags:
``` python
all_clusters
```
%% Cell type:code id: tags:
``` python
sns.scatterplot(all_encodings[:,0], all_encodings[:,1], hue=np.argmax(all_clusters,axis=1))
```
%% Cell type:code id: tags:
``` python
# import umap
# from sklearn.manifold import TSNE
# from sklearn.decomposition import PCA
# from sklearn.discriminant_analysis import LinearDiscriminantAnalysis
# import plotly.express as px
```
%% Cell type:code id: tags:
``` python
def plot_encodings(data, samples, n, clusters, threshold, highlight=None):
reducer = LinearDiscriminantAnalysis(n_components=n)
clusters = clusters[:samples, :]
# filter = np.max(np.mean(clusters, axis=0), axis=1) > threshold
clusters = np.argmax(clusters, axis=1)#[filter]
rep = reducer.fit_transform(data[:samples], clusters)
if n == 2:
df = pd.DataFrame({"encoding-1":rep[:,0],"encoding-2":rep[:,1],"clusters":["A"+str(i) for i in clusters]})
enc = px.scatter(data_frame=df, x="encoding-1", y="encoding-2",
color="clusters", width=600, height=600,
color_discrete_sequence=px.colors.qualitative.T10)
#if highlight:
# ig.add_trace(go.Scatter(x=, y=)
elif n == 3:
df3d = pd.DataFrame({"encoding-1":rep[:,0],"encoding-2":rep[:,1],"encoding-3":rep[:,2],
"clusters":["A"+str(i) for i in clusters]})
enc = px.scatter_3d(data_frame=df3d, x="encoding-1", y="encoding-2", z="encoding-3",
color="clusters", width=600, height=600,
color_discrete_sequence=px.colors.qualitative.T10)
return enc
```
%% Cell type:code id: tags:
``` python
plot_encodings(all_encodings, 10000, 2, all_clusters, 1, 10)
```
%%%% Output: error
---------------------------------------------------------------------------
NameError Traceback (most recent call last)
<ipython-input-53-299d75cbfe33> in <module>
----> 1 plot_encodings(all_encodings, 10000, 2, all_clusters, 1, 10)
<ipython-input-52-f274dd41a1f9> in plot_encodings(data, samples, n, clusters, threshold, highlight)
1 def plot_encodings(data, samples, n, clusters, threshold, highlight=None):
2
----> 3 reducer = LinearDiscriminantAnalysis(n_components=n)
4 clusters = clusters[:samples, :]
5
NameError: name 'LinearDiscriminantAnalysis' is not defined
%% Cell type:markdown id: tags:
# Preprocessing
%% Cell type:code id: tags: