{ "cells": [ { "cell_type": "markdown", "id": "aaa3548c-c138-467b-90a2-efa7134c2747", "metadata": {}, "source": [ "# WPLSC (within-participants PLSC)\n", "\n", "Regular PLSC captures multivariate patterns of association when the unit of observation is a participant (or a participant average for a within-participants condition). However, we often want to examine associations within participants, where the unit of observation is a single trial. To do this we can use a within-participants variant of PLSC ([Roberts et al., 2016](https://doi.org/10.1016/j.neuroimage.2016.04.028), which is implemented in `pyplsc` in the `WPLSC` class. Briefly, the method is to compute cross-correlation matrices within participants, average these matrices across participants, and apply singular value decomposition to the average matrices." ] }, { "cell_type": "markdown", "id": "a7f2dcf2-88c0-4160-ad2c-9d31f53151eb", "metadata": {}, "source": [ "## Setting up simulated data\n", "\n", "We'll simulate data with 10 participants, in which a single latent variable explains the relatinoship between a set of 50 observed variables in the data array and a set of 4 covariates. The latent variable will have a sinusoidal pattern of loadings onto the observed variables, and positive loadings onto all of the covariates. First, we will import the necessary libraries and seed the random number generator for reproducibility:" ] }, { "cell_type": "code", "execution_count": 1, "id": "7e051aea-75b1-497c-9b9a-b014de2be5a4", "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import pandas as pd\n", "from pyplsc import WPLSC\n", "from matplotlib import pyplot as plt\n", "\n", "np.random.seed(123)" ] }, { "cell_type": "markdown", "id": "3b7750cd-44be-4f09-a964-cd40f01f204f", "metadata": {}, "source": [ "Next, we will create a set of sinusoidal loadings from the latent variable onto the 50 observed variables. Importantly, these loadings are shared by all the participants---between-participants variation in loadings is not explicitly modeled:" ] }, { "cell_type": "code", "execution_count": 2, "id": "ffe99140-1f85-455e-b0c8-5f974ad45caf", "metadata": {}, "outputs": [], "source": [ "n_var = 50\n", "data_loadings = np.sin(0.2*np.arange(n_var)).reshape(1, n_var) # Sinusoidal loading pattern" ] }, { "cell_type": "markdown", "id": "1874b3d4-2740-4f44-a14b-f4e88d00acf2", "metadata": {}, "source": [ "Next, for each participant, we'll simulate a normally distributed latent variable, compute the data array as the latent variable multiplied by the loadings plus noise, and compute each covariate as simply## Fitting and evaluating the model\n", "\n", "Since the data are not stratified by any experimental condition, we can fit the model by simply providing the data array and the covariates withut a design matrix specifying conditions: the latent variable plus noise." ] }, { "cell_type": "code", "execution_count": 3, "id": "26b24799-bb3e-4bb4-96f9-4117f1ef51ef", "metadata": {}, "outputs": [], "source": [ "n_subj = 10\n", "n_trials = 40\n", "n_cov = 4\n", "subjwise_data = []\n", "subjwise_covs = []\n", "for subj in range(n_subj):\n", " # Simulate latent variable\n", " latent_var = np.random.normal(size=(n_trials, 1))\n", " # Simulate data array\n", " data_noise = 0.5*np.random.normal(size=(n_trials, n_var))\n", " data = latent_var @ data_loadings + data_noise\n", " subjwise_data.append(data)\n", " # Simulate covariates\n", " cov_noise = np.random.normal(size=(n_trials, n_cov))\n", " covariates = latent_var + cov_noise # Effectively all loadings of 1\n", " covariates = pd.DataFrame(covariates, columns=['Cov. %s' % cov for cov in range(1, n_cov + 1)])\n", " subjwise_covs.append(covariates)" ] }, { "cell_type": "markdown", "id": "fc15e13e-1503-4a48-8e70-86ca6dccf7ce", "metadata": {}, "source": [ "## Fitting and evaluating the model\n", "\n", "To fit the model, we simply provide the lists of participant-wise data arrays and covariate arrays to `WPLSC.fit()`:" ] }, { "cell_type": "code", "execution_count": 4, "id": "ff95fb13-c2b0-4e83-aece-03e5d53c3861", "metadata": {}, "outputs": [], "source": [ "mod = WPLSC()\n", "mod.fit(data=subjwise_data, covariates=subjwise_covs)" ] }, { "cell_type": "markdown", "id": "4e3d9d46-7162-449e-adc7-8fb4d3ca7a7f", "metadata": {}, "source": [ "We can use permutation testing to evaluate the significance of the latent variables identified by the model. This yields one $p$ value per singular value:" ] }, { "cell_type": "code", "execution_count": 5, "id": "e5291189-f1a8-4dde-a559-c268e6c89e91", "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "Getting permutations: 100%|███████████████████████████████████████████████████████████| 10/10 [00:00<00:00, 117.39it/s]\n", "Permuting: 100%|██████████████████████████████████████████████████████████████████| 1000/1000 [00:02<00:00, 355.13it/s]" ] }, { "name": "stdout", "output_type": "stream", "text": [ "[0.000999 0.82317682 0.44355644 0.10689311]\n", "[0]\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "\n" ] } ], "source": [ "mod.permute(1000, return_null_dist=False)\n", "print(mod.pvals_)\n", "is_sig = mod.pvals_ < 0.05\n", "sig_lvs = np.where(is_sig)[0]\n", "print(sig_lvs)" ] }, { "cell_type": "markdown", "id": "e8a4a105-3637-4833-ad35-e8461d98c40e", "metadata": {}, "source": [ "Thus there is one significant pair of latent variables, in line with how we simulated the data. Next, we can perform bootstrap resampling to assess the reliability of the data saliences and to evaluate how reliably the design scores correlate with each covariate." ] }, { "cell_type": "code", "execution_count": 7, "id": "76a685cb-eb1f-40e9-a0cb-6a4ead9bfe19", "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "Getting participant resamples: 100%|████████████████████████████████████████████| 1000/1000 [00:00<00:00, 18619.42it/s]\n", "Getting trial resamples: 100%|█████████████████████████████████████████████████████████| 10/10 [00:02<00:00, 4.41it/s]\n", "Resampling: 100%|█████████████████████████████████████████████████████████████████| 1000/1000 [00:07<00:00, 139.28it/s]\n" ] } ], "source": [ "mod.bootstrap(1000, return_boot_stat_dist=False)" ] }, { "cell_type": "code", "execution_count": 11, "id": "c1f42f2a-3e7c-4e80-ac72-8a8883a2f690", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAosAAAHrCAYAAACn9tfQAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjgsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvwVt1zgAAAAlwSFlzAAAPYQAAD2EBqD+naQAAkPxJREFUeJzt3Qd8U1X7B/Cne7dAJ7tsZCMqgoqKKIoLJ+IAFcGFIo5XEAXRP6KoCCqK4n5duLcMURAFxBdE9t6ldFBoS/fI//M7yQ1Jm7RNSZqbm9/38wkkN+tk9Oa555znOQEmk8kkREREREQOBDraSERERETEYJGIiIiIasSeRSIiIiJyisEiERERETnFYJGIiIiInGKwSEREREROMVgkIiIiIqcYLBIRERGRUwwWiYiIiMgpBotkWEOGDJHRo0d7uxkkInv37pWAgAB57733+H7UU2pqqtx6660+997Xt93kPnPnzpVWrVpJSUkJ31aqFwaL5DL86ODH53//+5/d9q+++kptf+utt5zed/Hixeo2L7/8co3PgZ3ao48+Ks2aNZOIiAjp27evum9d/fnnn7Jo0SL1GLZ27twp1157rTRu3FgiIyPl7LPPlt9++83hY7z66qtyyimnSFhYmDRv3lwefPBBKSgoqNPzFxcXy/Tp06VLly7qeXD/6667TjZt2mR3u99//12uuOIKadmypYSHh0tKSopcfPHFqv1VvfHGG9KmTRtp0qSJ3HLLLZKXl2d3fWVlpfTu3VueeeaZOrWRiFyzZcsW9fcZHR1t/TvMysqq8/2/++47OfXUU9XfOoK3KVOmSHl5ebXbHTt2TMaMGSOJiYkSFRUl559/vqxdu7ba7ebPny8333yzdOjQQe1XzzvvPIfPi2C9tLRU7UOI6gVrQxO54t1338V64qa///7bbntxcbEpLi7OdP755zu976233moKCgoyZWRk1PgcN9xwgyk4ONj08MMPm9544w1Tv3791OXly5fXqY1XXnml6aKLLrLbtn//flNCQoIpOTnZNG3aNNOsWbNMPXv2VI+7bNkyu9v+5z//Ua/x2muvNb3++uum++67T92u6mM6c/XVV6vb33333aZ58+aZpk6dakpKSjLFxMSY9u7da70drkNb/+///s/01ltvmZ5//nnVpsDAQNPPP/9svR1ed0BAgGncuHGm2bNnm1JSUkxjxoyxe865c+ea2rRpoz4HvamsrDQVFRWZysvLvd0Un9W6dWvTyJEjfe69r2+79ebAgQNq/9GuXTv1N4h9SOPGjdXfa0lJSa33/+mnn9TfMPaPb775ptqn4O/8rrvusrtdRUWFqX///qaoqCjTk08+aXr11VdNXbp0UfuO7du329323HPPNUVHR6vHRFtw2Rns0/BZ4PtA5CoGi+S2YBFGjRqldoBpaWnVrsMPFoLJiy++uMbH/+uvv9TjI3CyvS920ggaa4NAFIEagi9b99xzj9q+detW67aCggJTy5YtTaeeeqp126FDh9TtbrnlFrv7v/LKK6pd3333XY3Pf/DgQXU7BLq2fv31V7V95syZNd4fbUJAO3jwYOu2Rx991C4Ix2eAgFFz9OhR9UP25ZdfmvSkrKysTj+kZNygy1fbXRUO/CIiIkz79u2zblu8eLH6m8YBbW0Q8CGwxN+EZtKkSSqA3LJli3Xb/Pnz1WN+/vnn1m2ZmZmmRo0amYYPH17tABjBJXTt2rXGYPF///ufetwlS5a48KqJzDgMTW6FIREMh3766afVrvvxxx8lNzdXbrrpphof44svvpCgoCA1DKPBsM2oUaNk5cqVcuDAgRrvj+fB0M6gQYPsti9fvlwN03bq1Mm6DUPEGAbGEM+OHTvUNjwH7n/DDTfY3V+77Oi12crPz1f/Jycn221v2rSp+h/D6jVBmzD8hKEoTVFRkRo612AIrLCw0Hr5ySeflO7du8vVV18tdVFWVqYe47bbbqt2HYa38X4//PDD6jKGryZPnix9+vSRuLg4NSx2zjnnVBu+1+bGvfDCCzJr1ixp166dGsLfvHmzw3lz69evV8Njbdu2tQ7B33777XLkyBG7x8Vrw30xhQC3b9SokWoH2m77Hmg+/PBDOeOMM9T7iPdswIABakqCrZ9//lm9BryWmJgYufTSS6tNEXAGn8sDDzygpg7g9bVv316ee+459b0HHIRj2BCfYWZmpvV+eB/xGeF90aYzaK9t69atcv3110tsbKzEx8fLuHHj1FSGmuTk5KjPCI+JYVHc95JLLpF///3X4edi+97jfcR90tLSZOjQoeo82ovHq6iosLs/Xhc+z65du6rPCd/rO++8U44ePWp3O7zu//u//5MWLVqo9x7vQV3fU2cwrNqtWzdZs2aN9O/fX/3tYCoG5uA1tC+//FIuu+wyNXyswT6mY8eO8tlnn9V4X/wN4IR9WnBwsHX7Pffco9437PM0OI/32PZvGZ8Nvh/ffvut3bxDfAcDA+v2M46/X/zN4zGIXMVgkdwKP8z4sfj444+rXYdt+BHBj1NN/vnnH7UDxo+fLQQAsG7duhrvv2LFCvWD27p1a7vt2Mk6CtTQJsAPknY7qHrbqrdzBsEA3oMXX3xRvv/+ezl48KCsXr1a7rrrLvVDVzUI1QK07OxsFTQ89thjsnHjRrngggus159++umyYMECFfQgqMVja+8HfoTw44kf9LoKCQmRq666Sr755hsVxNjCNrwHWjvRNsxDxQ83giIEOJinNXjwYIefxbvvviuvvPKK+mFEO/ED5QjmoO7evVsFfbg9ng+BOBKT8ANaFX4sEYhjLijOI/iZOnWq3W1wGfPI8PqeeuopdRk/qL/++qv1Nv/9739VcIgACa/niSeeUO8h5q8isKoJgtNzzz1XBaQjRoxQc2/POussmThxoprTCgjM3nnnHRXs4TPXYH4agie8PwhSq742bZ4rXj8e1/ZgyRG8d/isEMDMnDlTHnnkEdmwYYNq36FDh6Q2CArxGeJvBQE+7ofP680337S7HQJDPDZe5+zZs9Xn9dFHH6n74qBDgwMKvJc9e/aU559/Xh0EXHTRRXWe5+sMglK8Jwh2ZsyYof627r77bvUe1wYHp/i7qu10/PjxGh8HQTUC/9NOO63adfg7xD6rJtr1Ve+POdl4Pbb3x3nMa6waBOJ58P3bvn271Bce19F8aKJaWXoYidwyDA2PPPKIun7btm3Wbbm5uabw8PBqwyiOYDhl4MCB1bZv2rRJPS7m5tXk7LPPNvXp06fa9ssvv1wN5eTl5dltx9A2HveFF15Ql9esWaMuP/3003a3W7BggdqOOUK1wVA6hs1xe+2ENqWnpzu8PYactduFhoaa7rzzTjX0rsF8M8yD1G6DofP169er6zCPsuq8p7pYuHCheqzvv//ebvuQIUNMbdu2tXvuqkPJGPbGUPntt99u3bZnzx71eLGxsWrYzJZ2Hb47msLCwmpt+uSTT9Ttfv/9d+u2KVOmqG22zwVXXXWVKT4+3np5x44dagoEtmtDcxptnlZ+fr76DowePdru+sOHD6spElW3V4XvBOaSVZ07NmHCBDUXF8OCGgxNot0ffvihadWqVer6Bx54wO5+2mu74oorqk2ZwPZ///3X6XAu5qZWfZ14n8PCwkxPPfVUje89HgfbbG8HvXv3tvvbwVxZ3O6jjz5y+Legbcfnje/tpZdeajcn7rHHHlO3q+8wNIZVcf8XX3zRug3fxV69eqk5wKWlpXW6f22n2tqHfR1u98EHHzjd39U0VxhTanAb2++H5vTTTzedeeaZ1sv4flX9rsOPP/6oHgPvvSO1DUMD5jljKJ3IVexZJI8MRYNt7yKGcNBzUtsQtDbkiuG9qjAEpl1fEwxj2g7ZatAbgSHEYcOGqaN3HKFjOFHL6tYeF0ffyL5GrxN6gdDbhGFL9LCgx6q25wc8f69evWTChAmq9wc9N3gcZEQ7Gl589tlnVa/h22+/LWeeeabq7bPNksSwPN5D9CqivWg7hh+RXYley6efflr1flx++eWqtwL/19a7NHDgQElISFAZlba9OOjxw3tk+9yhoaHWIUkMf6Jt6CVxlKF5zTXXqGGz2tj23OI9QQ8PXjs4elzbXjrAMDI+ay0rHO8z2ocerqq9MujtA7w2fAeGDx9u17OE14jP3FlmvObzzz9Xz4vP1/b+GI5ETx2y2zXoGUTv23333ad6O9Hj7CxT/d5777W7jPvATz/95LQt+BvRXieeG+8FeksxzcLR++eIo/cUPZa2rxdD/hdeeKHd60UvH55Le79++eUX9Z1Fu7X3GvD3dbIwbIu/PQ2+i7iMnr7aevnRU4rPvLbTf/7znxofR/ubr+9+qbb72973ZPd/NcH3Fvd3NH2DqCYnJk8QuUmPHj3UPKNPPvlEDVlqgSMCE/x41iWIcFQPTAuyapvzB46GMTGfC8OdCOAQEALmm02bNk39WODHT4PADAET5tABggkMMy5btky2bdtW69AXfnQxdPfQQw9ZtyO4wlAuAlAErrYQWNoG22gf5pXZzmXS2qvBjzMeH8ObeG/xnJgXiaFvBJ833nijLF26tMYfYQR2+GzwfuMHCuWPMLRoGyzC+++/r354MUxuO/SIYfWqHG1zBEEnhokx9Gw7t097D6uynSsG2gEBAlxMWdi1a5cKnlCuyBltXioCZUeqTn1wdH/MtXQWDFd9HQj+ESTifpge4ey7i9IntnAfvJaahsURGGNY+LXXXpM9e/bYzTXE0HJtEHxUfR14T23nIqLd+CySkpJqfL379u1z+Drw+I4O3FyBg5+qw/aYpgJ4f7QDDEcQ1LqD9rnVd79U2/1t7+uO/V9t+0XbgJ6oLhgskkcg4EFQhl4wzMlBDwR6A2wndzuDgAe9ZFWlp6dbfzxqgh/KqpPvNWPHjlVzrvCDjx4KBGn4Qbf9AQLURfzjjz/Uj+Xhw4fVjyASMPDctrdzBIFmRkaGSpyxhTlhCEYwZ6hqsGgL7cJ9EfChF8DZj8NLL72k3k+8JiT9oL0IGlAEGXO7MGcM8yXx/juDeYKovYaeU8wlxUT9zp07q3lnGszPQ+CK6xEAI3BA8Iz5dQjQqqrrjxnm6SGAwmPic0CwjgAIdey0ZBFbeM66Hhg4oz0u5i3i86yqtu8n7o9eNmc9UVW/GwjWtR9+zCfs169fndpZlx9z9FJijiAOaNCzjLmhCDDRm+fo/avr+2kLj4PPG3MUHalLD7I34YCk6pxcR/CdRQ+qM1pymrYPsoVteO8d9QY6uj/m0Fa9vzb/WLuts+epy/6vJtgvYu71yQSc5J8YLJJHYJgPk/7Ra4VEE/R61GUIGhA4ILjE8KJtT89ff/1lvb4mCHYQsDmDXgrbH20MoWHniQn8VSFI1HpLkASBHXZtq1EgUISqWaUIarDNURHeqhAk4vZI6HC0Y0c7kHmKYUIEONqQs/ZDov2PoLumYBEJSfhxwlA0EjyQCDJp0iS726B3E4GnVnRdgx7Nk/nRWrJkiepZxLBx1Z6/+kBvHIIbfE7OviO4DSAAqpotX9fnQDJEXe6LzwjDskjywAEAMo3Rs1418Up73bY9ssj8xmtB4O8MPhdkHGsHOxoMs6On2R3wevH3gb+NmgIM7TXhdeC7okEilLMDt7rCdxtJMra9i1qSR03vDyCjGKMBtRk5cmSNK9zg4BGBcdWFCADTQGrbJ2nX4/62gSFeGw7obJOZcFtUbsDnbzudAvs/BHq1HazWBAeTWGiAyFWcs0gegSFDDIsiCEHPFH4IUfqiLrDCCoIq26xM9M5g+BbzyqoemVeFQBA/ULZzr5xBzxaCIJTlqalnATtu9CZhZ207zwtDshiate0J0HbmVUvsYH4hfvRQvsfZsKX2Y49gF6/T2fAfem0R6KEXzrZMD9qirTQBjnrPbOHHCO83hq7R24ZAtuoQtNYDZduDhx8ulBiqL0ePCa5kdFeFnk+8HmRBV+1Z054HwRoOQNArZzucrqltNQ70huJ1L1y40OHnZnsggKUm0Q4Ec/guI6jH98xRT+icOXPsLmO6hDZ1oqb3sOpj4eDBUa98feH14m8RPZdV4bVq5Z0QPGM+L9pt26aT+Txtn8d25RFtJRIEb7UNM7trziJgysYPP/xgV7oLBzwIXDEXuaZ9AsoO4SAW3wPbg8jXX39dHYDhb1CD8zjgxH5Jg3mi+GwxF7mmHszaYC5rXffDRLbYs0j1htIVKOdSFWrEoXYdhqJxxIyj56q9VTVBQIidL3omEUxhnh7mzGF+UtVeFEdQFgU/zOgRsT1ix7wq/PhhiBdBFMqYoOQM5lhWTTzQ6tzhKB87f/SQogcB7bCdO4cfZhyp2/ZMYIeOHwcELXhOzKlCTxGWD0QvHgIGDYIB9PzhNSMw3L9/vwqK8Z7ZJp7YQjtwHYbSNehhwZxI9Hri8VHqBo/pqBerKgSH+JFHTyGSZqr2PKA0C364UGoH7y16J/C+YW5gbSVHnEHAhmAXw+V4f9FzgwQfPHZ94XuC7xkCGxyooFcJP6x///236mnFsDmeFz/QSDjBvFAMwyPowPuO+pzoQcPn5AyGzBH04z3Be41gBQcAGGJGTx++o+jVw2eIx8N3QuvZxXuMvwk8P+rr2cLrxvcSwT+CURxgYc6p7XSAqtAGfMcwrQIBANqA4WLbnr2ThakTmD6C9w5lktBLiqAQPYgIXjBnEsGNVqMRt0O7UOoGSWSY3uCol1PrEaytVBHgs0OyGW6LAzF899EWBF5oS03cNWcRUNIKrxm9udg/4LuPEkH4m7GtV+ponwC4LT5jvIf43qE8Fr5rd9xxh93fHN5P7DPwmOglx/uHeakIMquWikJClZZUhQMdfBcx4gD4+8JJg2QgDMtfeeWVbntPyI+4nD9Nfk8rnePshGWxICcnR5XxwLbNmze79L6hbAxWQMEqJXgMlJdwVjLCEZQiueCCC+y2oT1YWg+PiTIfWBoPK6NULaWjvUastoAyFlhmC4+FFViq0sqSVC29gecaP368qWPHjqr9WF0FSxju3r3b7nZYygulfnA9Vo1JTExUJX5sS8fYQlmSvn37mh588MFq1+3cudM0YMAAVdoH/+/atatO7xUeE6V48Dqw7KCj65955hlVugWvBeVVfvjhB/Wasa3qe2G78k5N5Vuw0g3K3KCUDcrWXHfddWr1HNwOJWWqlpfJyspy+D3EY9t65513VBvRVm0JNKy0Yeu3335T5YrwvCjphDJHWIoSq1zUBuV3Jk6caGrfvr36HuGzw/JsKL2EUi74/uNx8TlWhdeL75T2PdBeG/4+sLQkvmto89ixY+1KJzkrnfPQQw+ZmjZtqsqhnHXWWaaVK1eq12tbQsVZ6Ry0oyqtPVVheTqU1MHzoI3du3dXy8fh89KgjA+WtdTac95555k2btzocAUXvGe25WKcwetASRh8Lihxhc8Kj4e/G2/A60GpqsjISPW9vemmm1TZpbrsE+Drr79WZX/w3WzRooXp8ccfd1j+B/sPrIaF0lB4LrwPjkqVaZ+Xo5Pt3xBgX9eqVSsu90f1EoB/vB2wErkb5vwg8xjDQVUzNIn0AtUC0FuEXiF3zTPUO/SWoecdQ7roqa4J/oYxBIteOKo/TONBby6mr6BXlMhVnLNIhoRhSAz3YJiTiPQDyWuYV1xboEjug2kRGLKvWleTqK7Ys0hE5CX+2LPoCvYsEukDexaJiIiIyCn2LBIRERGRU+xZJCIiIiKnGCwSERERkVMMFomIiIjIKQaLREREROQUg0UiIiIicorBIhERERE5FSx+prKyUg4dOiQxMTESEBDg7eYQUQPC6qb5+fnSrFkzCQz072Nl7guJ/JfJxX2h3wWLCBRbtmzp7WYQkRcdOHBAWrRo4defAfeFRHSgjvtCvwsW0aOovUGxsbHebg4RNaC8vDx1sKjtB/wZ94VE/ivPxX2h3wWL2tAzAkUGi0T+iVNQuC8kIqnzvtC/J+0QERERUY0YLBIRERGRUwwWiYiIiMgpBotERERE5BSDRSIiIiJyisEiERERETnFYJGIiIiInGKwSEREREROMVgkIiIiIqcYLBIRERGRUwwWiYiIiMgpBotERERE5BSDRSIiIiJyisEiERER+b2S8gq/fw+cYbBIREREfm17Rr70eHKRTPtxs7eboksMFomIiMiv/bUnR0rKK2XBpsPeboouMVgkIo8oKCiQgIAAdcJ5IiK9Sj9WpP4/kFMkuUVl3m6O7jBYJCIiIr+WnltsPb81Pc+rbdEjBotERETk1w5ZehZhM4PFahgsEhERkV+z7VncwmBRf8HinDlzJDU1VcLDw6Vv376yevXqGm8/a9Ys6dSpk0REREjLli1l/PjxUlx84kMmIiIiqqvKSpMctgkW2bOos2Bx/vz58uCDD8qUKVNk7dq10rNnTxk8eLBkZmY6vP3HH38sEyZMULffsmWLvP322+oxHnvssQZvO5ErmOxBRKRPRwpKpbSi0np5++HjUmZzmbwcLM6cOVNGjx4tt912m3Tp0kXmzp0rkZGR8s477zi8/YoVK+Sss86SG2+8UfVGXnTRRTJ8+PBaeyOJiIiIHEnPNc9XTIwJk5iwYBU47so6zjdLD8FiaWmprFmzRgYNGnSiMYGB6vLKlSsd3qd///7qPlpwuHv3bvnpp59kyJAhTp+npKRE8vLy7E5EREREcOiYeQi6WaMI6dw0Rp3nvEWdBIvZ2dlSUVEhycnJdttx+fBhx0Ux0aP41FNPydlnny0hISHSrl07Oe+882ochp4+fbrExcVZT5jnSERERGTbs9gsLly6NI1V5zcfYseSrhJcXLF06VJ55pln5LXXXlNzHL/66iv58ccf5emnn3Z6n4kTJ0pubq71dODAgQZtMxEREek/E7ppXIR0aWYJFpkRbSdYvCQhIUGCgoIkIyPDbjsup6SkOLzPE088Ibfccovccccd6nL37t1V4sCYMWNk0qRJahi7qrCwMHUiIiLfVlpeKYEBIsFBPtXPQT5SY7FZI/Qsxll7Fk0mk1qBirzYsxgaGip9+vSRJUuWWLdVVlaqy/369XN4n8LCwmoBIQJOwIdKRETGlJlfLH3+b7Hc8cH/uL8nj/UsdkiOlqDAADlaWCYZeSV8py28eniGsjnz5s2T999/X5XCufvuu1VPIbKjYcSIEWoYWXP55ZfL66+/Lp9++qns2bNHFi9erHobsV0LGn0Jy6kQEdXNX7tzJL+4XJZuy5Jftzour0Z0MutCN20ULuEhQdIuMUpd3pyeyzfU28PQMGzYMMnKypLJkyerpJZevXrJggULrEkv+/fvt+tJfPzxx1WXMP5PS0uTxMREFShOmzbNi6+CiIg8bUdGvvX8i4u2y/mdkiQQY9JEJ6Gi0iQZ+eYexGZxEep/JLlszziuhqIHdrZPwvVXXg0WYezYserkLKHFVnBwsCrIjRMREfkP/HhrkHzw88bDcmmPpl5tExljegMCxuDAAFVnEZDk8s26Q0xyscFZwtTgOPxORK7anmnuWTytdWP1/8zF29SPPJE7aiwmx4aruYqgJblsST/Rm+3vGCwSEZGulZRXyL4jher89Ku7S6PIENmVVSDfrkvzdtPIIDUWm8aFW7edYinMvfdIgRwvKfda2/SEwSIREena7qwC1YuIpdjaJ0XLnQPaqe2zftnBNXzppKRbehabNjLPV4T46DBJjg0TFFnZdpjFuYHBIhER6dp2S3ILypogyXFk/9aSEB0m+3MK5fP/HbTeLq+4TL5Yc1BmLt6uzhPV5pDN6i22uJKLzhJciIiIarIz05zc0jHZPDwYGRos957fTqZ+v1le+XWHRIcHy4/rD8lv27JU4W74bWum/HfUGdIoMpRvLtXes1g1WGwWq75Pm6vMWzxyvEQO5xVL12bmeY3+gj2LRETkIz2L5mARbuzbSvUGoaDy/Z/8Iws3ZahAEcPUjSNDZENargyf95f6cSeqdc6izTA0WFdysVn277etmTLwxWVy6ct/yLbD/pX8wmCRiIh0bYelbE7H5GjrtrDgIHn0ks7qfIvGEXL3ee3k53HnyOLxA+TTMf3UMPWW9DwZPm+VKo9C/mt31nFZtOmww+vSLD2LzasEi1qSy9b0PHUQMnPRNrntvb8lt8g8veHPndniTzgMTUREulVcVqGyUm2HoTVX9mou53dOUokvtmv4dkqJkU/HnCk3vbVK1We84c1V8vEdZ0pKlaFG8g8PzF8n6w/myiejz5R+7eLtsuyzLT3PVYehW8dHSWRokBSWVsg1r69QPdWQGh8pe48Uyj8Hjok/Yc8iERHpOhMa5RRjw4MlyVI02VZseIhdoKjBcPT8Mf3UUDUeA4GjNp+R/EdZRaVaiQWWbrdfJjIj1xwohgUHSpMo+7mtqLnYOcV8cLIhLVciQoJk1rBeMu2q7mrb2n1HxZ8wWCQiIt3aYSnGjV5FR0FhTVITomT+nf0kJjxY1WW0nX9G/mHfkQIptxRvX7nriMNMaPQqOvpu9W5lLgDfNjFKvh17lgzt3Vx6tmwkqN2ddqxIMvP8Z3oDg0UiIvKp5BZXtGwSKX0sq75sOOhfQ4d0IpNe6yHMLSxzUJDbfr6i5v4LOsgrw3vLd2PPtk6BiA4Ltp5fu99/vk8MFomISPdrQtsmt7iqR3NzZivmrZH/Bososr1qz5FqS/01beR4LmtcRIhc3rOZChAd9Tj+c8B/hqIZLBIRkW7tsPQsVk1ucUU3S7CoJSmQ/wWLoUHmcGeFTRaz1rPYzEnPojOntmqk/v9nH3sWiYiIvJ4JvS+n0Lp6S331aGH+cd+ReVyKSivc1j7Sv51Z5mDxsh5N1f8rbOYtnljqz7Us+d6WnsX1acf8ZrlJ9iwSEZFue4UwdNgoMkQSo6tnQtcV1vlNjAlT60szycV/VFaarD2LN53ZWpDDggMGLTHlUG5xvXoW2yZEqSHq4rJK2VplhRejYrBIRET6zoROcj0T2hbuq81bZJKL/0DGMgK6kKAA6dkiTro2i7XrXTyxeotrPYuBgQHSq2Ujv5q3yGCRiIh0ndxyMkPQmu4tLEkunLfod0PQbRKiJDgoUPq3S1CXV+zKVtMRjlkyo51lQ9fkVMtQtL/UW2SwSEREhk1u0fSwBIsbmBHtN3ZZhqBRoB36W1Zv+XPnEWuNxajQIFXw3VW9LUkujsrnYGWYd/7YY/3+GgGDRSIiMnzPopYRjd6mgpLyk3480j9tvmL7RPP354w2TSQ4MEANT/+1O0dta9oool5THHq1aqTmQO7PKbQuGaiZ8+tOeeqHzepkFAwWiYhIdzBMeOBoodt6FpNiwtVKHUiY2WRZ/o38I1hsZ+lZjAwNtvYIfrn2oMM1oesqNjzEGoT+Y9O7iOSZecv32D2/ETBYJCIi3WZCY83ehJPIhLbV3Vqc23/q4/krk8mkMp9th6FBm7e4xjLX0NVMaEfzFv/Zf2Le4qwlO6SozFyeKT23WJV/MgIGi0REpN9l/mx+6N0VLLI4t/FlHy+V3KIyNVTcztIDCGe1NweLGlczoR3PWzQHi7uyjsv8vw+o8xjuhn1HzL3jvo7BIhER6c52rWyOG4agq2ZEM8nF+LQh4JaNIyU8JMi6HSVvImwun1TPYmtLce6DuVJeUSkzFmxVtTwHnZIsXSxlevYeKRAjcD0FiMigUif86LHHriw1F3+FU55YIIGh9T+arcneZy/1yOMSNbQdblgT2lnP4u7sAskvLpOY8BC3PTbps2yO7RA0hAYHyultmsjv27NOumexfWK0xIQFS35JuXzy9wFZuClD0KH46MWd5OVfd6ogcp9BgkX2LBIRkX6Hod3YsxgfHSbNG5l7kjamMcnFn8rm2DrLUkKnvjUW7YpzW4ain/7enPl8XZ+W6jubGh+pLu/JNsYwNHsWvdjj1FC9TcAeJzLq95vfbePBXLODR8118Dq5MVjU6i2idMqGtGPSzyZoIGOXzbGlJbmcTDa07TrRy3dkS2lFpYQFB8r4Czuq7anxUep/o/QsMlgkhzgkS0TessmyykqLxhHSOCrUrY+NeYs/bzyshgjJ+EtFamVzbGHZv1vObK3WHI8KO7kwqLelZxFuP7uNpFiCz9SESEMluDBYJCIiXdGylbU5hu7Uo7n5x50Z0caVV1wmGXklToehMXz89NBubnmuU1s1lriIENWreNe57azbW1t6FrFSDMrn2CbZ+CIGi0REpCva+s1a9rI7aQEoenxyC8skLpJJLkadr5gUE6YCOU+KiwiRhQ8MkOCgALvnio8KtSa/HMgpdOvcW29gggsREenKRg/2LCI4bNXEPETI3kWDz1d0Y43OmmDouWrheCwh2NoyFL3XAEPRDBaJiHRm+vTpcvrpp0tMTIwkJSXJ0KFDZdu2bXa3KS4ulnvvvVfi4+MlOjparrnmGsnIyBBfh94+bZ6XJ4JF2x7L9WlcycWfyuY0tNaWoei92b6f5MJgkYhIZ5YtW6YCwVWrVsnixYulrKxMLrroIikoOPGjM378ePn+++/l888/V7c/dOiQXH311eLrtN4+9P41inRvcoumh7aSC5NcDGlnhj6CxTZasGiAjGjOWSQi0pkFCxbYXX7vvfdUD+OaNWtkwIABkpubK2+//bZ8/PHHMnDgQHWbd999V0455RQVYJ555pniqzyZ3FKtZ5HBorF7Fh2UzWlIreONkxHNnkUiIp1DcAhNmjRR/yNoRG/joEGDrLfp3LmztGrVSlauXOnwMUpKSiQvL8/upEeof+ip5BZN12bmx0a9xaMFpR57Hmp4yDxGQgm0d+PqP/WRmmDuWdzDYWgiIvKkyspKeeCBB+Sss86Sbt3M5T4OHz4soaGh0qjRiRpvkJycrK5zNg8yLi7OemrZsqUuPzitt08bKvYEZK1qvT6bDukzaKb6QWBWaRKJDQ+WxCpJJw0t1aZ8Tkl5hfgy9iwSEekY5i5u3LhRPv3005N6nIkTJ6oeSu104MAB0Rv08mkrt3T1YLAI3Sy9i8yINm4mNDKSvSkhOlSiQoPEZBI5kGP+XvsqXQSLc+bMkdTUVAkPD5e+ffvK6tWrnd72vPPOU1+AqqdLL720QdtMRORpY8eOlR9++EF+++03adGihXV7SkqKlJaWyrFj9tm8yIbGdY6EhYVJbGys3UlvtMAN6+p6uj5eN0swuvEQV3IxioKScvlh/SFdJLdYy+cYJCPa68Hi/Pnz5cEHH5QpU6bI2rVrpWfPnjJ48GDJzMx0ePuvvvpK0tPTrScccQcFBcl1113X4G0nIvIEk8mkAsWvv/5afv31V2nTpo3d9X369JGQkBBZsmSJdRtK6+zfv1/69evnsx+KFixqgZwndWsea7e0IPm2X7dmyEUv/S4LN5nLR13YxfFBU0Nrk2CMjGivZ0PPnDlTRo8eLbfddpu6PHfuXPnxxx/lnXfekQkTJlS7vTbBW4OhmcjISAaLRGSooWdkOn/77beq1qI2DxFzDSMiItT/o0aNUgfa2Ceil/C+++5TgaJPZ0Jr8xU9mNxSdRgaBZOxPFxsOFdy8UWZ+cUy9fvN8uP6dOt64tOu6i7ndkwUPWhtkIxorwaLGEZBVh/m0mgCAwNVhp+zjL6qUD7ihhtukKgoc/TuKAMQJ41eMwCJiDSvv/66ddqNLZTHufXWW9X5l156Se0vUYwb+ziMyLz22ms+/SaeKJtjn7jjCY2jQqV5owiVEb0pLU/6tYv3+HOS+wPFwS/9LkcLyyQoMEDuOLuNjBvUQSJDvd4PVi0jmj2LJyE7O1sqKipUBp8tXN66dWut98fcRgxDI2B0BhmAU6dOPZlmEhE1+DB0bTDHG/O9cTKCI8dLVOAGXS1DxJ6GWo54TiwvyGDR9yzblqUCRQT9b9zSp0GmL9Q3I9rXg0Wvz1k8GQgSu3fvLmeccYZPZwASEfk7rVexbUJUgw0Ja/MWmeTimzLzzaOGCPT1GChqyVqQdrRISssrxVd5NVhMSEhQySlV1zOtKaNPg2WvMF8R83Zq4gsZgERE/k6br+jJYtxVaeV50LNIvicjr1j9nxzr3XqKNUmMCZPI0CBV+/HAUd+dt+jVYBFFZZHVZ5vRhwK0uFxbRh/WQ8U8nZtvvrkBWkpERL6+zJ+zJJfd2QVyvKS8wZ6X3B0shuv2LQ2wKZ+zz4eHor0+DI1svnnz5sn7778vW7Zskbvvvlv1GmrZ0SNGjLBLgLEdgh46dKjEx3NSMhGRr/NGsIhen5TYcFU0eUs6kx99TUaeeRg6KUa/waLtUPTebN/tWfR6ytCwYcMkKytLJk+erMpD9OrVSxYsWGBNekHdMGT82UI9sT/++EMWLVrkpVYTEZG7ZOWXSHpusWDBDU+v3FIV5rodzitWw+Cnp9qXZiN9y/SBYWijZER7PVgEFJ/FyZGlS5dW29apU6c6ZQsSEZH+bbRJbokOa9ifJSS5/LIlg0kuPqay0mRNcEmJ85GexSO+27Po9WFoIiLyb+utxbg9X1+xKm3YG7UWyXccKSiV8kqT6o1OiNZ3z2JrzlkkIiI6OX/uylb/e6P8ifacOzLzpai0osGfn04uuSU+KkxCgvTd79XGMgx98GiRlFX4Zvkcfb/DRERkaKt2H5HVe3IkJChABne1X6ChISTFhKmeKZQ22XKYvYu+tHqLL8xX1L5j4SGBUlFpUgGjL2KwSEREXoG55zMXbVfnrz+tpbRobJ7b1dClTbprxblZb9HnMqH1XDbH9jvm6yu5MFgkIiKv+HPnEVm9N0dCgwNl7MD2XvsUtKFoBou+wxcKctuyBovZvhks6iIbmvxLYGi4tH70B283g4i83Kv44uJt6vyNZ7SSpnERXg8WNzDJxWf4Us8itLZkRO/z0YxoBotexKCJiPzV0m1Z8s/+Y2ou1z3nt/NqW6xJLhn5UlxWIeEhQV5tD7lSY9FXgsUon17FhcPQRETU8HMVF5vnKo7ol+r1FTiaxYVL48gQVYple0a+V9tCdYNC6r40DN1a61nMYc8iETnBXmSiExZtzlDL+0WGBsmdA9p6/a1BAgJ6F5fvyFbt8ka9RzLmUn9Vg8WDOUUqKzooMEB8CXsWiYioQVfeeMnSq3hr/1SJ10lBZa04N8r4kL6hVuGRAt+as9g0LkKVhyqtqJT0XN8rn8NgkYg8msiEE84TwXf/HpKth/MlJixYxuigV1FzYRdzjcdFmzLkeEm5t5tDNcg+XiJY8Re9c/FRoT7xXgUFBkhLS2mo/T6Y5MJgkYiIGkRWfolM/X6TOo9AsVGkfn7oe7VsJG0To6SorEJ+2pDu7eZQnYagwyTQh4ZzW/vwGtEMFomIqEGSWh7/ZoMcLSyTU5rGyp3nejcD2tG8xWtObaHOf7HmoLebQ3WqsehbIxattYzoHN/LiGawSEREDTL8vHBThgQHBsgL1/VQhbj15upTm0tAgHne4gEfzVr1r7I5+pjvWletmnAYmoiIyOk6vlO+Mw8/Y6WWrs3MySR6TEI4u32COv/lWvYu6r9sjm/1LKYmcBiaiIjI8fDz1xvlWGGZdGkaK/ee771l/epCG4pGsIjMbdIfX1u9RdOqiXkYev+RAvV34Uv0Nw5ARESGGn5GXUWUDXnx+p4SEqTvn53BXVMkOixYDuQUyd97WUZHz3MWkeDiS1o2iVDTHApKK+RIQan4En3/1RIRkc/KLSqTyd+ah5/vG9hBJbboXURokAzpnqLOcyhanzJ9tGcxLDhImlnWQPe1Zf8YLBIRkUdsTMtVASOW07v7PH1lP9fk2j4t1f8/rk+XwlLWXNSbjHzfnLNom+Syz8fK5zBYJCIij0g7Zl6pol1StO6Hn22dntpY/ahjuHDhpsPebg7ZKC6rUPNfIcUHg8VUH01y8Z2/XiIi8imHLMFi80bmoTdfgZqLKKMDX65J83ZzyMEQdFhwoMRGBPvce9PKJsnFlzBYJCIijwaLzXwsWLTNiv5zV7a1h5T0NQSNoN7XpFpWcdnnY3U8GSwSEZFHpOcW+2yw2LJJpPRu1UitQfzX7iPebg5VW73FtzKhNa20YJHD0ERERCfmLCLBxRe1TYi2C3pJR+tC++B8Rdsl/3IKSiWv2Dz30hewZ5GIiNwORYd9eRgaUuLMvVeHGSzqb6m/GN8MFqPDgiUhOlSd3+9DvYsMFomIyO2OFpZJcVmlOp/ioz2LKZaaeNrycuR9vj4M7avlcxgsEhGR22m9ignRYRIeEuST77BWmoU9i/obhvbVAxDboeh9Ob6TEc1gkYiI3O7EELTv/qg3tQQk7FnU41J/vvu9aq0luWQbvGdx165d8vjjj8vw4cMlMzNTbfv5559l0ybzsk5EROTfrMGiZSjXF2krhGQfL5GyCvOQOnmXEYahW1vL5xi4Z3HZsmXSvXt3+euvv+Srr76S48ePq+3//vuvTJkyxRNtJCIiH+PLZXM08VGhEhIUoMrnZOabhz/Je46XlKtVdXw5G9puGNrIcxYnTJgg//d//yeLFy+W0FBzRg8MHDhQVq1a5e72ERGRL5fN8eFh6MDAAOtw5+FcFubWS68iMopx8lWtLQkumN6A5QsNGSxu2LBBrrrqqmrbk5KSJDs7213tIiIiH+brZXOqzVvMZc+ibuYr+vAQNDSJClXBLnqsDx4tNGaw2KhRI0lPT6+2/Z9//pHmzc1raRIRkX87dMz3h6Eh2RIsprNnUTfrQmtZ6r4qICDAOm9xb7ZBg8UbbrhBHn30UTl8+LB6wZWVlfLnn3/Kww8/LCNGjPBMK4mIyGcgGSTTsoavLw9DQ1NLYKL1apH3nEhu8e3vlH2Si0GDxWeeeUY6d+4sLVu2VMktXbp0kQEDBkj//v1VhjQREfk3/KhXmkRCgwIlIcq3hwy1en5c8s/7DhtkGNo2yWX/Ed/IiA52dfkm9Ci+/PLLMnnyZDV/EQFj7969pUOHDp5rJRER+dwQNAItJIkYIVhkz6J+hqF9dak/R0kue48YsGcRwWL79u3l4MGDqmdxyJAhcv31159UoDhnzhxJTU2V8PBw6du3r6xevbrG2x87dkzuvfdeadq0qYSFhUnHjh3lp59+qvfzExGRe2nz+3x9CNp2fhx7Fr3PSMPQrSzD0PuNOAwdGBioAsMjR4645cnnz58vDz74oKrPuHbtWunZs6cMHjzYWui7qtLSUrnwwgtl79698sUXX8i2bdtk3rx5TKwhItJl2RzfTm6x7VlEr1YlxtbJazLyfb8gtybVMgx9IKdQyn2g4LvLhYqeffZZeeSRR+T111+Xbt26ndSTz5w5U0aPHi233Xabujx37lz58ccf5Z133lH1HKvC9pycHFmxYoWEhISobeiVJCJyJyTuYQGC5cuXy759+6SwsFASExPVlJtBgwapkRWqvWxOcwMEi1qdxdKKSskpLFVrXVPDw8imti60EXoWU2LDJTQ4UErLK1WvdUvLsLRhElyQ8YyhYvQCRkRESJMmTexOdYVewjVr1qgdr7UxgYHq8sqVKx3e57vvvpN+/fqpYejk5GQVrCLhpqLCeVHLkpISycvLszsRETlSVFSkFh3QptlgGVNMfQkKCpKdO3eqUZA2bdqo67gIQe1zFpv68FJ/GvygawHiYcuqNNTwDh4tUoFVQIBIYozvB+yBgQHSwnIwhddmuJ7FWbNmueWJUcAbQR6CPlu4vHXrVof32b17t/z6669y0003qXmK2Hnfc889UlZW5nSpwenTp8vUqVPd0mYiMjbMgcYBKaa3YMqLNoJhCz2NH3/8sSojNmnSJDU6Qs4Kcvt+DxCkxIWp9aERLHZrHuft5vilnzaY6zufkdpEwkOCxAgSYsJkd3aB+m4ZLlgcOXKkeHNoCCvFvPnmm+pIv0+fPpKWlibPP/+802Bx4sSJal6kBj2LHEIiIkcWLVokp5xySo1vTuvWrdV+BbVl9+/fzzfS4MPQkBIbIRvT8qylW6jh/bDeHCxe1rOZYd7+REuPtSGDRUCP4DfffCNbtmxRl7t27SpXXHGFCuDqKiEhQd0+IyPDbjsup6SkOLwPMqBxpG/7PNixo5wPhrVt16rWIGMaJyKi2tQWKNrCvqhdu3Z8U6vILy6TvOJy8z7bKMFiHIehvWlvdoFsSMsVVGG6pJvj+MAXJUSH+kyw6PKcRQz9YoeKuYtfffWVOt18880qYNy1a1edHweBHXoGlyxZYtdziMsYBnLkrLPOUs+P22m2b9+ugkhHgSIR0clAggv2b9gnYRQD/vvf/8off/zBN9YJrcRMbHiwWv/WCLS5l+xZ9I4fLUPQ/dslGCrBKMHyWrLyDRgs3n///epo+sCBA6rcDU4YisGkb1znCgwPY27Q+++/r3op7777bikoKLBmRyMgxXCPBtcjG3rcuHEqSETmNBJckPBCROROX375pSrlhUS+f/75RyXLQW5urtrvUG3zFY3Rq2ibfcsEF+/4/t9D6v/LezYVI0mwJOpkHy8VvXP5sA/lJJAFaJv5HB8fr0rqoOfPFcOGDZOsrCy1GgyGknv16iULFiywJr0gCEWGtAZzDRcuXCjjx4+XHj16qPqKCByxVjURkTshKxrlvHDQ+umnn1q3Yz+H66jmTGijzFeEppZai+xZbHg7M4/L1sP5EhwYIIO7GmcI2rZn0ReGoV0OFjH/Lz8/v9p2LPtXn6HgsWPHqpMjS5curbYNw0EsWUFEnoai/1j3vqq4uDhVToccY88iudMP6829iud0SJBGkcaabpagzVk04jD0ZZddJmPGjJG//vpLFcnECcHbXXfdpZJciIiMAIl2mCNdFeYrtm3b1itt8gWHLEv9NTVI2RzbVVyOl5SrBB5qGIgvrFnQPYyTBa1JtBmGxms1VLD48ssvqzmL6OHDes44YVgGa0bPnj3bM60kImpgqJ+IaS44MA4ICJBDhw7JRx99pErmYP40+UfZHECiTowlWUdbn5g8b1tGvhqGDg0KlAu72tdkNtIwdGlFpbWCgGGGoRs1aiTffvutOuLWSucgOxrBIhGRUWDJUVReuOCCC9RyfxiSxjQcBIv33Xeft5un+zmLRkpw0XoX8zOPy+HcEmmfFOPt5vhVYsu5nRIlNrx6gXxfFx4SpA5C8kvK1bzFuAj9vsZ61zVAcMgAkYiMCLVk//zzT1Vp4ZFHHlEHx5iX3aVLF4mOjvZ283SrstIk6ZZhaCMGizsyj1tfHzXkELSxsqCrZkSrYDG/RNolRhtnGPqaa66R5557rtr2GTNmyHXXXeeudhEReQ0K/1900UVy9OhRlbiHIPGMM85goFiL7IISKaswqeLJyQZYv9dWiqV8DoehGwZWzNl3pFDCQwJl0CnGG4KuXphb3+VzXA4Wf//9dxkyZEi17Zdccom6jojICLp166bWoyfXh6BRlzA4yOWfF58on6MVHaeGyYIe2DlJogxS3N2Xy+e4/NfsrEQOlr7CustEREaAWoqYn/jDDz9Ienq62r/Znsg/yuZoki3BInsWG8byHdnq/0u6GXcI2peCRZfD9e7du8v8+fNVIW1bKFqLoRoiIiPQRlBQEgzZ0LZzqXAZ8xrJf4JF9iw2rCMF5uApNT5KjCzBqMHiE088IVdffbVaB3rgwIFqG9Zz/uSTT+Tzzz/3RBuJiBrcb7/9xne9vpnQll44I9GW/GPPoufhgOxoobmeZaNI/WYIu0NCTKhPrA/tcrB4+eWXyzfffKPWRv3iiy/UuqlYeu+XX36Rc8891zOtJCJqYNyfuc7YPYsR1kSEkvIKCQsO8naTDKu4rFJKyyv9I1iMNvcsZhktwQUuvfRSVVaioKBAsrOz5ddff+WOlYgMB8v6vfjii3LHHXeo00svvSS5ubkef14kC+LAvFmzZmrIGwfoVXteMBWoadOm6oB90KBBsmPHDtHL6i1GDBYbR4ZIaLD5JzMzT9+9QL7uaKE5cMJ60CiIbmSJ2iouOu9ZdDlYPHDggBw8eNB6efXq1fLAAw/Im2++6e62ERF5zf/+9z+1WhUCxJycHHWaOXOm2rZ27VqPPjcOxHv27Clz5sxxeD1KlWE1rblz56oVZqKiomTw4MFSXFysk55F4w1DI2jXyucc5iouDRIsYi1o2/nCRpRoM2dRz0v+uRws3njjjda5PIcPH1ZHtAgYJ02aJE899ZQn2khE1ODGjx+vklv27t0rX331lTrt2bNHLrvsMnWA7EkoRYZs7KuuuqradfhBmTVrljz++ONy5ZVXqmlAH3zwgVqOsGoPZEPKLSqz1opr2SRSjMgaLLJ8jkfl+sl8Rdth6JLySrX2uGGCxY0bN6ritPDZZ5+p7OgVK1aoNVPfe+89T7SRiMgrPYuPPvqoBAefGAbD+f/85z/qOm9BwKodqGvi4uKkb9++snLlSqf3Kykp8Wj5n52Z+dasYSMuzaat4gIMFj1LS27B0L/RRYQGSVRokO4Lc7scLJaVlan1UQFJLTjyhs6dO6taZERERhAbGyv79+93OBUnJsZ7awMjUITkZPtVLXBZu86R6dOnq6BSO7Vs2dKt7dqecVz93yHZuOsma8EiC3N71rGiE8PQ/iAhRv/lc1wOFrt27armySxfvlwWL14sF198sdqOIZD4+HhPtJGIqMENGzZMRo0aperKIkDECfVkkegyfPhwn/tEJk6cqJJztBNejzvt0ILFJP2ub3uyuORfwzimDUNHGL9n0a7Woo6TXFxOM8K60JhH8/zzz8vIkSPVJGz47rvvrMPTRES+7oUXXlCT60eMGCHl5eXWlaruvvtuefbZZ73WrpSUFPV/RkaGyobW4HKvXr2c3g8jQtqokCfssAxDd0yO9oOeRXMiD3nG0QJzz2LjKD/pWYzW1oc2ULB43nnnqXI5mO/SuHFj6/YxY8ZIZKQxJzUTkf/BsqazZ89Ww7dYhACQCe3t/VybNm1UwIjFELTgEPtjZEUjkPUWrWexfZLxh6EzWDrHo44VmXsW4/ysZzFLx3MW61XAKCgoyC5QhNTUVHe1iYjI6zBUiyX9mjRpohL5NCihg0QXzGn0lOPHj8vOnTvtklrWrVun2tKqVSuVjY1s6Q4dOqjgEStroSbj0KFDxVuZ0Fo5mQ5G7lm0WcWlstIkgYHGLuviLccspXMa+8ucxWgDzlkkIvIHN9xwg5qjWBWqQOA6T0K2de/evdUJHnzwQXUehbgBGdn33XefGtE5/fTTVXC5YMECCQ/3Tn1DLRMawZRRM6EhKSZMEB+WV5okS8c/7L7On7Kh7RJcdDxnkcEiEZEDGNY9//zzHU7FwXWehOdAPcWqJ608GeZSoq4tsp9RiBuVKTp27Oi1z9Ga3GLgXkUIDgq0Lvu3P6fQ280xfM9inJ8Ei4mWOYt6PgBhsEhE5KQuoZbYUrV8WFERExwclc3paOCyOZrW8eY5q/uOMFj0dDa0vwxDJxqxdA4RkT9AdQdHy5iidFifPn280ia90jKhjVw2p2qwyJ5Fz0APupbg4i/BYoK1dI7BElywbumyZctUwdrSUvsXd//997urbUREXoMEEqyS8u+//8oFF1ygtiED+e+//5ZFixbxk3E4DG38nsVWTaLU//uPFHi7KYaUX1IuFZUmv1nuzzZYLCqrkIKScokKq1do5lEut+iff/6RIUOGSGFhoQoakZ2HUjooJ5GUlMRgkYgM4ayzzlLL56GmLJJaIiIi1DrMb7/9tspCpuqZ0O39oGexlWXd632cs+gRxwrMvYrhIYESHmJeBs/oosKCJSIkSAWLGIo2RLA4fvx4ufzyy9VQDJaMWrVqlSpUe/PNN8u4ceM800oiIi9AHUOse0/O7cw8bs2E9oe6eNZhaM5Z9OhSf/4yBK1JiAmVAzlFKlhsHW/uvfbpOYuo9fXQQw9JYGCgqreISeBYY3TGjBny2GOPeaaVREQNbO3atbJhwwbr5W+//VbVMcR+rur0G3+2IyPfLzKhNa0sweKRglI5XlI9AYrcUzbHHw48HBbm1um8RZeDRfQiIlAEDDtj3iKgl9Hda40SEXnLnXfeKdu3b1fnd+/erdaKxnSbzz//XNU5JPtM6A4GXrnFFupIavX/2Lvofv5WkNtXCnO7HCyiMCwmeMO5556risRimAYrCnTr1s0TbSQianAIFLXl9BAgYn/38ccfq1qHX375JT8RP1oTuqpWlmHC/TlMcvFY2Zwo/+xZzDZKsPjMM89YF6+fNm2aWvYP65FmZWXJG2+84Yk2EhF5pYRHZWWlOo+i10jsA0y7QVIf+VdBbodJLpy36HZHtYLcEf7Vs5hoKcyt12DR5QSX0047zXoew9BYYoqIyGiwr9PK56BU2Ouvv25dpzk5OdnbzdNhJrR/DENDa2ZEN0BBbj/rWYzRd61Fl3sWBw4cKMeOHau2PS8vT11HRGQEs2bNUkkuY8eOlUmTJkn79u3V9i+++EL69+/v7ebpKhM6OTbMrxIStCQXzll0P85ZLBFD9CwuXbrUYSYg1iddvny5u9pFRORVqKlomw2tQd1FVIKgE5nQ/rDMn6OeRa7i4sFsaD/rWUy09Cxm+XqwuH79euv5zZs3qwXsNRUVFWo4unnz5u5vIRGRjoSHh3u7CbqxI9O/MqE1Wh28tGNFUlZRKSFBXDnXXfy+ZzHfx4NFZAUGBASok6PhZqxu8Morr7i7fUREpFPb/azGoiYpJkxCgwOltLxSDh0r0mURZV+lrQvtL0v9aRIsCS4FpRVSVFohEaFBvhksYlI3sgPbtm0rq1evlsTEROt1oaGhKtmFQzNERP6XCe1PZXMgMDBAZURjziYyohksus/RAq3Oon8Fi9FhwRIWHCgl5ZUqI7qlZaqDXtS577x169aSmpqqSkkgSxCXtRNK6ZxMoDhnzhz12Bje6du3rwpGnUGNM62HUztxWIiIqGHlFftnJrSG8xbdr6LSJHnF5lVxGvlZUe6AgIATq7jocN5ivVerxrxFrN5SNdnliiuucOlx5s+fLw8++KBaaxqBIjIQBw8eLNu2bVO9lY7Exsaq623fZCIid/r999+lc+fO1fZDZWVlsnLlShkwYIBfv+Far6K/ZUJXy4jOKfR2UwxViknjj9+phJgwNQ9Wj/MWXQ4WsezVVVddpbIEEaRhaNo2YEOyiytmzpwpo0ePlttuu01dRtD4448/yjvvvCMTJkxweB88V0pKiqtNJyKqs/POO0/VU/z666/lzDPPtG7PycmR888/3+V9ndH4ayZ0tVqLR7iKi7sLcseEBftl0lCitTC3/motuvxpjBs3Ttq0aSOZmZlqndRNmzapI3AMTaOsjivQK7lmzRpV9NbaoMBAdRlH7s4cP35cDX9jJYUrr7xStcGZkpISVQPS9kREVBc33HCDXHDBBWr6iy3tINmfaZnQ7ZP8a75i1Z5FruLi/oLcjfxsqT9fWPLP5WARQdxTTz0lCQkJKrDD6eyzz5bp06fL/fff79JjYcksHJ1XXQ0Bl21L89jq1KmT6nX89ttv5cMPP1RzKFEg9+DBgw5vj3bFxcVZTwgwiYhqgxGMiRMnyn//+19VmBvTZaqOpPiz3Vl+Hiw20daHLuTBg5vL5jTys6X+DBksIriLiTEPOyBgPHTokDqPnj7beYSe0q9fPxkxYoQq5XPuuefKV199pTKzna1LjZ19bm6u9XTgwAGPt5GIfJ8WGF599dVqwQGs3HLJJZc4XMHKH2Va5lU1jfPPupMtm0QIjhkKSyvkiCWDl9xTkNvfyuZULZ9jiGCxW7du8u+//6rzSEiZMWOG/Pnnn6q3EWV1XIFgE1nUGRkZdttxua5zEkNCQqR3796yc+dOh9eHhYWphBjbExGRK7CPQZUGBIoYliaRLEuwmBjtn8FiWHCQNI01v3YORbuHvxbk9oX1oV0OFh9//HE19AsIEFF/8ZxzzpGffvpJXn75ZZceC/UZ+/TpI0uWLLFuw2PjMnoQ69rTiWQblO8hInKXkSNHqsUGNDiAXbZsmQoWW7VqJf5e4kTrTUuKNf/A+aMTGdFMcnHrnEU/7VlM1PEwtMvZ0Chro2nfvr1s3bpVZQc2bty4XvN4MA8IO2UkyJxxxhmqdE5BQYE1OxpDzlhGEHMPtQAVmYl4bhzlY53Wffv2yR133OHycxMROfPuu+86HKl4//33/f5NQ9YqAkbs8ptE+WcvEKAw96rdOexZdHM2tL/VWKzas6j12huizqKtJk2a1Pu+w4YNk6ysLJk8ebJKasFcRKwzrSW9oJYjkmg0R48eVaV2cFsEqOiZXLFihXTp0sUdL4WI/Bj2N670GqalpamDWX+TmWf+MWsSGeqXJU402sotrLXo3qX+/G31Fk2yZVpDfkm5HC8pV6u66EWdWoIJ3nWFhBNXIdMQJ0eqluN56aWX1ImIyN1OP/10GTp0qBqpwHlHkCj32WefyezZs2XMmDEuV4EwAm2FiURLT4g/9yzC/iMszO3WbGg/DRajw4JVMXIUJ8ea43qqYVqnYBElZ2wzBFGkFtswdAyolYghYVeCSiIivcHKVNOmTZMLL7xQLSOKkYtmzZqp8xjVwPWo63rqqaeq5L4hQ4aIXye3+Hmw2FqrtchVXNw8Z9E/h6GhWaMIFSym+WKwaDt359FHH5Xrr79erbSirQeNJJN77rmHmcZE5NPi4+PVqlIIGLGS1B9//KHmRBcVFanqDTfddJOat42qEP6MwaJZa0utRbwfhaXlEhmqn2FDnw4W/XCpP03zRhGyJT1P0o4WiZ64/M1GQWzsQLVAEXAeiSoojo2EEyIiX4Ys6GuvvVadqDoGi2ZxkSESGx4secXlat5i5xSWZnNHgou/ls6B5o3M8xYxDK0nLs9MLi8vVxnQVWGbVlKHiIiMKzO/2K7Uhz+zJrlw3uJJKSmvUAXO/T1YbNYoQpfBoss9iyhpM2rUKNm1a5cqdQN//fWXPPvss9ZyN0REZFzsWbSvtbghLZcZ0Scp1zIEHRggEhPuv8P5zRubg0XMWdQTlz+RF154QRWnffHFFyU9PV1tQ0HsRx55RB566CFPtJGIiHSYDZ0U45+rt9hqbcmI5iou7lnqD9nAgYgYxd97FovFp4NF1Dz8z3/+o055eXlqG5fQIyLyH+xZPIEZ0e7h70v92Sa4wOG8YimvqJRgndQxPalWcK1lIiL/UlxWIfnF5eq8v5fOgZaWnsWDR1lr0S09i35aY1GDecAhQQFqhaQMHa3k4r8TA4iIaoG52ViCdMuWLeoyVooaN26ctGvXTvy9VzE0OFBlAvs7bdWNLMuqNlQ/7Fk0wxB807gINQcWSS5aT6O36aN/k4hIZxYuXKiCw9WrV0uPHj3UCcl8Xbt2lcWLF4u/ytQKckeHSQAWh/ZzWu8qlmgrsmTzUv2X+vPX1VtsNdNh+RweFhIROTBhwgQZP368qvRQdTsWJ8AqL/7cs5gUyyFoiAkLlvCQQCkuq1TvDbKjqf41FhtF+PecRWjeCN+hHDmoo8Lc7FkkInIAQ88oE1bV7bffrpb981fWdaFZY1FB76rWu6jVnyTXHSsw9yw2Zs+iNDdKz+KSJUvUKTMzs1ohbqzwQkTk6xITE2XdunXSoUMHu+3YlpSUJP4qK89SkJvJLVYoIXQgp8ja60quO1Zk6VmMYs9iMx0W5nY5WJw6dao89dRTctppp6n6ipyzQkRGNHr0aBkzZozs3r1bLWUKf/75pzz33HNqeVPx955FBotWWi+rNp+T6p8N7c/rQuu5MLfLweLcuXPlvffek1tuucUzLSIi0oEnnnhCYmJi1AIEEydOVNuaNWsmTz75pNx///3ir6xzFlmQ20qbv8mexZNfwcXf6yza9iymHS0Sk8mki045l4PF0tJS61E2EZFRYQeNBBec8vPz1TYEj/6OBbmrS+KcRfcluHDOojSLMweLBaUVkldcrla18bkElzvuuEM+/vhjz7SGiEiHECQyUDRjsFjdiQQXDkPXB3rPjmnD0AwWJSI0SOItczfRu6gHdepZtJ2fg4SWN998U3755RdVdywkxD7inTlzpvtbSUTUAE499VSVvNe4cWPp3bt3jcM/a9eu9csfdc5ZrE4bkucwdP0UlVVIaYU5WZbD0CeGoo8UlKokly7NYsUngsV//vnH7nKvXr3U/xs3bvRMq4iIvODKK6+UsLAw63k9zBXSE/T+lFWY1PmEaM4t07Bn0T3JLVjmLjI06CQfzTiFuTek5cqhXB/qWfztt9883xIiIi+bMmWK9TwSWcie1quIocKwYP6oV52zeOR4iVrTNyiQBxmuOFqgzVcM5QGaXWFu/QxDuzxnEQVptcnetgoKCtR1RERG0LZtWzly5Ei17ceOHVPX+fV8RRbkthMfHSaIDytNIkcKOG/RVbmWpf5YkLv6kn96KZ/jcrD4/vvvS1FR9cZj2wcffOCudhERedXevXuloqL6Wr8lJSVy8OBB8UdMbnEMPYlNolg+p7641F91zXVWmLvOpXPy8vLU5Gac0LMYHm6OegE71J9++smvVzUgImP47rvvrOcXLlwocXFxdvs6JMC0adNG/JG2nB0Lcjseis4+XqIyors2+CdjkILczITWbWHuOgeLjRo1UnMJcOrYsWO167Edq7sQEfmyoUOHWvdpI0eOtLsO1R9SU1NVoW7/Lsht7kWjE1QAnc6M6PrItdRYZCZ09cLcOPgoLa+U0GCXB4K9EywiyQW9igMHDpQvv/xSmjRpYr0uNDRUWrdurVY3ICLyZdp69+g9/PvvvyUhIcHbTdINDkM7pwXQLJ/jupwC9ixWhTqLYcGBUlJeKRl5xdKyiTnhRffB4rnnnqv+37Nnj7Rq1YoZS0RkaNjXkT3WWKxD+Zw881A91V3asUL1f9O4E9Pb/F1AQICat7g7u0AOHi3yjWBx/fr10q1bNwkMDJTc3FzZsGGD09uiUDcRkRGgysOyZctk//79aqlTW/64PnRmnpYNzR91pz2LlvJCVHcHcszz8rwdEOlxKHp3doEuklzqFCyiCPfhw4dVAgvOI+LFkHRV2O4oe5CIyNdgMYIhQ4ZIYWGhChox9SY7O1siIyPVvtAfg0UtEEqK5ZzFqpJiw+0Caqq7A0fNPYsMFh2Xz/GZYBHDMYmJidbzRERGN378eLn88stl7ty5KiN61apVKsHl5ptvlnHjxom/KSmvsK7fyzqLzoeh2bPoeo3F/OJydb6FJQOYqhTm9pVgEckrjs4TERnVunXr5I033lDTb4KCglR9RRTjnjFjhsqSvvrqq8WfHDleal2SLS4ixNvN0e0wNHoWMfLGpSLr5kBOoXX5yMjQOqdR+IVmOirM7XIuNpJbRowYIW+//bbs2rXLM60iIvIy9CIiUAQMO2PeIqCX8cCBA+JvtCzfBKxWwuXsnPYsFpVVSEEpp2PV1UHLEHTzxpyv6Kwwt08Gi88884wqyP3cc89Jhw4dpGXLlmpYZt68ebJjxw7PtJKIqIH17t1blc7RqkFMnjxZPvroI3nggQdUwp+/Qb03YEFux9ArFh1m7hljRnQ9kls4BO20MDfmLDrKE9F1sIjA8M0335Tt27dLWlqaPP/882r7PffcI507d/ZEG4mIGhwOjJs2barOT5s2TRo3bix33323ZGVlqX2gv2FBbhfmLVoCa6odk1ucS7GUEiouq7SucuMt9ZoggOzAP/74Q5YuXaqKdSNrEEfa5513nvtbSETUwHAUj6FnrQcR5xcsWODXnwMLctctWNyTXWDthaXaoYYgtOQwdDVhwUHqO4W/vbSjRdIkKlR8Jljs37+/Cg5POeUUFRxOmDBBBgwYoI66iYiMEiy2b99eNm3apKbbELJ8LetCR7NsTq1JLgwWXU5wadmEmdDO5i2qYPFYkXRvcWKdet0PQ2/dulWioqLUkDNOCBoZKBKRkSCxBUHikSNHvN0U3WDPYu04DO36QZnWs9iCPYs1Jrl4u9aiy8Eidp6//vqrnHnmmbJw4UI566yzpHnz5nLjjTeqJJf6mDNnjqSmpqrEmb59+8rq1avrdL9PP/1UlScYOnRovZ6XiMiZZ599Vh555BHZuHEj3yS7BBeu3uJMkuW9ycznkn91kX28VGWPBwScKBND+iyf43KwiOAMS/ph9YIvvvhCfv75Z7nwwgvl888/l7vuusvlBsyfP18efPBBmTJliqxdu1Z69uwpgwcPlszMzBrvt3fvXnn44YflnHPOcfk5iYhqgxJhOHDFPikiIkKt4GJ70oP6HmjXB3sWa8eexfolt6TEhqv5eVRD+RxLD6zPzFlEQIfEFpyQ5JKfny/du3eX++67T5WXcNXMmTNl9OjRctttt6nLWC3hxx9/lHfeeUfNh3QESwredNNNMnXqVFm+fLkcO3bM5eclIqrJrFmzdP0GaQfa2GciUER7caC9bds2lZDj7uFCZkO7sD405yy6Nl+RQ9BOtYo315/ce6RAfCpYPOOMM1T9MQSGCPKQ3IIitfVRWloqa9askYkTJ9rNFRo0aJCsXLnS6f2eeuoptTMcNWqUChZrglUXcNLk5eXVq61E5F+wSoue1edAu77yS8qlpLzSWpSbHGPPomus8xWZ3OJUh6QY9f/u7AKpqDRJkJcK4rscLObk5EhsbKxbnjw7O1v1EiYnJ9ttx2Uk0jiC3kysHoOluOpi+vTpqgeSiMgo6nOgfTIHzlpPWUxYsESEcriwtp7FIwWlUlZRKSFBLs/08svVW5jcUvMwdHhIoKq1iJ7Y1IQo8QaXv8nuChTrA0Pet9xyi0qkSUhIqNN9sDPNzc21nvxxmS4iMpaaDrQPHz7s9MAZo0DaCatv1RXWO4bEWPYq1qRxZKgEW3p+tLW0yTmu3lI7LK3ZNiFand+ReVy8xaurdiPgCwoKkoyMDLvtuJySklLt9liLGoktl19+uXVbZaV5aCQ4OFjN1WnXrp3dfcLCwtSJiMif4cAZcxxtexbrGjBmHbcEixyCrvWHHcP0h/OKVUa0tgIHOcbVW+qmfVK0bE7Pk52Zx+XCLvYHiA3Fq33koaGh0qdPH1myZIld8IfL/fr1q3Z71HXcsGGDGoLWTldccYWcf/756rwrR8pERL7K1QNtwEEzRoZsT3XFTOi6S7L0vmq9seQY5t9ptQNbNjEncZDzYBEQLPplzyLgSBcTyU877TSVPIOMvoKCAuukbZSvQB1HDKGgPIS2/JamUaNG6v+q24mI3EWbvqKXA1LbA22tzqx2oD127Fi3P1/zRuEy6JRk6dXSvL8l57TeV603lhzLyCuWsgqTGrZH6RyqQ7CY5cfB4rBhwyQrK0smT56s5tr06tVLrcGqzcXZv3+/mrhNRNSQysvLVXLcyy+/LMePm3fS0dHRqkwY6sKGhITo+kDbnS7u1lSdqHbsWXStbE6zRhFey/D1tWBxV+ZxVcYK9a51HyxiUvV7772njmBROFubM6jB6i6uwpGws6Nh1HOsCdpCRORuCAq/+uormTFjhnVaDDKNn3zySbWS1euvv67rA23yds8iV3GpyQFL2RyuCV271PgoFVAfLymXjLwSr8yFdTlYHDdunArQLr30UjX0640Il4jI0z7++GO1pOgll1xi3YbVqzAUPXz4cK8Hi7UdaJN3JFqGVDlnsWYsyF13ocGB0rpJpKq1iHmLPhEsYuf52WefyZAhQzzTIiIiHUBCCJbSq6pNmzZqziBRjau4cM5ijZgJ7fpQtDlYzJezO9StdKA7uTwZEDvJ9u3be6Y1REQ6gR67p59+2q6QNc5PmzaNvXlU6you7Fms4+otjc1rH1Pd5i16q9aiyz2LDz30kMyePVteffVVDkETkWH9888/am52ixYtpGfPnmrbv//+q1ZPueCCC+Tqq6+23hZzG4mqrg/trWQEX3DQkuDC1Vt8o3yOy8Eiltv77bff5Oeff5auXbtWywjkTpOIjABlua655hq7bXopnUP6pa2dXVpRKXlF5RIX6d2seT0qLa+U9DxzAhATXFzMiM7ykWARO9CrrrrKM60hItKJd99919tNIB8UHhIkcREhkltUplZxYbBYHYpxm0x4rwK5KlAdtUs0B4vZx0vlWGGpNIoM1XewyB0oERFRzfMWESxiKLpDcgzfKifJLRiC5jB93USFBUuzuHA5lFushqJPS20iPlGUG/W9sBYzdOrUSRITE93ZLiIir/viiy9U9QcsDoC5irbWrl3rtXaR/uct4gc9M5+ruDhyIMdSY5HJLS5plxTttWDR5WxorBBw++23S9OmTWXAgAHq1KxZMxk1apQUFpqPFoiIfB1WbsFqKChyjWQXrJISHx8vu3fvtqu9SFRTkgtVd9CmZ5F8I8klsD5LTC1btky+//57OXbsmDp9++23ahsypYmIjOC1116TN998U1555RVVMuw///mPLF68WO6//37Jzc31dvPIF8rn5HMVF0e4eovvrRHtcrD45Zdfyttvv62OrGNjY9UJBbrnzZunhmyIiIwAQ8/9+/dX5yMiIiQ/P1+dv+WWW+STTz7xcutIz5Itq7gczmPPoiNcvaV+2if6UM8ihpodrT2alJTEYWgiMoyUlBTJyclR51u1aiWrVq1S5/fs2aPq5xE50zTOXGj6cK55bh45HoZu2YTD0K7QkqVQ0LywtFx0HSz269dPpkyZIsXFJ7rXi4qKZOrUqeo6IiIjGDhwoHz33XfqPOYujh8/Xi688EIZNmwYy4dRjbS1e9NzOQxdFYIclH+Blpyz6JImUaHqBLuzCkTX2dBYvWXw4MHVVjUIDw+XhQsXeqKNREQNDvMVKysr1fl7771XJbesWLFCrrjiCrnzzjv5iZBTTS3BYkZesVRWmiQwkKu4aNIsy/zFhAezBmU9h6JXF+SooehuzeNEt8Fit27dZMeOHfLRRx/J1q1b1bbhw4fLTTfdpOb1EBH5uvLycnnmmWdU5QccGMMNN9ygTkR1SXBBfFhWYZIjBaXWhBeyXROaQ9D1LZ+zeq85WNR9ncXIyEgZPXq0+1tDRKQDwcHBMmPGDBkxYoS3m0I+KCQoUAWIGXklcji3mMGijUOWeZzNG5l7X8k3yufUKVjEvB1kP2MdaG0OjzMYoiEi8nUXXHCBKgmWmprq7aaQD0qJi1DBYnpukXRv0XDDhXqH4Nl2Xif5RvmcOgWLQ4cOlcOHD6uMZ5x3Bsv2VFRUuLN9RERegQPkCRMmyIYNG6RPnz4SFRVldz0PjKkmTWPD5V9VPodJLg6DRUt5IapfsLg3u0DKKipVL7ZugkVtknfV80RERnXPPfeo/2fOnFntOh4YU22YEe2YFjyj55Vch/WhI0ODpLC0QvYdKbQGj57mckj6wQcfSElJ9UKjWDcV1xERGQEOjJ2dOIJCtWlmmZOXfoy1Fm1p5YS0jHFyDQ5U21mLc5sXCtBlsIh6Y46WusLqBriOiMgIeGBMJ0PrOWOtRXucs3jytN7E3dkF+g0WsXIBItuqDh48KHFxnMRLRMbAA2M6GVrPGecsnpBfXCbHS8wrj3DO4sn3WmuBt65K5/Tu3VsFiTghSxClJTQYksESWBdffLGn2klE1KB4YEwnQwuG0LPo7Lvkb1CkXCvIHRVWr8p9JCe+W7oMFrUs6HXr1qkVXKKjT0yqDA0NVeUlrrnmGs+0koiogfDAmNwh2fKDXlpeKUcLy6zLtPkzzld07xSHhuy1rnOwiPWgAUEh1kbF8n5EREbDA2Nyh9DgQEmIDpPs4+ZaiwwWT/SEaYE0GbBnUTNy5EjPtISISAdsD4yxvF9YGJdqo/rPW0SwiB/1rs04p18LbpgJ7Z6yTFnHSxqs1qLLz4D5iS+88IKcccYZkpKSIk2aNLE7EREZwcCBAyUrK8t6efXq1fLAAw/Im2++6dV2ke9grUV76ayx6BbxUaESEhQgJpNIVn71Uoa6CBanTp2qitRiKBoldB588EG5+uqrJTAwUJ588knPtJKIqIHdeOON8ttvv6nzWMFq0KBBKmCcNGmSPPXUU/w8qO4Z0Q04XKhnGVy9xS0CAwMkKeZEApUug8WPPvpI5s2bJw899JDKiB4+fLi89dZbMnnyZFm1apVnWklE1MA2btyoRlDgs88+k+7du8uKFSvUPvC9997j50G1Ys+iPSa4uP+7pWWY6y5YxBE2dpqAjGitQPdll10mP/74o/tbSETkBWVlZdb5ir/88ot1LejOnTtLeno6PxNyodYiV3GxX+qPCS6+luTicrDYokUL646yXbt2smjRInX+77//5kRwIjKMrl27yty5c2X58uWyePFiax3ZQ4cOSXx8vLebRz4gJdayissxDkMXl1VITkGpej+Y4HLytIC7ocrnuBwsXnXVVbJkyRJ1/r777pMnnnhCOnToICNGjJDbb7/dE20kImpwzz33nLzxxhty3nnnqek2PXv2VNu/++476/A0UZ3Wh7YU5vZnmXnmRIyw4ECJiwjxdnN8XkoD9yy6XDrn2WeftZ5HkkurVq1k5cqVKmC8/PLL3d0+IiKvQJCYnZ0teXl50rhxY+v2MWPGSFRUFD8VqpVWT7CorELyisolLtJ/gyTUmtR6FbmajRt7FvU6DF1Vv379VEY0A0UiMlrpnPz8fLtAEVAiDAfKRLUJDwmyFuNO9/N5i5yv6NvD0HXqWcSwS11pk8CJiHzZ0qVLpbTUPMfKVnFxsZrHSFTX4ULM1cNQdOeUWL9907QeMG34lNw0DJ3XMGuPB7uy/FVt0FgU7XbVnDlz5Pnnn1eZ1pgX9MorrzidE/TVV1/JM888Izt37lTZihj+RhmfW265xeXnJSKqav369dbzmzdvVvslDfZvCxYskObNm/ONozrBsOvm9Dy/r7Wolc3R1jWmk5MUG9aga4/XKVisrKz0WAPmz5+vhrGRddi3b1+ZNWuWDB48WLZt2yZJSUnVbo8hIBTFRfmK0NBQ+eGHH+S2225Tt8X9iIhORq9evdSBL04Yiq4qIiJCHdAS1QVrLZpxqT/3CgsOUiu5HCkoVe+tLoJFZzAcEx5+cl3KWA1m9OjRKuADBI2o1/jOO+/IhAkTHE46tzVu3Dh5//335Y8//mCwSEQnbc+ePWpYp23btmrFlsTEROt1OEDFgWlQUBDfaXJxFRfOWQTWWHQfvJcqWMwrki7NYvW3NvTTTz+thmFQlHv37t1qO0rovP322y49FuYDrVmzRi2jZW1QYKC6jAzr2mCHjjI+6IUcMGCAw9uUlJSobEbbExGRM61bt5bU1FQ1onLaaaepy9qpadOmDBTJJdqwa0Mty6ZXnLPoyfI5nl8f2uVgcdq0aWqpqxkzZqijbE23bt3Usn+uQFkKBJ/Jycl223HZdp5QVVg1BoEqnv/SSy9VQ0IXXnihw9tOnz5d4uLirKeWLVu61EYi8l+7du1S9WRxAIvT/fffr7YR1RXXhxYpr6iUzHxzsMyC3L6ZEe1ysPjBBx/Im2++KTfddJPdETYSU7Zu3SoNISYmRtatW6dWjUHwijmPyFx0ZOLEiSq41E4HDhxokDYSkW9buHChdOnSRQ1F9+jRQ53++usvtbILVnQhqgvOWRTJOl4ilSaR4MAAiY82J2aQO3sWPT/FweU5i2lpadK+fftq2zFkg+xkVyQkJKiAMyMjw247LqekpDi9H4aqtTZgMvqWLVtUD2LV+YyAtV219V2JiOoKc6bHjx9vtxCBtv3RRx91OppBZEvrSTteUi75xWUSEx7it0PQSTFhEhTo2RIv/iTZ2rOow2FoHGk7qjH2xRdfSO/evV16LAwj9+nTx7p8oBZ04jKKfdcV7oO5iURE7oKD0FGjRlXbjmVNUVKHqC4iQ4Oty9s11Gobup2vaAluyPeSp1zuWZw8ebKMHDlS9TAiSEPdQySYYHgaZWxchSFkPB4mkqO2IkrnFBQUWLOjseY0kmnQcwj4H7dt166dChB/+ukn+e9//yuvv/66y89NROQMsqAx3QW1XG1hm6OyXkQ1/ajnFpWpJJcOyTF+90ZpyT1NWWPRZ9eHdjlYvPLKK+X777+Xp556Sq2PiuDx1FNPVdvqMyyDZbOysrLU4yCpBcPKKHqrJb3s379fDTtrEEjec889cvDgQVXvDPUWP/zwQy6/RURuhZJeWAcaFR/69++vtv3555/y3HPPqYNcorpCj9rWw/l+27OYYUnA0NbKJvfQemrzisulsLRc9WJ7ikuPXF5erlZPwTCMOyd4jx07Vp0cqZq48n//93/qRETkSSgHhmS6F198USXKQbNmzeTJJ59UWdFErg4X+mv5nBM9iwwW3QnzX6NCg6SgtEIdiLRNjBZdzFkMDg5WJXMQNBIRGRlWcEGCC0YxtGoKOI+FADy9DisZS0qsudYiiif7I85ZbIgkl2J9JbhccMEFsmzZMs+0hohIZzBN5t9//1Un1IYlcpW/9yxqgQx7Fn23jqfLA9yXXHKJKh2xYcMGlcmMeYu2rrjiCne2j4jIKzA/GgW5kbyHZD5AqS8k3WEhgMjISH4y5FrxZD8MFrHSmva6OWfR/bT31NM9iy4Hi0gu0dZ0rgpDM1iRhYjI1yGJBaMoSN4766yz1DasQY/5ig899BArMFCd+XPPYk5BqZRWmA+2GCz6Uc+idoRNRGRkX375paofa1vsf8iQIaoKw/XXX89gkVzuWUT5HE9nreqNFiAnRIdJaLDLM99IJ+VzXPrksEILklw2btzouRYREelAYWFhtXXrATUWcR2RK1mr0WHBftm7qJXNSYnjSmqekGKpXam9z7oIFkNCQqRVq1YcaiYiw8MqUlOmTJHi4hM74aKiIpk6dapLK0wR2a0Rfcy/gkUtONYywskzPYuePghxuS980qRJ8thjj6lVU5o0aeKZVhERedns2bNl8ODB0qJFC+nZs6fahozo8PBwWbhwobebRz6mReMI2Zl5XPblFMjZkiD+QhseZSa0ZyRbemyzjpdIWUWlhAQF6iNYfPXVV2Xnzp2qOG3r1q2rZUOvXbvWne0jIvKKbt26yY4dO+Sjjz6SrVu3qm3Dhw+Xm266Sc1bJHJFx+QYWbotS3ZkHPerN07L0uW60J6REBUmwYEBUl5pkqz8EmnWKEIfweLQoUM90hAiIr1BeRws+0d0sjokmVfX2HY43z8LcnOpP48IDAxQWeZpx4pUYK6bYBFzeIiIjO7IkSMSHx+vzh84cEDmzZun5ixefvnlMmDAAG83j3xMp5QY9f+OTP8KFtNzzavWcBjac9Bri2Axw4PzFuudv79mzRrZsmWLOt+1a1fp3bu3O9tFROQVWHAAASECxA4dOsinn34qF198sSrSHRgYKC+99JIqqcNRFnJFe0vPYvbxUjlyvETio/0jOzgjr0T9z2Fo305ycXkmZGZmpgwcOFBOP/10VZwWJ6zkgmUAsSwWEZEv+89//iPdu3eX33//XdVYvOyyy+TSSy9Va0MfPXpU7rzzTnn22We93UzyMait2KqJedWf7X4yb/F4Sbk6AQtye44WiHuyfI7LwSKWv8rPz5dNmzZJTk6OOqHuYl5engociYh82d9//y3Tpk1Tq7a88MILcujQIbVyFXoVccI+UEt4IXJFx2Rz7+L2DP8YikYPKkSGBkmUpc4k+UnP4oIFC+S1116TU045xbqtS5cuMmfOHPn555/d3T4iogaFA+CUlBR1Pjo6WlV8aNy4sfV6nMcBM1F9MqL9KVjEkDs0iQr1dlMMLVlb8k9PPYtY7g/FuavCNi4FSERGgHXua7pMVB/+FixiXWiIZ7DoUU0bYBja5X5hzFccN26cfPLJJ6rWIqSlpcn48ePVvEUiIl936623SliYOQEBK7jcdddd1pqyJSXmoTWi+geLx8VkMhn+IEQbhvaXZB49DEN76ntVr6LcV1xxhaSmpkrLli3VNmQNooDthx9+6PYGEhE1pJEjR9pdvvnmm6vdZsSIEQ3YIjKKtolREhggkltUJpn5JYZP+jhi6VnkMLRnJcWag/HS8ko5VlgmjT3Qk+tysIgAEau0/PLLL9ZJ3pi/OGjQILc3joioob377rtefdORXPPjjz/KunXrJDQ0VI4dO1btNvv375e7775bfvvtNzWvEgHu9OnTJTiYSQR6Fh4SJKkJUbI7q0ANRRs+WLTMWYyP5pxFTwoLDlJD/QjO0buoi2AR0MV54YUXqhMREblPaWmpXHfdddKvXz95++23q11fUVGhSvkgCWfFihWSnp6uejoxb/yZZ57hR6FzHZNiLMHicTmnQ6IYWU6BZRiacxYbpI5nk4JSKa2o9G6Cy6+//qqynlEipyrUH0Nh7uXLl7u7fUREfmXq1KlqDjhqPTqyaNEi2bx5s5r206tXL7nkkkvk6aefVhUpEGiSj5TP8YNl/7Rh6Pgozln0tPl39pPFD54rvVo28m6wOGvWLLVGamxsbLXr4uLiVKHamTNnurt9RERkY+XKlSqQTE5Otm4bPHiwOpBH/VtnkJiD29ieqOF1tCz7t90Plv3ThqGbcBja59U5WPz333/VklfOXHTRRWoJQCIi8pzDhw/bBYqgXcZ1zmBOIw7stZOWoEjeyYjeYcmINjKWzjGOOgeLGRkZDusrajCxmsv9ERFVN2HCBDXXu6aTp1eFmThxopoypJ1QxYIaXmp8lIQEBahl8A55cMUNb0MgfESbs8jSOT6vzgkuzZs3V8v6tW/f3uH169evl6ZNm7qzbUREhvDQQw+p2o01adu2bZ0eC4ktq1evrnYwr13nDOpGarUjyXtCgwOlTUKUSnDBvMXmjSIM+XHkl5RLWYW555QJLn4ULA4ZMkSeeOIJNRQdHm6f7l9UVCRTpkyRyy67zBNtJCLyaYmJierkDsiSRnmdzMxMSUpKUtsWL16s5pMjCZF8YyhaBYsZ+XJ+Z/NnaNT5ilGhQapkEPlJsPj444/LV199JR07dpSxY8dKp06d1HYMnSALD+UcJk2a5Mm2EhEZHmooYn1q/I/9KuotAkZ1UFMR88MRFN5yyy0yY8YMNU8R++d7772XPYc+NW8xXbYZeNk/rWwOk1v8LFjEBGrU9EIhWMx90SbmYq4NMvEQMFaddE1ERK6ZPHmyvP/++9bLvXv3Vv+jAPd5550nQUFB8sMPP6h9MXoZsQwhinI/9dRTfKt9MMnF8AW5WTbHEFwqyt26dWv56aef5OjRo7Jz504VMHbo0EEaN27suRYSEfmR9957T53qsi8m3661uCMzXyorTRKINQANW2ORq7cYQb1WcEFwePrpp7u/NURERAbXOj5KJboUl1XKgaOF6rJhy+awxqJ/lc4hIiKikxcUGCDtE829i9sMupJL9nHLnEUOQxsCg0UiIiKvDUUbc94iC3IbC4NFIiIiLy37Z9SeRWuCC4ehDYHBIhERUQPrmGRZI9qg5XO0BJcmTHAxBAaLREREDayTpWdxd1aBlFVUGrbOYgKX+jMEBotEREQNDMv8hQYFSmlFpRw22BrRKKunzVlkz6Ix6CJYREHv1NRUtYxg3759q617amvevHlyzjnnqPI9OA0aNKjG2xMREekNaismxZrX6s7IM1awmFd8Yl1oBovG4PVgcf78+fLggw+qtaXXrl0rPXv2VCvCYN1TR5YuXSrDhw9XqxmsXLlSWrZsqZa/SktLa/C2ExER1VfTuHD1/2GDBYtHLGVzosOCuS60QXg9WJw5c6aMHj1abrvtNrXe6dy5cyUyMlLeeecdh7f/6KOP5J577pFevXpJ586d5a233pLKykpZsmRJg7ediIiovpJjLcGiwYahOQRtPF4NFktLS2XNmjVqKNnaoMBAdRm9hnVRWFgoZWVl0qRJE4fXl5SUSF5ent2JiIhINz2LBgsWs1k2x3C8GixmZ2dLRUWFJCcn223H5cOHD9fpMR599FFp1qyZXcBpa/r06RIXF2c9YdiaiIhINz2LBhuGZkFu4/H6MPTJePbZZ+XTTz+Vr7/+WiXHODJx4kTJzc21ng4cONDg7SQiIqoqxaA9i1rZnHgu9WcYwd588oSEBAkKCpKMjAy77bickpJS431feOEFFSz+8ssv0qNHD6e3CwsLUyciIiI9MWqCizYM3YSrtxiGV3sWQ0NDpU+fPnbJKVqySr9+/Zzeb8aMGfL000/LggUL5LTTTmug1hIREbl/GDozr0QqK82lZoyAw9DG49WeRUDZnJEjR6qg74wzzpBZs2ZJQUGByo6GESNGSPPmzdXcQ3juuedk8uTJ8vHHH6vajNrcxujoaHUiIiLyBUkx5mARhblzCksNs9rJEW0Ymj2LhuH1YHHYsGGSlZWlAkAEfiiJgx5DLell//79KkNa8/rrr6ss6muvvdbucVCn8cknn2zw9hMREdVHaHCgChCzj5eoeYuGCRa1YWjOWTQMrweLMHbsWHVyVoTb1t69exuoVURERJ6VEmcOFrGKS7fmcYZ4u49YlvqLjwr1dlPITXw6G5qIiMiXpVjmLaYbJCMa60If1YJFDkMbBoNFIiIiL5fPMcr60HlF5VJuSdbhutDGwWCRiIjIyz2LRqm1mG1JbokJC5aw4CBvN4fchMEiERGRlxhtFRdr2RwOQRsKg0UiIiIvaRoXYaiexSPHzT2LHII2FgaLREREXsyGNlLPojUT2iBlgMiMwSIREZGXh6Hzi8uloKTcMDUWWTbHWBgsEhEReUlMeIhEhwUbpndRm7PIYWhjYbBIRETkRcmx5iHbDAPMW+QwtDExWCQiItJBrUUjFObWElw4DG0sDBaJiIi8KCU2wnDD0CydYywMFomIiHSQEW2EVVyyLQkunLNoLAwWiYiIvMgo60NXVprkaKE5WExg6RxDYbBIRETkRSmWwty+3rOYW1QmFZZ1oRtHhnq7OeRGDBaJiIi8yCjrQ2uZ0DHhwRIazPDCSPhpEhER6SAbOut4iZRVVPp8cguHoI2HwSIREZEXocxMSFCAmEwiWfnm0jO+iOtCGxeDRSIiIi8KDAyQpBjLULQPz1u0FuSO4nxFo2GwSEREpJOhaF+et2hdFzqawaLRMFgkIiLyMiMkueQUaKu3mOtGknEwWCQiIvKyZEuw6Mvlc7Itw9AsyG08DBaJiIi8rKkB1ofO4TC0YTFYJCIi8rJkbc6iD/cs7s8ptOslJeNgsEhERKSTOYu+OgyNsjlpx4okIECka7NYbzeH3IzBIhERkY6GoU0ouGiRW1gm7/65x1rDUK82pOWq/9smRElMeIi3m0NuFuzuByQiIiLXJMWaM4hLyyvlWGGZNI4KVUHj2E/WyvId2ZJ2tEgev6yLbt/W9QfNwWKPFo283RTyAPYsEhEReVlYcJA1i1ibt/jx6v0qUITN6XmiZ1qw2L15nLebQh7AnkUiIiKdzFvE+sqotRgVGizTftxivW5n5nHRsw1px9T/PVowWDQi9iwSERHpaBWXQ7lF8sgX/0phaYX0bGke1s3ML5G84jLRIyTlZOSVSGCASBcmtxgSg0UiIiId0ErOzF22S/7akyORoUHyyg29Jdkyn1GvvYvaEHTH5BiJDOWApRExWCQiItJRRvSBnCL1/2NDTpFW8ZHSPila18HihoPmIWjOVzQuBotEREQ6qrUI53RIkJv6tlLn2yeag8VdOg0W11vK5nC+onExWCQiItIB9CJCTFiwPHdNDwlAhWsEizruWUR5H2smNMvmGBYnFxAREelA3zZN5MnLu0jvVo2lWaMI6/Z2WrCYpb9gEau2IIM7JChATmka4+3mkIcwWCQiItIB9CTeelabatu1nsUDOYVSXFYh4SFBohcbLL2KnVJiVK1IMiYOQxMREelYYnSYxIYHS6VJZE92gejJv9Zi3Fy5xci8HizOmTNHUlNTJTw8XPr27SurV692ettNmzbJNddco26PI7BZs2Y1aFuJiIgaGn7v9DpvkcW4/YNXg8X58+fLgw8+KFOmTJG1a9dKz549ZfDgwZKZmenw9oWFhdK2bVt59tlnJSUlpcHbS0RE5A16DBZtk1uYCW1sXg0WZ86cKaNHj5bbbrtNunTpInPnzpXIyEh55513HN7+9NNPl+eff15uuOEGCQszFyklIiLym2CxAZNcyisq5eDRQlm9J0e+//eQpOea6z9q9h0plPzicgkNDlQFucm4vJbgUlpaKmvWrJGJEydatwUGBsqgQYNk5cqVbnuekpISddLk5el7MXYiIiJnwaKnay2it/DJ7zbJL1sy5XBesVRgoqRF80YR8u3YsyQh2txZ86+lGHeXprESEuT1WW3kQV77dLOzs6WiokKSk5PttuPy4cOH3fY806dPl7i4OOupZcuWbntsIiKihtA+0dxztzu7wC6Aczf0Fr6/cp8qiYPnQUmc1vGREh8Vqrbd8+FaKS2vtMuE5hC08Rn+UAA9l7m5udbTgQMHvN0kIiIilzRvHCFhwYEqUEMJHU/5c1e2+v/UVo1k9WMXyLanL5Flj5wv8+88U6LDgmX13hx56odNVVZuYSa00XktWExISJCgoCDJyMiw247L7kxewdzG2NhYuxMREZEvCQoMkLaJnk9yWbHriPr/vE5JkhQbLoGB2ioyMTL7hl6CRWU+XLVf/rtqn2zkMn9+w2vBYmhoqPTp00eWLFli3VZZWaku9+vXz1vNIiIi8sskl8pKk6yyBIv928VXu/6CU5Ll4Ys6qfOTv90ohaUVEhESJO0sQSwZl1dXcEHZnJEjR8ppp50mZ5xxhqqbWFBQoLKjYcSIEdK8eXM171BLitm8ebP1fFpamqxbt06io6Olffv23nwpREREHtXewz2L2zPz5UhBqUSGBjkdWr7nvHayJT1Pflifri53ax6rej3J2Lw6Z3HYsGHywgsvyOTJk6VXr14q8FuwYIE16WX//v2Snm7+QsKhQ4ekd+/e6oTtuC/O33HHHV58FURE7rN3714ZNWqUtGnTRiIiIqRdu3aqFi0OkG2tX79ezjnnHLWgARL3ZsyYwY/B4Dxda3HFTnOv4umpTVQ5HGcFwmdc20NlQMOprRp7pC2kL15fG3rs2LHq5MjSpUvtLmPlFqT1ExEZ1datW9WUnDfeeEONmGzcuFHVo8WoCw6QtRJgF110kSo1hvq0GzZskNtvv10aNWokY8aM8fZLoAYon4PfQgRu7rTCktziaAjaVmRosHww6gz55p80ufrUFm5tA+mT14NFIiI64eKLL1YnDVat2rZtm7z++uvWYPGjjz5SPY1YwADzv7t27apGZrDQAYNF40pNiBSM+OaXlEtmfokkx4a7tQD3X7tz1Pn+7RJqvT1qLd5xTlu3PT/pm+FL5xAR+TqU/WrSpIn1MhYuGDBggAoUNVgqFUHl0aNHHT4GFidAj6TtiXxLWHCQtI6P8shQ9MZDeSoIjQ0Pli7NWDWE7DFYJCLSsZ07d8orr7wid955p3UbFi5wtKCBdp0jXKDAGLTMY3cHi9oQ9Jlt45mwQtUwWCQiagATJkxQc8xqOmG+oi1UfMCQ9HXXXafmLZ4MLlBgDJ5KcllZQ8kcIs5ZJCJqAA899JDceuutNd4G8xNtqz+cf/750r9/f3nzzTftboeFCxwtaKBd52yBApzIt3kiWCwpr5C/95rnK57Vvvb5iuR/GCwSETWAxMREdaoL9CgiUMTCBe+++64EBtoPAmHhgkmTJklZWZmEhISobYsXL5ZOnTpJ48YsZWJknijM/c/+Y1JcVqmSVrTHJ7LFYWgiIh1BoHjeeedJq1atVPZzVlaWmodoOxfxxhtvVMktqMe4adMmmT9/vsyePVstdEDG1i7RnOCSlV8iuUVlbl3iD0PQ7i7HQ8bAnkUiIh1BDyGSWnBq0cK+hp1WZzYuLk4WLVok9957r+p9TEhIUIsbsGyO8cWEh0hKbLgcziuWnZn50qf1iSz5+lpZx/qK5L8YLBIR6QjmNdY2txF69Oghy5cvb5A2kb50ax6ngsVl27JOOlgsLC1Xw9B1ra9I/onD0ERERD7k8p5N1f9fr0s76VXN/t57VMorTdK8UYS0bBLhphaS0TBYJCIi8iEXdUmRqNAgOZBTJGv3Oy7CXp8l/jhfkZxhsEhERORDIkKDZHA3c4mkr/9Jq/fjVFSaZNEmc8ml/u05X5GcY7BIRETkY67ubU5++mF9upSWV9brMX7akC57sgskLiJEBp1ivyIQkS0Gi0RERD6mX7t4SYoJk2OFZbJse5bL96+sNMmrv+5U5287K1VlWRM5w2CRiIjIxwQFBsiVvZqp89/UYyj6ly0Zsi0jX6LDguW2/m080EIyEgaLREREPmho7+bq/8VbMlwq0I0M6ld/M/cqjujXWuIi2atINWOwSERE5IO6NI2VjsnRas7igo3pdb7f7zuyZf3BXAkPCZRRZ7NXkWrHYJGIiMgHodSN1rtY16xo9Cq+smSHOn9T39YSHx3m0TaSMTBYJCIi8lFX9jIHi6t250jasaJab//Xnhz5376jEhoUKGMGtG2AFpIRMFgkIiLyUVh5pW8b85J/3607VOvttQzo609vIcmx4R5vHxkD14YmIiLyYVf1bq56DN/+Y7f8e+CYhAQHSkhQgOo91IpvV5hMUlJeKX/szJbgwAC5c0A7bzebfAiDRSIiIh92Sfem8vQPmyX7eKks2HS41ttffWpzadkkskHaRsbAYJGIiMiHYQWWL+/pr3oVSytMUlZeKWUV5hMEBgZIUECAqs2IpQKv6Gmuz0hUVwwWiYiIfFznlFh1IvIEJrgQERERkVMMFomIiIjIKQaLREREROQUg0UiIiIicorBIhERERE5xWCRiIiIiJxisEhERERETjFYJCIiIiKnGCwSERERkVMMFomIiIjIKQaLREREROQUg0UiIiIi0newOGfOHElNTZXw8HDp27evrF69usbbf/7559K5c2d1++7du8tPP/3UYG0lIiIi8ideDxbnz58vDz74oEyZMkXWrl0rPXv2lMGDB0tmZqbD269YsUKGDx8uo0aNkn/++UeGDh2qThs3bmzwthMREREZndeDxZkzZ8ro0aPltttuky5dusjcuXMlMjJS3nnnHYe3nz17tlx88cXyyCOPyCmnnCJPP/20nHrqqfLqq682eNuJiIiIjC7Ym09eWloqa9askYkTJ1q3BQYGyqBBg2TlypUO74Pt6Im0hZ7Ib775xuHtS0pK1EmTm5ur/s/Ly6tzOytLCsXXufJ6ga/ZPz5nI3zWrrxm7bYmk0n8nfYe1Oc7Q0S+zdV9oVeDxezsbKmoqJDk5GS77bi8detWh/c5fPiww9tjuyPTp0+XqVOnVtvesmVL8Sdxs8Tv8DX7h/p8zvn5+RIXFyf+DO+BP+4Licj1faFXg8WGgF5L257IyspKycnJkfj4eAkICBA9RPfYWR84cEBiY2PFX/jj6+Zr9v7njKNo7BybNWsm/g7vAf7+YmJi6rQv9OXvry+33dfbz7br8713dV/o1WAxISFBgoKCJCMjw247LqekpDi8D7a7cvuwsDB1stWoUSPRG3yQvrYTcAd/fN18zd7l7z2KtlN+WrRo4VffX19uu6+3n23X33vvyr7QqwkuoaGh0qdPH1myZIldzx8u9+vXz+F9sN329rB48WKntyciIiKi+vP6MDSGiEeOHCmnnXaanHHGGTJr1iwpKChQ2dEwYsQIad68uZp7COPGjZNzzz1XXnzxRbn00kvl008/lf/973/y5ptvevmVEBERERmP14PFYcOGSVZWlkyePFklqfTq1UsWLFhgTWLZv3+/Gi7R9O/fXz7++GN5/PHH5bHHHpMOHTqoTOhu3bqJL8IQOWpMVh0qNzp/fN18zeTLfPn768tt9/X2s+3GeO8DTKwhQURERER6LcpNRERERPrFYJGIiIiInGKwSEREREROMVgkIiIiIqcYLNYTMrfvu+8+adu2rco0QpX0yy+/vFoNSE/atGmTXHPNNZKamqpWYEDZIaO/5nnz5sk555wjjRs3ViesI7569WpDv+avvvpKlZZCMfmoqChVMeC///2vR59TD6/bFkpk4Ts+dOhQrzy/v5szZ47az4SHh0vfvn09+jd3Mn7//Xf1PcWqFPi+oFKGLeRzovJG06ZNJSIiQu0/duzYIXqA8nCnn366WlEnKSlJfde3bdtmd5vi4mK599571Qpk0dHRav9fdZEKb3j99delR48e1uLPqHv8888/677djjz77LPqu/PAAw/4RPuffPJJ1V7bU+fOnd3edgaL9bB3715VTPzXX3+V559/XjZs2KDK/Zx//vnqQ2kohYWF6sccX25nK9gY7TUvXbpUhg8fLr/99pusXLlSBTEXXXSRpKWlGfY1N2nSRCZNmqRe7/r161UNUpwWLlzokefTy+u2bc/DDz+sDhKo4c2fP1/Vw0UJjrVr10rPnj1l8ODBkpmZqbuPAzV60T4Et47MmDFDXn75ZZk7d6789ddf6uALrwU/qN62bNky9fe1atUqtdBEWVmZ2rfhNWnGjx8v33//vXz++efq9ocOHZKrr75avA0rAeF3aM2aNaru8cCBA+XKK69UHRp6bndVf//9t7zxxhsq8LU1Xuft79q1q6Snp1tPf/zxh/vbjtI55JpLLrnE1Lx5c9Px48erXXf06FHr+X379pmuuOIKU1RUlCkmJsZ03XXXmQ4fPqyu27Ztmwlv/5YtW+zuP3PmTFPbtm1d/khat25teumll/zqNUN5ebl6nvfff9/kL68ZevfubXr88cdNnqCn143Pt3///qa33nrLNHLkSNOVV155Uq+NXHfGGWeY7r33XuvliooKU7NmzUzTp0/X9duJ79/XX39tvVxZWWlKSUkxPf/889Ztx44dM4WFhZk++eQTk95kZmaq17Bs2TJrW0NCQkyff/659Tb4+8JtVq5cadKbxo0bq79bX2l3fn6+qUOHDqbFixebzj33XNO4cePUdr23f8qUKaaePXs6vM6dbWfPootycnJULwuOAHFUWpW27jSWLcSRFW6PaB5Hirt371ZFyKFjx45qaPGjjz6yuz8u33jjjaInen7N6F3FETh63/zhNeM3EEPBGJ4aMGCAuJveXvdTTz2lhuRGjRp10q+NXFdaWqp6izBcq8EiCbiMnm5fsmfPHjW9wva1YG1cDKvr8bXk5uaq/7V9Gz4H7Ots24/hxlatWumq/RUVFWraCHpEMRztK+3GPg+rwtm2E3yh/Tt27FBTLzDSeNNNN6nFTNzddgaLLtq5c6f6wbadE+AIftAxfIfVZjCkhx3SBx98oH5Y0dUN+FA/+eQT6322b9+uPlxs1xM9v+ZHH31U/ZFU/QM32mvGDwfmm2A9dezQXnnlFbnwwgvF3fT0ujGU8vbbb6t5quQd2dnZ6sdfW1FLg8sIvHyJ1l5feC04GMOcubPOOsu6OhnaiL9/7YBNb+3H/gD7KMxxvuuuu+Trr7+WLl266L7dgOAWUyy0ZYVt6b39ffv2lffee08d5GPuKA6KMGUnPz/frW1nsOiiui54s2XLFjWfDicN/nDwoeE6uOGGG9R8LMxR0XpdTj311Fp/qBuaXl8z5sjgjxw7JUy8N/JrxqT3devWqUBs2rRpag4Z5m+6m15eN3Z0t9xyiwoUExIS6v16iHwRerk2btyo9m++olOnTmofhbmgd999t4wcOVI2b94senfgwAEZN26c2j+5+3ekIVxyySVy3XXXqXmWmH/7008/ybFjx+Szzz5z6/MwWHQR1qJGttHWrVtP+s1HUgomAqN3BvC/3noV9fqaX3jhBRUsLlq0qNpkZCO+Zgz9tW/fXmVCP/TQQ3Lttdc6PAo2yuvetWuXCjSR2RocHKxO6Ln87rvv1HlcT56HQD0oKKha9iQuezqpzt209ur9tYwdO1Z++OEHlcSHxBEN2ohpAQgE9Nh+9GBhH4WRBuybkGg0e/Zs3bcbox1I1sKBrLavwQgJEqFwHr1wpTpuf1U4YMc0IIwSufW9P+nZlX7o4osvrjUBYNGiRaagoCDT/v37rddt2rRJTSz9+++/rdvee+89U1JSkmnFihWmwMBAU1paWr3a5OkEFz295ueee84UGxvr8cnFenrNVd12221qErZRX3dRUZFpw4YNdicktwwcOFCdLykpcctrpboluIwdO9YuwQXfD19NcHnhhRes23Jzc3WT4IL2IZEIyUPbt293mqzwxRdfWLdt3bpVN4kWVZ1//vkqKU3v7c7Ly6u2rznttNNMN998szqv9/Y7StRBctHs2bPd2nYGi/Wwa9cutdPp0qWL+hDwh71582b14XTu3Nn6h9+rVy/TOeecY1qzZo3pr7/+MvXp06faDzy+qBERESqb6YILLrBuP3jwoKlTp07qfs7gB/Off/5Rp6ZNm5oefvhhdX7Hjh2Gfc3PPvusKTQ0VLUhPT3desIfiFFf8zPPPKOCM7QHz48fu+DgYNO8efPc/pr19LqrYja0d3z66acqoELgj+/BmDFjTI0aNbJmvusJ9gPaPhE/iMi+x3lk7mv7D7T922+/Na1fv14dgLRp00YdnHjb3XffbYqLizMtXbrUbt9WWFhovc1dd91latWqlenXX381/e9//zP169dPnbxtwoQJKmt7z5496n3F5YCAALXf0nO7nbHNhtZ7+x966CH1ncF7/+eff5oGDRpkSkhIUNn07mw7g8V6OnTokDoKRI8eghccaaOMyG+//Van0iK2rr/+erVje+edd6zb8MFjm+3jVaXdpurJUz1OenjNeG5HrxnlA4z6midNmmRq3769KTw8XB0x4g8dP+CepIfXXRWDRe955ZVX1A8OvgvoaVy1apVJj/B9crR/wHdHO8h54oknTMnJySoAxgEMyjzpgaN24/Tuu+9ab4Og9p577lH7gcjISNNVV12lAkpvu/322637isTERPW+aoGinttd12CxSMftHzZsmOos0vbTuLxz5063tz0A/7h3xJyIiIiIjIIJLkRERETkFINFIiIiInKKwSIREREROcVgkYiIiIicYrBIRERERE4xWCQiIiIipxgsEhEREZFTDBaJiIiIyCkGi0RERHWQmpoqs2bNMsx75c3XExAQIN98802db//kk09Kr169arzNrbfeKkOHDnVD66gqBotEROTXDhw4ILfffrs0a9ZMQkNDpXXr1jJu3Dg5cuSIt5tmWOnp6XLJJZd4uxlURwwWiYjIb+3evVtOO+002bFjh3zyySeyc+dOmTt3rixZskT69esnOTk5XmtbRUWFVFZWipGUlpaq/1NSUiQsLMzbzaE6YrBIRER+695771W9iYsWLZJzzz1XWrVqpXq8fvnlF0lLS5NJkybZ3T4/P1+GDx8uUVFR0rx5c5kzZ471OpPJpIZL8RgIhNBTef/991uvLykpkYcffljdD/fv27evLF261Hr9e++9J40aNZLvvvtOunTpoh7jrbfekvDwcDl27JhdO9DzOXDgQOvlP/74Q8455xyJiIiQli1bquctKCiwXp+ZmSmXX365ur5Nmzby0Ucf1fi+4P2o7XnR84r3Aq8nMjJSunfvrgJuW+edd56MHTtWHnjgAUlISJDBgwc7HIZ+9NFHpWPHjupx2rZtK0888YSUlZVVa9cbb7yhXh9ud/3110tubq7T14BAe/r06er14nX37NlTvvjiixpfNznGYJGIiPwSeg0XLlwo99xzjwombKHn66abbpL58+erIFDz/PPPq6Djn3/+kQkTJqjgafHixeq6L7/8Ul566SUV0KCnEsEQAigNgqaVK1fKp59+KuvXr5frrrtOLr74YnVbTWFhoTz33HMqSNy0aZNqAwJIPLZtjyPahetg165d6nGuueYa9bi4DsEjns92Ph+G23/77TcVML322msqgHTmggsuqPV5i4uLpU+fPvLjjz/Kxo0bZcyYMXLLLbfI6tWr7R7r/fffVwH5n3/+qXptHYmJiVHB8ubNm2X27Nkyb9489V7aQq/vZ599Jt9//70sWLBAfQb47JxBoPjBBx+o58R7OX78eLn55ptl2bJlTu9DTpiIiIj80KpVqxAFmr7++muH18+cOVNdn5GRoS63bt3adPHFF9vdZtiwYaZLLrlEnX/xxRdNHTt2NJWWllZ7rH379pmCgoJMaWlpdtsvuOAC08SJE9X5d999Vz3funXr7G4zbtw408CBA62XFy5caAoLCzMdPXpUXR41apRpzJgxdvdZvny5KTAw0FRUVGTatm2betzVq1dbr9+yZYva9tJLLzl9f2p7XkcuvfRS00MPPWS9fO6555p69+5d7XY1ve/w/PPPm/r06WO9PGXKFPX+HTx40Lrt559/Vq8xPT1dXR45cqTpyiuvVOeLi4tNkZGRphUrVtg9Lt6r4cOHO31ecizYWRBJRETkD2x7DmuDeYxVL2sZxegpxHkMo6Knb8iQIWroNzg4WDZs2KB65jDUagtD0/Hx8dbL6IHr0aOH3W3Qk3fmmWfKoUOH1NA2hpAvvfRS1fMH//77r+pRtB1axmvCMOyePXtk+/btqg3oBdR07tzZen9nantevJ5nnnlG9fZhyB7zEfF6MERsy/Z5nUGP5csvv6x6SY8fPy7l5eUSGxtrdxsM72PI2/a9x2vctm2b6gmu2guJXtoLL7zQbjva2Lt371rbQ/YYLBIRkV9q3769mju3ZcsWueqqq6pdj+2NGzeWxMTEOj0e5tIhcMF8RwxNY4gUw9YY9kQAFBQUJGvWrFH/24qOjraex3A42mTr9NNPl3bt2qnh67vvvlu+/vprNWSrwWPfeeeddvMjbQMsBIv1Udvz4rVhyBgBMobbMQ8TcxO1JBYNttcEQ/MITKdOnarmNMbFxannfPHFF6W+8J4AhshtA0xgYo3rGCwSEZFfQo8eep4wfw/z2WznLR4+fFj1pI0YMcIueFu1apXdY+DyKaecYr2Mx0BvIk5InkEPHnoV0ZuFnjjME0QiiqsQTKE9LVq0kMDAQNXDpzn11FPVXD8Ev46gDeipQ6CKABAQ1FZNXnH1eTEH8corr1TzAAG9fAhMkZzjihUrVqhyRbbJRPv27at2u/3791t7ObX3Hm3q1KlTtdtqCUK4DxKX6OQwwYWIiPzWq6++qoZO0aP1+++/qyQQJE8giESP1LRp0+xujwBpxowZKihCJvTnn3+uklwAvW5vv/22SvZASZ4PP/xQBY8IhDD8jMALwedXX32lhoeRCIIkDPR+1Qb3Xbt2rWrPtddea9c7hkxiBFxIaFm3bp1KmPn222+tCS4IpjAsjt7Hv/76SwWNd9xxR7WkHleft0OHDqoHFc+NXlg8fkZGhkvvv/Y4COrQm4hhaAxHoxezKmRnjxw5Ug27L1++XPWkIiO66hC0ljCDzHMcBCDBBo+L1/HKK6+oy+QiJ3MZiYiI/MLevXtVckRycrIpJCTE1LJlS9N9991nys7OtrsdElymTp1quu6661TyREpKimn27NnW65Gw0bdvX1NsbKwpKirKdOaZZ5p++eUX6/VIfJk8ebIpNTVVPU/Tpk1NV111lWn9+vXWBJe4uDin7TzjjDNUYsivv/5a7Tokr1x44YWm6Oho9dw9evQwTZs2zXo9kkCQfIIElVatWpk++OAD9XpqSnCp7XmPHDmiEkrwnElJSabHH3/cNGLECGuSiZbggkSZ2hJcHnnkEVN8fLx6LCQNoV227wUSXHr27Gl67bXXTM2aNTOFh4ebrr32WlNOTo71NrYJLlBZWWmaNWuWqVOnTur9TkxMNA0ePNi0bNmyWl8z2QvAP4ywiYiIiMgRDkMTERERkVMMFomIiIjIKQaLREREROQUg0UiIiIicorBIhERERE5xWCRiIiIiJxisEhERERETjFYJCIiIiKnGCwSERERkVMMFomIiIjIKQaLRERERCTO/D9vX7epeZG6AAAAAABJRU5ErkJggg==", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "labels = mod.design_sal_labels_['covariate']\n", "n_sig = len(sig_lvs)\n", "\n", "fig = plt.figure(constrained_layout=True)\n", "subfigs = fig.subfigures(nrows=n_sig)\n", "if n_sig == 1:\n", " subfigs = [subfigs]\n", "\n", "for plot_idx, lv_idx in enumerate(sig_lvs):\n", " subfig = subfigs[plot_idx]\n", " subfig.suptitle('LV %s (%.2f%% variance explained, p = %.3f)' % (\n", " lv_idx,\n", " 100*mod.variance_explained_[lv_idx],\n", " mod.pvals_[lv_idx]))\n", " ax = subfig.subplots(ncols=2)\n", " ax[0].bar(x=labels,\n", " height=mod.boot_stat_val_[:, lv_idx],\n", " yerr=mod.get_boot_stat_yerr(lv_idx))\n", " ax[0].set_ylabel('Correlation with data score')\n", " ax[1].plot(mod.data_sals_z_[:, lv_idx])\n", " ax[1].set_ylabel('Bootstrap ratio (z score)')\n", " ax[1].set_xlabel('Observed variable')\n" ] }, { "cell_type": "markdown", "id": "1cf88077-bf52-45b1-ac7a-56e59d553437", "metadata": {}, "source": [ "As we can see, the model recovers the roughly sinusoidal pattern of loadings. `mod.flip_signs(lv_idx=0)` can be used to obtain an equivalent solution in which the covariates have positive correlations with the latent variable." ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.13.2" } }, "nbformat": 4, "nbformat_minor": 5 }