{ "cells": [ { "cell_type": "code", "execution_count": 67, "id": "geological-comfort", "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "/home/shellytrigg/bin/anaconda3/lib/python3.8/site-packages/pandas/core/frame.py:4441: SettingWithCopyWarning: \n", "A value is trying to be set on a copy of a slice from a DataFrame\n", "\n", "See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n", " return super().rename(\n" ] } ], "source": [ "# Import necessary libraries\n", "import pandas as pd\n", "import matplotlib.pyplot as plt\n", "import seaborn as sns\n", "import numpy as np\n", "from sklearn.decomposition import PCA\n", "from sklearn.preprocessing import StandardScaler\n", "from pycombat import Combat\n", "\n", "# Load the gene counts data\n", "gene_counts_path = '/mnt/c/Users/strigg/Downloads/salmon.merged.gene_counts (2).tsv'\n", "metadata_path = '../data/SraRunTable (3).csv'\n", "\n", "# Replace with the correct paths to your files\n", "gene_counts = pd.read_csv(gene_counts_path, sep='\\t')\n", "metadata = pd.read_csv(metadata_path)\n", "\n", "# Simplify metadata for alignment\n", "metadata_subset = metadata[['Experiment', 'treatment', 'Collection_Date', 'batch']]\n", "metadata_subset.rename(columns={'Experiment': 'Sample'}, inplace=True)\n", "\n", "# Align gene counts with metadata\n", "gene_count_samples = set(gene_counts.columns[2:])\n", "metadata_samples = set(metadata_subset['Sample'])\n", "aligned_gene_counts = gene_counts[['gene_id'] + list(gene_count_samples.intersection(metadata_samples))]\n", "\n", "# Filter genes with low expression (e.g., total counts < 10 across all samples)\n", "filtered_gene_counts = aligned_gene_counts.loc[\n", " aligned_gene_counts.iloc[:, 1:].sum(axis=1) > 10\n", "]\n", "\n", "# Calculate CPM normalization\n", "total_counts_per_sample = filtered_gene_counts.iloc[:, 1:].sum(axis=0)\n", "cpm_counts = filtered_gene_counts.iloc[:, 1:].div(total_counts_per_sample, axis=1) * 1e6\n", "cpm_counts.insert(0, 'gene_id', filtered_gene_counts['gene_id'])\n", "\n", "# Transpose CPM data for visualization\n", "cpm_values = cpm_counts.drop(columns='gene_id').T\n", "\n" ] }, { "cell_type": "code", "execution_count": 33, "id": "stunning-direction", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
0145679101112...38253382543825638257382583825938260382613826238263
SRX98455340.36961458.7686558.4663820.3271098.7480289.07698517.5566740.73922896.46967113.017257...11.0884251.65568722.3616583.3265280.00000025.5220441.1088435.35940636.961418225.279843
SRX98455410.1937691.60730340.2191771.74088710.6712408.87695416.7873830.53576857.5057188.825163...5.8934430.00000020.1805781.7858920.000000162.7611831.4287134.28614043.575761185.732753
SRX98455260.17178310.5202056.6020581.59494716.92063212.06323614.6063242.063937121.13722913.548794...7.3031620.31752935.8807522.4368750.317529121.2745602.6989955.08046028.577590197.502900
SRX98455300.0000001.2988834.4482401.83214617.5884598.68476314.5763501.010242113.0026505.647253...9.0921793.92969731.6061463.3424580.0000002.0023000.4329615.33985127.70949877.500001
SRX98455430.6737372.91952671.6624721.79663215.38702722.5706341.5720530.00000097.6922963.401922...8.3094211.97270247.6107391.5720530.000000286.7628502.9195264.04242185.789162375.271438
\n", "

5 rows × 32707 columns

\n", "
" ], "text/plain": [ " 0 1 4 5 6 7 \\\n", "SRX9845534 0.369614 58.768655 8.466382 0.327109 8.748028 9.076985 \n", "SRX9845541 0.193769 1.607303 40.219177 1.740887 10.671240 8.876954 \n", "SRX9845526 0.171783 10.520205 6.602058 1.594947 16.920632 12.063236 \n", "SRX9845530 0.000000 1.298883 4.448240 1.832146 17.588459 8.684763 \n", "SRX9845543 0.673737 2.919526 71.662472 1.796632 15.387027 22.570634 \n", "\n", " 9 10 11 12 ... 38253 \\\n", "SRX9845534 17.556674 0.739228 96.469671 13.017257 ... 11.088425 \n", "SRX9845541 16.787383 0.535768 57.505718 8.825163 ... 5.893443 \n", "SRX9845526 14.606324 2.063937 121.137229 13.548794 ... 7.303162 \n", "SRX9845530 14.576350 1.010242 113.002650 5.647253 ... 9.092179 \n", "SRX9845543 1.572053 0.000000 97.692296 3.401922 ... 8.309421 \n", "\n", " 38254 38256 38257 38258 38259 38260 \\\n", "SRX9845534 1.655687 22.361658 3.326528 0.000000 25.522044 1.108843 \n", "SRX9845541 0.000000 20.180578 1.785892 0.000000 162.761183 1.428713 \n", "SRX9845526 0.317529 35.880752 2.436875 0.317529 121.274560 2.698995 \n", "SRX9845530 3.929697 31.606146 3.342458 0.000000 2.002300 0.432961 \n", "SRX9845543 1.972702 47.610739 1.572053 0.000000 286.762850 2.919526 \n", "\n", " 38261 38262 38263 \n", "SRX9845534 5.359406 36.961418 225.279843 \n", "SRX9845541 4.286140 43.575761 185.732753 \n", "SRX9845526 5.080460 28.577590 197.502900 \n", "SRX9845530 5.339851 27.709498 77.500001 \n", "SRX9845543 4.042421 85.789162 375.271438 \n", "\n", "[5 rows x 32707 columns]" ] }, "execution_count": 33, "metadata": {}, "output_type": "execute_result" } ], "source": [ "cpm_values.head()" ] }, { "cell_type": "code", "execution_count": 35, "id": "chemical-seminar", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "( SRX9845534 SRX9845541 SRX9845526 SRX9845530 SRX9845543 \\\n", " count 32707.000000 32707.000000 32707.000000 32707.000000 32707.000000 \n", " mean 165.440053 171.200147 192.577787 211.851658 136.141405 \n", " std 1075.332847 957.385064 1048.747322 1345.628396 686.896627 \n", " min 0.000000 0.000000 0.000000 0.000000 0.000000 \n", " 25% 2.000000 2.000000 4.000000 2.489000 2.000000 \n", " 50% 23.895000 25.000000 33.037000 28.248000 20.406000 \n", " 75% 103.660000 110.792000 131.232500 126.000000 96.000000 \n", " max 91474.000000 85682.000000 87475.000000 137044.000000 69707.000000 \n", " \n", " SRX9845527 SRX9845506 SRX9845518 SRX9845501 \\\n", " count 32707.000000 32707.000000 32707.000000 32707.000000 \n", " mean 179.994930 183.099091 180.849091 159.277506 \n", " std 1453.555755 1201.931238 974.852911 911.736745 \n", " min 0.000000 0.000000 0.000000 0.000000 \n", " 25% 2.000000 2.000000 5.581500 2.403500 \n", " 50% 22.000000 25.318000 32.227000 25.640000 \n", " 75% 98.000000 112.560500 125.000000 111.000000 \n", " max 147753.000000 112650.000000 86589.000000 86268.000000 \n", " \n", " SRX9845536 ... SRX9845505 SRX9845531 SRX9845510 \\\n", " count 32707.000000 ... 32707.000000 32707.000000 32707.000000 \n", " mean 236.052346 ... 171.201700 206.463199 165.741786 \n", " std 1271.652215 ... 986.912211 1131.480546 1038.359313 \n", " min 0.000000 ... 0.000000 0.000000 0.000000 \n", " 25% 5.000000 ... 3.000000 2.000000 3.000000 \n", " 50% 43.000000 ... 27.522000 23.721000 28.000000 \n", " 75% 172.000000 ... 113.682000 125.000000 119.000000 \n", " max 100197.000000 ... 97256.000000 98612.000000 93615.000000 \n", " \n", " SRX9845535 SRX9845513 SRX9845512 SRX9845525 \\\n", " count 32707.000000 32707.000000 32707.000000 32707.000000 \n", " mean 199.503417 191.538372 131.001250 230.778746 \n", " std 1311.955681 1232.982014 721.962074 1613.978983 \n", " min 0.000000 0.000000 0.000000 0.000000 \n", " 25% 3.000000 4.000000 2.000000 3.000000 \n", " 50% 28.601000 30.000000 20.700000 30.670000 \n", " 75% 118.922500 121.000000 93.943500 139.000000 \n", " max 130487.000000 113471.000000 77585.000000 145177.000000 \n", " \n", " SRX9845521 SRX9845500 SRX9845511 \n", " count 32707.000000 32707.000000 32707.000000 \n", " mean 210.007663 148.738319 129.222067 \n", " std 1196.141346 959.965734 688.428395 \n", " min 0.000000 0.000000 0.000000 \n", " 25% 3.000000 1.000000 2.000000 \n", " 50% 31.306000 18.320000 20.000000 \n", " 75% 134.000000 88.000000 90.000000 \n", " max 116651.000000 95149.000000 67648.000000 \n", " \n", " [8 rows x 44 columns],\n", " gene_id LOC111099029LOC111099030LOC111099033LOC1110990...\n", " SRX9845534 5411047.807\n", " SRX9845541 5599443.202\n", " SRX9845526 6298641.692\n", " SRX9845530 6929032.179\n", " SRX9845543 4452776.927\n", " SRX9845527 5887094.163\n", " SRX9845506 5988621.973\n", " SRX9845518 5915031.218\n", " SRX9845501 5209489.38\n", " SRX9845536 7720564.083\n", " SRX9845529 6904563.668\n", " SRX9845524 8228094.16\n", " SRX9845507 5451326.526\n", " SRX9845540 5058340.811\n", " SRX9845517 5033662.657\n", " SRX9845532 7787956.536\n", " SRX9845503 5090347.989\n", " SRX9845508 5380343.456\n", " SRX9845537 7915664.803\n", " SRX9845520 6881819.531\n", " SRX9845523 4540587.389\n", " SRX9845528 6565149.428\n", " SRX9845516 6858787.739\n", " SRX9845514 8076393.261\n", " SRX9845522 9031287.071\n", " SRX9845502 4523632.462\n", " SRX9845509 4345607.019\n", " SRX9845539 5521862.433\n", " SRX9845542 4701381.722\n", " SRX9845504 4731008.048\n", " SRX9845538 7748085.775\n", " SRX9845515 8087706.403\n", " SRX9845519 5850139.193\n", " SRX9845533 8105079.858\n", " SRX9845505 5599493.994\n", " SRX9845531 6752791.859\n", " SRX9845510 5420916.597\n", " SRX9845535 6525158.253\n", " SRX9845513 6264645.549\n", " SRX9845512 4284657.883\n", " SRX9845525 7548080.453\n", " SRX9845521 6868720.618\n", " SRX9845500 4864784.191\n", " SRX9845511 4226466.132\n", " dtype: object)" ] }, "execution_count": 35, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Summary statistics for the raw (unnormalized) gene count data\n", "raw_summary = filtered_gene_counts.describe()\n", "\n", "# Calculate library sizes (total counts per sample)\n", "library_sizes = filtered_gene_counts.sum(axis=0)\n", "\n", "# Display key statistics\n", "raw_summary, library_sizes" ] }, { "cell_type": "code", "execution_count": 41, "id": "divided-commander", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAABDAAAAFgCAYAAABNIolGAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAAAuGUlEQVR4nO3de7xcZX33/c8XAoKEkxJBY0IsHhAtCE9AK1SLKBVbi4dWQY1g9abe9QC3+FRLrQ8eeiu9lVtbtYiChwiCCiieoVSraEVCGgkQUEQUDEKAhCTKweDv+WOtDcM4e2eHZPasJJ/367VfmVnrWtf6zZoDzHeuda1UFZIkSZIkSV22xagLkCRJkiRJWhsDDEmSJEmS1HkGGJIkSZIkqfMMMCRJkiRJUucZYEiSJEmSpM4zwJAkSZIkSZ1ngCFJm6EkpyT5xw3U1+wkq5Ns2d7/dpLXbIi+2/6+nuSoDdXfOuz33UluTfKrqd73pqDLx29Dv0alXkn+JMmNo65DkjZFBhiStIlJcn2SO5OsSrIiyfeTvDbJfZ/5VfXaqnrXJPt69kRtquoXVTW9qu7dALWfmOQzff0fVlWfWt++17GOWcDxwF5VtduA9UP5gpLkXUkWJ1mT5MQB61+W5OdJfp3ki0keNqDNKW2gtDrJPUl+23P/6xu65nEex4THTw9OkkcmOS3JTe37++ok70iyXbu+2tfG6iS/THJyT7B4fft62KWvz0XtdnPWsu9Ptq/LRw3tAa6nJI9Ock4bnN3RvpeOHnVdkqQNxwBDkjZNz6+q7YHdgfcCbwFO29A7STJtQ/fZEbsDt1XVLVO832uBvwO+2r8iyZOAjwLzgF2B3wAf6W/XhlPTq2o68L+Bs8fuV9VhPf0N87l70MdvE35NrZc2rPovYFvgj9r393OAnYA9epru0z73hwAvA/5Hz7qfAUf29PmHbX9r2/d2wIuBO4CXP8j6p+J5nQ/cQPP6ezjwSuDmKdivJGmKGGBI0iasqu6oqvOBlwJHJXky3Pdr6rvb27sk+Uo7WuP2JN9NskWS+cBs4MvtL7p/l2RO+2vtq5P8AviPnmW9X1D2SPLD9lfQL42NFBg0cmFslEeS5wInAC9t9/ejdv19w/3but7WjkK4Jcmnk+zYrhur46gkv2h/hf2H8Y5Nkh3b7Ze1/b2t7f/ZwIXAo9o6PrkuxzzJE9uaVyS5Mslf9Kx7eJIvJ1mZ5NI0p1lc3PN8faqqvg6sGtD1y4EvV9V3qmo18I/Ai5Jsvw61XZ/kLUkuB36dZFqStyb5afuL/lVJXtjT/ugkFyd5X5LlSX6W5LC+9de12/4sycvHO35J/qI9Hiva4/PECep6bPtcvirJDe2+X5tk/ySXt318qO+x/XWSJW3bbybZvWfdc9KMVrij3S4THKMTk3w+yWfax7U4yeOT/H37mrshyaE97XfM/aMiftk+p2OjHvZI8h9Jbmtfj2ck2anvcb+5fUx3JDk7yTbjlPYmmtfFK6rqeoCquqGqjq2qy/sbV9XVwHeBJ/csnk/zpX7MUcCnxzsWPV4MrADe2W7Te7weluQTSZa2x/6L7fI/SXJj+7z+CvhEkock+UDbdml7+yFt+4GfQ+26t7THdlWSa5IcMk6d+wOfrKpfV9Waqvrv9v00Vuvnk/yqPdbfSRMKjq37ZJKPpDllbXWS7yXZra1xefv62ben/fXta+Kqdv0nxnvukjwqzciQZe375I096w5IsiDNZ8LNSU6exPMhSZstAwxJ2gxU1Q+BG4E/HrD6+HbdDJpf9k9oNql5wC9oRnNMr6p/7tnmmcATgT8dZ5evBP4aeBSwBviXSdT4DR44YmCfAc2Obv8OBv4AmA58qK/NQcATaH6BfnvvF+U+/wrs2PbzzLbmV1XVvwOHAUvbOo5eW+1jkmwFfBm4AHgE8AbgjCRPaJt8GPg1sBvNF8F1mdvjScCPxu5U1U+Be4DHr0Mf0PwC/2fATlW1BvgpzetiR+AdwGeSPLKn/VOBa4BdgH8GTktjO5rn9bB2NMDTgUWDjl+SxwOfBY6jeZ19jSYY23pQXTSvmbF9P44mgPsA8A/As9tj8ZIkzwRI8gKa1+2L2v6/2+6PNKdMnAO8rX0MPwUOXMsxej7Nl/2dgf8Gvknz/0wzab7Ef7Sn7afaeh8L7AscCozNrxHgPTTvgycCs4AT+/b1EuC5wGOAvWle34M8Gzi3qn63ltqbHSd70Tyv/92z+AfADmlCti1pjutnBm3f5yia43kWsGeS/XrWzQceSvOcPAL4vz3rdgMeRjMi4hia5+9pwFOAfYADaJ4XGOdzqH3vvB7Yv32d/Slw/Th1/gD4cJIjkswesP7rNK+nRwALgTP61r+E+18nd9OMeFnY3v8C0B8uvLytZw+a9+Hb+tbThjBfpnnvzqT5XDouydhn5weBD1bVDm0/nxvnsUmSMMCQpM3JUpovE/1+CzwS2L2qfltV362qWktfJ7a/ct45zvr5VXVFVf2aZqTAS8Z+lV5PLwdOrqrr2lEIfw8ckQeO/nhHVd1ZVT+i+dLwe0FIz5e3v6+qVe0v2u+nOT1jfTyNJlR5b1XdU1X/AXwFOLLd54uB/6+qflNVV9F8+Z2s6TRD+HvdAUx6BEbrX9pf7u8EqKrPV9XSqvpdVZ0N/ITmi+WYn1fVx9o5Tj5F81rZtV33O+DJSbatqpuq6spx9vlS4KtVdWFV/RZ4H82pC08fr67Wu6rqrqq6gCb4+WxV3VJVv6QJKcZ+Ef8b4D1VtaQNZf438JR2FMbzgKuq6gvtvj8ArG1i0e9W1Tfbvj5P86X6ve32ZwFzkuyUZFeasOa49v1wC80X+CPaY3tt+5jvrqplNF+An9m3r39pj//tNF90nzJOTQ8HblpL3QALkyxv+/o48Im+9WOjMJ4DXA38cqLO2iDgYODMqroZuIg2eGuDrsOA11bV8vbz4z97Nv8dzev97vZ5fTnwzvY5XEYTmI2958b7HLoXeAiwV5Ktqur6Nrwb5K9oXhf/CPwszfwe+4+trKrT2/f73TRB0j5pR3C1zquqy6rqLuA84K6q+nT72j+b+19vYz7UvmZvB/6JntNzeuwPzKiqd7afCdcBH6N9jbSP+7FJdqmq1VX1g3EemyQJAwxJ2pzMBG4fsPz/0My9cEGa0wHeOom+bliH9T8HtqL5FXN9Partr7fvadz/hRoe+OX0NzRf/PvtAmw9oK+ZG6C+G/p+JR/rd0Zba++xWdtx7LUa2KFv2Q7AqiR/nPsn6hwvRBi4zySvbL/orUiyguaUg97n6r7jWVW/aW9Ob8OplwKvBW5K8tUke46zzwc8b+3xuYEHHu9Bx6J3/oI7B9wfe253Bz7Y8xhupxn9MLPd9319t1+K13bc+/dza90/Se1YwDK93e9WNI9/bN8fpfmFnySPSHJWe/rDSprRDv3vg8m8XgFuo/mCvzb7VdXOVbVHVb1twIiN+TRzYxzN5E4fmQcsqapF7f0zgJe1o41mAbdX1fJxtl3WhgFjBr1/xyYFHfg5VFXX0ozcORG4pT2eAycSbUOUt1bVk2g+ExYBX2xHDG2Z5L1pTpdayf2jOHqfj8m+3sb0f84Nqmt3mtOpVvS8Rk7g/s+sV9OM3rg6zWllfz7osUmSGgYYkrQZaH+FnAlc3L+u/UXy+Kr6A5qh82/qOcd8vJEYaxuhMavn9myaXxlvpfkV/aE9dW1J88V+sv0upflC0Nv3GtZ9or5b25r6+5rw1+hJWArMSs8VX3r6XUZT66N71vUep7W5kp7RJEn+gOaX6R+3v1aPTdT5pHF7aNx3jNsRCh+jGaL/8KraCbiCCeaIeEBHzSiF59B8sb667WuQBzxvSULz2HuP99qe+4ncAPxNVe3U87dtVX2fZtTCfce5Z98bwg00pxrs0rPfHXqeg/fQPK6921MEXsEkj+0A/w68sO+1tc6q6uc0k3k+Dzh3Epu8EviDdu6IX9GMItmFZuTFDcDD0jOvR//u+u4Pev8ubesa93Ooqs6sqoPabQs4aRKP81aakT6Pohl59jLgcJpTcXYE5rRNH+zzAb//Obd0QJsbgJ/1vTa3r6rntXX+pKqOpAm9TgK+0J6eJUkawABDkjZhSXZof9E7C/hMVS0e0ObP00yaGGAlzZDtsV+bb6aZI2JdvSLJXkkeSjNnwBfaX7B/DGyT5M/aX3DfRvMlfMzNNMPzx/vv02eB/5XkMUl6r7KxZpz2A7W1fA74pyTbt1/k38Tk5gO4T5Jtev+AH9KENH+XZKskf0LzZeysdp/nAicmeWg7WuGVff1t1fazBTCt7Xfs1JszgOe3oy22ozmu51bVoAk/J2s7mi+Ey9r9v4oHTvo40WPfNc3EnNvRfIlfzf2vm36fA/4sySHt8358u83316P2XqcAf592UsY0E2v+Vbvuq8CTkryoPdXojTRzM6y3qrqJZr6T97fvtS3STNw5dprI9jTHZUWSmcD/ux67O5lmxM2n2tcrSWamuVTq3uvY16uBZ7WjaMaV5I9o5mU4gObUlqfQvD7OBI5qH//XgY8k2bl9/T5jgi4/C7wtyYw0c5O8nfY9N97nUJInJHlWmsk+76IZCTHwdZbkpCRPTjM57fbA/wSurarbaJ6Lu2lGsjyU5rNjfb0uzaVbH0YzquLsAW1+CKxMMxHptu1IkCe3oTJJXpFkRjtSZkW7zXpfklqSNlUGGJK0afpyklU0v/79A82Xn1eN0/ZxNL/urqaZtO4jVfXtdt17aL5wrEjy5nXY/3zgkzTD47eh+dJIVd0B/C3Nufm/pPmy33tVks+3/96WZOGAfk9v+/4Oza/Id9FMlPlgvKHd/3U0I1PObPufrJk0X6Z6/2YBf0Hz6/StNJc5fWU1V4SAZqTDjjTHZT7NF7q7e/r8WNvPkTTP2520cwS080u8libIuIXmC9nfrkO9v6eaeTjeT/O83wz8IfC9SW6+BU0QsZTmlI1njldPVV1DM/rgX2mOy/NpJoe9Z33q7+n/PJpfr89qTw+4guY5GPsl/q9oLid8G83rfbKPcTJeSXM60lXAcprJHsdO9XgHsB/NXCVfZXIjHgZq51l4Os3IoUva9/dFbd/XrmNfP62qBZNoehTwpapaXFW/GvujmXjyz9sv7vPamq6meV0eN0F/7wYWAJcDi2kmyHx3u268z6GH0Dx3t9K8bx5BExYM8lCauStW0Lyvd6d5P0JzuszPaT53rqKZ8HN9nUkTYF3X/r27v0EbXD6fJvz5Wfs4Pk7zOQDNBK5XJllNc1yP6DvtRpLUI7XWedokSdIwJDkJ2K2q1uVqJJJGLMn1wGuqueqOJGmKOAJDkqQpkmTPJHu3kwoeQDOU/7xR1yVJkrQxmLb2JpIkaQPZnua0kUfRDLd/P/ClkVYkSZK0kfAUEkmSJEmS1HmeQiJJkiRJkjpvkzqFZJdddqk5c+aMugxJkiRJkvQgXXbZZbdW1Yz+5ZtUgDFnzhwWLJjMVcEkSZIkSVIXJfn5oOWeQiJJkiRJkjrPAEOSJEmSJHXe0AKMJLOSfCvJkiRXJjl2grb7J7k3yV/2LLs+yeIki5J4XogkSZIkSZuxYc6BsQY4vqoWJtkeuCzJhVV1VW+jJFsCJwHfHNDHwVV16xBrlCRJkiRJG4GhjcCoqpuqamF7exWwBJg5oOkbgHOAW4ZViyRJkiRJ2rhNyRwYSeYA+wKX9C2fCbwQOGXAZgVckOSyJMdM0PcxSRYkWbBs2bINWLUkSZIkSeqKoQcYSabTjLA4rqpW9q3+APCWqrp3wKYHVtV+wGHA65I8Y1D/VXVqVc2tqrkzZvzeZWIlSZIkSdImYJhzYJBkK5rw4oyqOndAk7nAWUkAdgGel2RNVX2xqpYCVNUtSc4DDgC+M8x6JUmSJElSNw0twEiTSpwGLKmqkwe1qarH9LT/JPCVqvpiku2ALapqVXv7UOCdw6pVkiRJkiR12zBHYBwIzAMWJ1nULjsBmA1QVYPmvRizK3BeOzJjGnBmVX1jeKVKkiRJkqQuG1qAUVUXA1mH9kf33L4O2GcIZUlS51xwwQXcfPPNoy5D6rzly5cDsPPOO4+4Eqn7dt11Vw499NBRlyFJG9RQ58CQJEnaUO65555RlyBJkkbIAEOSRsxfyKTJmT9/PgDz5s0bcSWSJGkUhn4ZVUmSJEmSpPVlgCFJkiRJkjrPAEOSJEmSJHWeAYYkSZIkSeo8AwxJkiRJktR5BhiSJEmSJKnzDDAkSZIkSVLnGWBIkiRJkqTOM8CQJEmSJEmdZ4AhSZIkSZI6zwBDkiRJkiR1ngGGJEmSJEnqPAMMSZIkSZLUeQYYkiRJkiSp8wwwJEmSJElS5xlgSJIkSZKkzjPAkCRJkiRJnWeAIUmSJEmSOs8AQ5IkSZIkdZ4BhiRJkiRJ6jwDDEmSJEmS1HkGGJIkSZIkqfMMMCRJkiRJUucZYEiSJEmSpM4zwJAkSZIkSZ1ngCFJkiRJkjrPAEOSJEmSJHXe0AKMJLOSfCvJkiRXJjl2grb7J7k3yV/2LHtukmuSXJvkrcOqU5IkSZIkdd8wR2CsAY6vqicCTwNel2Sv/kZJtgROAr7Zt+zDwGHAXsCRg7aVJEmSJEmbh6EFGFV1U1UtbG+vApYAMwc0fQNwDnBLz7IDgGur6rqqugc4Czh8WLVKkiRJkqRum5I5MJLMAfYFLulbPhN4IXBK3yYzgRt67t/I4PCDJMckWZBkwbJlyzZYzZIkSZIkqTuGHmAkmU4zwuK4qlrZt/oDwFuq6t7+zQZ0VYP6r6pTq2puVc2dMWPGetcrSZIkSZK6Z9owO0+yFU14cUZVnTugyVzgrCQAuwDPS7KGZsTFrJ52jwaWDrNWSZIkSZLUXUMLMNKkEqcBS6rq5EFtquoxPe0/CXylqr6YZBrwuCSPAX4JHAG8bFi1SpIkSZKkbhvmCIwDgXnA4iSL2mUnALMBqqp/3ov7VNWaJK+nuTLJlsDpVXXlEGuVJEmSJEkdNrQAo6ouZvBcFuO1P7rv/teAr23gsiRJkiRJ0kZoSq5CIkmSJEmStD4MMCRJkiRJUucZYEiSJEmSpM4zwJAkSZIkSZ1ngCFJkiRJkjrPAEOSJEmSJHWeAYYkSZIkSeo8AwxJkiRJktR5BhiSJEmSJKnzDDAkSZIkSVLnGWBIkiRJkqTOM8CQJEmSJEmdZ4AhSZIkSZI6zwBDkiRJkiR1ngGGJEmSJEnqPAMMSZIkSZLUeQYYkiRJkiSp8wwwJEmSJElS5xlgSJIkSZKkzjPAkCRJkiRJnWeAIUmSJEmSOs8AQ5IkSZIkdZ4BhiRJkiRJ6jwDDEmSJEmS1HkGGJIkSZIkqfMMMCRJkiRJUucZYEiSJEmSpM4zwJAkSZIkSZ1ngCFJkiRJkjpvaAFGkllJvpVkSZIrkxw7oM3hSS5PsijJgiQH9ay7PsnisXXDqlOSJEmSJHXftCH2vQY4vqoWJtkeuCzJhVV1VU+bi4Dzq6qS7A18DtizZ/3BVXXrEGuUJEmSJEkbgaGNwKiqm6pqYXt7FbAEmNnXZnVVVXt3O6CQJEmSJEnqMyVzYCSZA+wLXDJg3QuTXA18FfjrnlUFXJDksiTHTND3Me3pJwuWLVu2gSuXJEmSJEldMPQAI8l04BzguKpa2b++qs6rqj2BFwDv6ll1YFXtBxwGvC7JMwb1X1WnVtXcqpo7Y8aMDf8AJEmSJEnSyA01wEiyFU14cUZVnTtR26r6DrBHkl3a+0vbf28BzgMOGGatkiRJkiSpu4Z5FZIApwFLqurkcdo8tm1Hkv2ArYHbkmzXTvxJku2AQ4ErhlWrJEmSJEnqtmFeheRAYB6wOMmidtkJwGyAqjoFeDHwyiS/Be4EXtpekWRX4Lw225gGnFlV3xhirZIkSZIkqcOGFmBU1cVA1tLmJOCkAcuvA/YZUmmSJEmSJGkjMyVXIZEkSZIkSVofBhiSJEmSJKnzDDAkSZIkSVLnGWBIkiRJkqTOM8CQJEmSJEmdZ4AhSZIkSZI6zwBDkiRJkiR1ngGGJEmSJEnqPAMMSZIkSZLUeQYYkiRJkiSp86ZNtDLJHwGvAP4YeCRwJ3AF8FXgM1V1x9Ar1Ebrggsu4Oabbx51GZKkTcTYf1Pmz58/4kokSZuKXXfdlUMPPXTUZWiSxg0wknwdWAp8Cfgn4BZgG+DxwMHAl5KcXFXnT0Wh2vjcfPPN3HTTTey8006jLkWStAnYYotm4Ohdd9454kokSZuC5StWjLoEraOJRmDMq6pb+5atBha2f+9PssvQKtMmYeedduKQQw4ZdRmSJEmS9AAXXXTRqEvQOhp3DowB4cWDaiNJkiRJkrS+HtQknkkWb+hCJEmSJEmSxjPRHBgvGm8VsNtwypEkSZIkSfp9E82BcTZwBlAD1m0znHIkSZIkSZJ+30QBxuXA+6rqiv4VSZ49vJIkSZIkSZIeaKI5MI4DVo6z7oUbvhRJkiRJkqTBxh2BUVXfnWDdguGUI0mSJEmS9Pse1FVIJEmSJEmSppIBhiRJkiRJ6jwDDEmSJEmS1HmTDjCSPGeYhUiSJEmSJI1nXUZgnDS0KiRJkiRJkibgKSSSJEmSJKnzxr2MKkCSTwAFBJid5PSxdVX110OuTZIkSZIkCVhLgAF8suf2QcCnhleKJEmSJEnSYBMGGFX1n2O3k6zqvS9JkiRJkjRV1mUOjHvWpeMks5J8K8mSJFcmOXZAm8OTXJ5kUZIFSQ7qWffcJNckuTbJW9dl35IkSZIkadOytlNI7lNVT1vHvtcAx1fVwiTbA5clubCqruppcxFwflVVkr2BzwF7JtkS+DDwHOBG4NIk5/dtK0mSJEmSNhNDuwpJVd1UVQvb26uAJcDMvjarq6rau9vRTBgKcABwbVVdV1X3AGcBhw+rVkmSJEmS1G1TchnVJHOAfYFLBqx7YZKrga8CY1c2mQnc0NPsRvrCj57tj2lPP1mwbNmyDVq3JEmSJEnqhqEHGEmmA+cAx1XVyv71VXVeVe0JvAB419hmA7qqAcuoqlOram5VzZ0xY8YGqlqSJEmSJHXJUAOMJFvRhBdnVNW5E7Wtqu8AeyTZhWbExaye1Y8Glg6tUEmSJEmS1GnjTuKZZBUPHPWQ9n6AqqodJuo4SYDTgCVVdfI4bR4L/LSdxHM/YGvgNmAF8LgkjwF+CRwBvGyyD0qSJEmSJG1aJroKyUXAbsC5wFlV9Yt17PtAYB6wOMmidtkJwGyAqjoFeDHwyiS/Be4EXtpO6rkmyeuBbwJbAqdX1ZXruH9JkiRJkrSJGDfAqKoXJNkReBHwsSTbAGfThBm3r63jqrqYwXNZ9LY5CThpnHVfA762tv1IkiRJkqRN34RzYFTVHVX1CeAw4BTgncDRU1CXJEmSJEnSfSY6hYQkTweOBP4YuBh4YVV9dyoKkyRJkiRJGjPRJJ7X00ymeRZwDLCmXb4fQFUtHH55kiRJkiRJE4/AuJ7mqiN/2v71KuBZQ6pJkiRJkiTpASaaxPNPprAOSZIkSZKkcY07iWeSxyX5YpIrknw2ycypLEySJEmSJGnMRFchOR34KvBiYCHwr1NSkSRJkiRJUp+J5sDYvqo+1t7+P0mctFOSJEmSJI3ERAHGNkn2BdLe37b3vlchkSRJkiRJU2WiAONXwMnj3PcqJJIkSZIkacp4FRJJkiRJktR5E12F5BVJ5g1Y/j+SvGy4ZUmSJEmSJN1voquQHA98ccDys9t1kiRJkiRJU2KiAGPLqlrVv7CqVgJbDa8kSZIkSZKkB5oowNgqyXb9C5NsD2w9vJIkSZIkSZIeaKIA4zTgC0nmjC1ob5/VrpMkSZIkSZoSE12F5H1JVgP/mWQ6zaVTfw28t6r+baoKlCRJkiRJGjfAAKiqU4BT2gAjg+bEkCRJkiRJGra1XUZ1C4CqWt0fXiTZI8lBwy5QkiRJkiRpohEYDwf+O8llwGXAMmAb4LHAM4FbgbcOvUJJkiRJkrTZm2gOjA8m+RDwLOBAYG/gTmAJMK+qfjE1JUqSJEmSpM3d2ubAuBe4sP2TJEmSJEkaiQkDjCR/CrwAmElzFZKlwJeq6hvDL02SJEmSJKkxboCR5APA44FPAze2ix8NvDHJYVV17PDLkyRJkiRJmngExvOq6vH9C5OcDfwYMMDQhJYvX85dd93FRRddNOpSJEmSJOkBlq9YwTZ33TXqMrQOxr2MKnBXkgMGLN8f8FmWJEmSJElTZqIRGEcD/5Zke+4/hWQWsLJdJ01o55135q477+SQQw4ZdSmSJEmS9AAXXXQR22y77ajL0DqY6DKqC4GnJtmNZhLPADdW1a+mqjhJkiRJkiRYy1VIANrA4gGhRZI9q+rqibZLMotmAtDdgN8Bp1bVB/vavBx4S3t3NfA/q+pH7brrgVXAvcCaqpo7mQckSZIkSZI2PWsNMMZxATB7LW3WAMdX1cL2NJTLklxYVVf1tPkZ8MyqWp7kMOBU4Kk96w+uqlsfZI2SJEmSJGkTMdFlVP9lvFXATmvruKpuAm5qb69KsoTmVJSretp8v2eTH9BcplWSJEmSJOkBJhqB8SrgeODuAeuOXJedJJkD7AtcMkGzVwNf77lfwAVJCvhoVZ26LvuUJEmSJEmbjokCjEuBK/pGSQCQ5MTJ7iDJdOAc4LiqWjlOm4NpAoyDehYfWFVLkzwCuDDJ1VX1nQHbHgMcAzB79trOapEkSZIkSRujLSZY95fAokErquoxk+k8yVY04cUZVXXuOG32Bj4OHF5Vt/XsY2n77y3AecAB49RyalXNraq5M2bMmExZkiRJkiRpIzNugFFVt1fVbx5sx0kCnAYsqaqTx2kzGzgXmFdVP+5Zvl078SdJtgMOBa54sLVIkiRJkqSN21qvQpJkMc18FL3uABYA7+4dNdHnQGAesDjJonbZCbRXL6mqU4C3Aw8HPtLkHfddLnVX4Lx22TTgzKr6xuQfliRJkiRJ2pRM5jKqXwfuBc5s7x/R/rsS+CTw/EEbVdXFNFcsGVdVvQZ4zYDl1wH7TKI2SZIkSZK0GZhMgHFgVR3Yc39xku9V1YFJXjGswiRJkiRJksZMNInnmOlJnjp2J8kBwPT27pqhVCVJkiRJktRjMiMwXgOc3l4ONTSnjry6nVzzPcMsTpIkSZIkCSYRYFTVpcAfJtkRSFWt6Fn9uWEVJkmSJEmSNGatp5Ak2THJycBFwL8neX8bZkiSJEmSJE2JycyBcTqwCnhJ+7cS+MQwi5IkSZIkSeo1mTkw9qiqF/fcf0eSRUOqR5IkSZIk6fdMZgTGnUkOGruT5EDgzuGVJEmSJEmS9ECTGYHxWuDTPfNeLAeOGl5JkiRJkiRJDzSZq5D8CNgnyQ7t/ZVJjgMuH3JtkiRJkiRJwOROIQGa4KKqVrZ33zSkeiRJkiRJkn7PpAOMPtmgVUiSJEmSJE3gwQYYtUGrkCRJkiRJmsC4c2AkWcXgoCLAtkOrSJIkSZIkqc+4AUZVbT+VhUiSJEmSJI3nwZ5CIkmSJEmSNGUMMCRJkiRJUucZYEiSJEmSpM4zwJAkSZIkSZ1ngCFJkiRJkjrPAEOSJEmSJHWeAYYkSZIkSeo8AwxJkiRJktR5BhiSJEmSJKnzDDAkSZIkSVLnGWBIkiRJkqTOM8CQJEmSJEmdZ4AhSZIkSZI6zwBDkiRJkiR13tACjCSzknwryZIkVyY5dkCblye5vP37fpJ9etY9N8k1Sa5N8tZh1SlJkiRJkrpv2hD7XgMcX1ULk2wPXJbkwqq6qqfNz4BnVtXyJIcBpwJPTbIl8GHgOcCNwKVJzu/bVpIkSZIkbSaGNgKjqm6qqoXt7VXAEmBmX5vvV9Xy9u4PgEe3tw8Arq2q66rqHuAs4PBh1SpJkiRJkrptSubASDIH2Be4ZIJmrwa+3t6eCdzQs+5G+sIPSZIkSZK0+RjmKSQAJJkOnAMcV1Urx2lzME2AcdDYogHNapxtjwGOAZg9e/Z61ytJkiRJkrpnqCMwkmxFE16cUVXnjtNmb+DjwOFVdVu7+EZgVk+zRwNLB21fVadW1dyqmjtjxowNV7wkSZIkSeqMYV6FJMBpwJKqOnmcNrOBc4F5VfXjnlWXAo9L8pgkWwNHAOcPq1ZJkiRJktRtwzyF5EBgHrA4yaJ22QnAbICqOgV4O/Bw4CNN3sGadjTFmiSvB74JbAmcXlVXDrFWDcnyFSu46KKLRl2GJGkTsGr1agC2nz59xJVIkjYFy1es4JHbbjvqMrQOhhZgVNXFDJ7LorfNa4DXjLPua8DXhlCapsiuu+466hIkSZuQO1Y2U2lt4/9sSpI2gEduu63fWTYyQ5/EU5uvQw89dNQlSJI2IfPnzwdg3rx5I65EkiSNwpRcRlWSJEmSJGl9GGBIkiRJkqTOM8CQJEmSJEmdZ4AhSZIkSZI6zwBDkiRJkiR1ngGGJEmSJEnqPAMMSZIkSZLUeQYYkiRJkiSp8wwwJEmSJElS5xlgSJIkSZKkzjPAkCRJkiRJnWeAIUmSJEmSOs8AQ5IkSZIkdZ4BhiRJkiRJ6jwDDEmSJEmS1HkGGJIkSZIkqfMMMCRJkiRJUucZYEiSJEmSpM4zwJAkSZIkSZ1ngCFJkiRJkjrPAEOSJEmSJHWeAYYkSZIkSeo8AwxJkiRJktR5BhiSJEmSJKnzDDAkSZIkSVLnGWBIkiRJkqTOM8CQJEmSJEmdZ4AhSZIkSZI6b2gBRpJZSb6VZEmSK5McO6DNnkn+K8ndSd7ct+76JIuTLEqyYFh1SpIkSZKk7ps2xL7XAMdX1cIk2wOXJbmwqq7qaXM78EbgBeP0cXBV3TrEGiVJkiRJ0kZgaCMwquqmqlrY3l4FLAFm9rW5paouBX47rDokSZIkSdLGb0rmwEgyB9gXuGQdNivggiSXJTlmgr6PSbIgyYJly5atZ6WSJEmSJKmLhh5gJJkOnAMcV1Ur12HTA6tqP+Aw4HVJnjGoUVWdWlVzq2rujBkzNkDFkiRJkiSpa4YaYCTZiia8OKOqzl2XbatqafvvLcB5wAEbvkJJkiRJkrQxGOZVSAKcBiypqpPXcdvt2ok/SbIdcChwxYavUpIkSZIkbQyGeRWSA4F5wOIki9plJwCzAarqlCS7AQuAHYDfJTkO2AvYBTivyUCYBpxZVd8YYq2SJEmSJKnDhhZgVNXFQNbS5lfAowesWgnsM4y6JEmSJEnSxmdKrkIiSZIkSZK0PgwwJEmSJElS5xlgSJIkSZKkzjPAkCRJkiRJnWeAIUmSJEmSOs8AQ5IkSZIkdZ4BhiRJkiRJ6jwDDEmSJEmS1HkGGJIkSZIkqfMMMCRJkiRJUucZYEiSJEmSpM4zwJAkSZIkSZ1ngCFJkiRJkjrPAEOSJEmSJHWeAYYkSZIkSeo8AwxJkiRJktR5BhiSJEmSJKnzDDAkSZIkSVLnGWBIkiRJkqTOM8CQJEmSJEmdN23UBUjS5u6CCy7g5ptvHnUZUueNvU/mz58/4kqk7tt111059NBDR12GJG1QBhiSJGmjsPXWW4+6BEmSNEIGGJI0Yv5CJkmSJK2dc2BIkiRJkqTOM8CQJEmSJEmdZ4AhSZIkSZI6zwBDkiRJkiR1ngGGJEmSJEnqPAMMSZIkSZLUeUMLMJLMSvKtJEuSXJnk2AFt9kzyX0nuTvLmvnXPTXJNkmuTvHVYdUqSJEmSpO6bNsS+1wDHV9XCJNsDlyW5sKqu6mlzO/BG4AW9GybZEvgw8BzgRuDSJOf3bStJkiRJkjYTQxuBUVU3VdXC9vYqYAkws6/NLVV1KfDbvs0PAK6tquuq6h7gLODwYdUqSZIkSZK6bUrmwEgyB9gXuGSSm8wEbui5fyN94UdP38ckWZBkwbJly9arTkmSJEmS1E1DDzCSTAfOAY6rqpWT3WzAshrUsKpOraq5VTV3xowZD7ZMSZIkSZLUYcOcA4MkW9GEF2dU1bnrsOmNwKye+48Glq5to8suu+zWJD9ftyolSdJGZBfg1lEXIUmShmr3QQuHFmAkCXAasKSqTl7HzS8FHpfkMcAvgSOAl61to6pyCIYkSZuwJAuqau6o65AkSVNvmCMwDgTmAYuTLGqXnQDMBqiqU5LsBiwAdgB+l+Q4YK+qWpnk9cA3gS2B06vqyiHWKkmSJEmSOixVA6eWkCRJ6hxHYEiStPmakquQSJIkbSCnjroASZI0Go7AkCRJkiRJnecIDEmSJEmS1HkGGJIkSZIkqfMMMCRJkiRJUucZYEiSJEmSpM6bNuoCJEmS1ibJ04E59Py/S1V9emQFSZKkKWeAIUmSOi3JfGAPYBFwb7u4AAMMSZI2I15GVZIkdVqSJcBe5f+0SJK0WXMODEmS1HVXALuNughJkjRankIiSZK6bhfgqiQ/BO4eW1hVfzG6kiRJ0lQzwJAkSV134qgLkCRJo+ccGJIkSZIkqfOcA0OSJHVakqcluTTJ6iT3JLk3ycpR1yVJkqaWAYYkSeq6DwFHAj8BtgVe0y6TJEmbEefAkCRJnVdV1ybZsqruBT6R5PujrkmSJE0tAwxJktR1v0myNbAoyT8DNwHbjbgmSZI0xTyFRJIkdd08mv9neT3wa2AW8OKRViRJkqacVyGRJEmdl2RbYHZVXTPqWiRJ0mg4AkOSJHVakucDi4BvtPefkuT8kRYlSZKmnAGGJEnquhOBA4AVAFW1CJgzsmokSdJIGGBIkqSuW1NVd4y6CEmSNFpehUSSJHXdFUleBmyZ5HHAGwEvoypJ0mbGERiSJKnr3gA8CbgbOBO4Azh2pBVJkqQpZ4AhSZK6bq/2bxqwDXA4cOlIK5IkSVPOy6hKkqROS3IN8GbgCuB3Y8ur6ucjK0qSJE0558CQJEldt6yqvjzqIiRJ0mg5AkOSJHVakkOAI4GLaObBAKCqzh1ZUZIkaco5AkOSJHXdq4A9ga24/xSSAgwwJEnajBhgSJKkrtunqv5w1EVIkqTR8iokkiSp636QZK9RFyFJkkbLOTAkSVKnJVkC7AH8jGYOjABVVXuPtDBJkjSlDDAkSVKnJdl90HIvoypJ0ubFAEOSJEmSJHWec2BIkiRJkqTOM8CQJEmSJEmdZ4AhSZKGJsk/JLkyyeVJFiV56hD39e0kc4fVvyRJGq1poy5AkiRtmpL8EfDnwH5VdXeSXYCtR1yWJEnaSDkCQ5IkDcsjgVur6m6Aqrq1qpYmeXuSS5NckeTUJIH7RlD83yTfSbIkyf5Jzk3ykyTvbtvMSXJ1kk+1ozq+kOSh/TtOcmiS/0qyMMnnk0xvl783yVXttu+bwmMhSZLWkwGGJEkalguAWUl+nOQjSZ7ZLv9QVe1fVU8GtqUZpTHmnqp6BnAK8CXgdcCTgaOTPLxt8wTg1KraG1gJ/G3vTtuRHm8Dnl1V+wELgDcleRjwQuBJ7bbvHsJjliRJQ2KAIUmShqKqVgP/D3AMsAw4O8nRwMFJLkmyGHgW8KSezc5v/10MXFlVN7UjOK4DZrXrbqiq77W3PwMc1LfrpwF7Ad9Lsgg4CtidJuy4C/h4khcBv9lQj1WSJA2fc2BIkqShqap7gW8D324Di78B9gbmVtUNSU4EtunZ5O7239/13B67P/b/LdW/m777AS6sqiP760lyAHAIcATwepoARZIkbQQcgSFJkoYiyROSPK5n0VOAa9rbt7bzUvzlg+h6djtBKMCRwMV9638AHJjksW0dD03y+HZ/O1bV14Dj2nokSdJGwhEYkiRpWKYD/5pkJ2ANcC3N6SQraE4RuR649EH0uwQ4KslHgZ8A/9a7sqqWtaeqfDbJQ9rFbwNWAV9Ksg3NKI3/9SD2LUmSRiRV/aMuJUmSuinJHOAr7QSgkiRpM+IpJJIkSZIkqfMcgSFJkiRJkjrPERiSJEmSJKnzDDAkSZIkSVLnGWBIkiRJkqTOM8CQJEmSJEmdZ4AhSZIkSZI67/8HBHD5zux9TakAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# Boxplot of CPM-normalized counts\n", "plt.figure(figsize=(15, 5))\n", "sns.boxplot(data=np.log10(raw_summary.T[['mean']] + 1), palette=\"coolwarm\")\n", "plt.title(\"Distribution of Log10-Transformed mean CPM Across Samples\")\n", "plt.xlabel(\"Samples\")\n", "plt.ylabel(\"Log10(CPM + 1)\")\n", "plt.xticks(rotation=90)\n", "plt.tight_layout()\n", "plt.show()\n" ] }, { "cell_type": "code", "execution_count": 42, "id": "accredited-lesbian", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAABDAAAAFgCAYAAABNIolGAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAAAnuElEQVR4nO3debgld1kn8O8bAiZkIYGEBgIksomAQZkmoEEQgsEgKMuobGERJmZcGXFc0QFkBnTQAUWIgGzNKgwREcEw8UFAROmwJUBUhEAgkHQgnQWSQOCdP6pucnJz7+3b6T73VHc+n+e5T59TVedX76lT5yT1rV/9qro7AAAAAFO2z6ILAAAAANgRAQYAAAAweQIMAAAAYPIEGAAAAMDkCTAAAACAyRNgAAAAAJMnwAC4gaiqU6rqd3dTW7evqsuq6kbj8/dW1dN2R9tje++qqiftrvZ2Yr3PraoLq+orG73uvcGUt9/u3kdhVlX9SFV9cdF1AOztBBgAe4GqOqeqLq+qS6tqe1V9sKpOrqqrf+e7++Tu/v11tvXgtZbp7i9094Hd/e3dUPuzqup1y9o/obtfs6tt72Qdt0vyjCR36+5brTB/LgcoVfX7VXVmVV1VVc9aYf7jqurzVfX1qvqrqrr5CsucMgZKl1XVN6vqWzPP37W7a17lfay5/bh+qurWVfUXVfXl8ft9dlU9u6oOGOf3uG9cVlVfqqo/ngkWzxn3h8OWtfmx8XVH7WDdrx73y9vM7Q3uoqq6bVX93zE4u3j8Lj150XUBMB8CDIC9x8O7+6AkRyZ5fpLfSPIXu3slVbXv7m5zIo5M8tXuvmCD1/uZJL+e5J3LZ1TV3ZP8eZITk2xK8o0kL1m+3BhOHdjdByb5X0nevPS8u0+YaW+en9313n578T61S8aw6p+S7J/kB8fv948mOSTJHWcWvef42R+X5HFJ/svMvM8leexMm983trejdR+Q5NFJLk7y+OtZ/0Z8rluSnJth/7tFkicmOX8D1gvAAggwAPYy3X1xd/91kp9J8qSqukdy9dnU546PD6uqvxl7a3ytqt5fVftU1ZYkt0/yjvGM7q9X1VHj2dqnVtUXkvz9zLTZA5Q7VtW/jGdB377UU2ClngtLvTyq6seS/HaSnxnX9/Fx/tXd/ce6njn2Qrigql5bVTcb5y3V8aSq+sJ4FvZ3Vts2VXWz8fXbxvaeObb/4CTvSXKbsY5X78w2r6rvHWveXlWfrKqfmJl3i6p6R1VdUlUfruEyiw/MfF6v6e53Jbl0haYfn+Qd3f2+7r4sye8meVRVHbQTtZ1TVb9RVZ9I8vWq2reqfrOq/mM8o/+pqnrkzPJPrqoPVNULquqiqvpcVZ2wbP5nx9d+rqoev9r2q6qfGLfH9nH7fO8add1p/CyfUlXnjus+uaruXVWfGNt48bL39rNV9elx2b+rqiNn5v1oDb0VLh5fV2tso2dV1Vuq6nXj+zqzqu5SVb817nPnVtXxM8vfrK7pFfGl8TNd6vVwx6r6+6r66rg/vr6qDln2vn9tfE8XV9Wbq2q/VUr71Qz7xRO6+5wk6e5zu/tXuvsTyxfu7rOTvD/JPWYmb8lwUL/kSUleu9q2mPHoJNuTPGd8zez2unlVvaqqzhu3/V+N03+kqr44fq5fSfKqqvquqnrhuOx54+PvGpdf8XdonPcb47a9tKr+taqOW6XOeyd5dXd/vbuv6u6Pjt+npVrfUlVfGbf1+2oIBZfmvbqqXlLDJWuXVdU/VtWtxhovGvefH5hZ/pxxn/jUOP9Vq312VXWbGnqGbBu/J788M++Yqtpaw2/C+VX1x+v4PACIAANgr9Xd/5Lki0l+eIXZzxjnHZ7hzP5vDy/pE5N8IUNvjgO7+w9nXvOAJN+b5CGrrPKJSX42yW2SXJXkT9ZR47tz7R4D91xhsSePfw9McockByZ58bJl7pfkezKcgf692QPlZf40yc3Gdh4w1vyU7v5/SU5Ict5Yx5N3VPuSqrpxknckOS3JLZP8UpLXV9X3jIv8WZKvJ7lVhgPBnRnb4+5JPr70pLv/I8k3k9xlJ9pIhjPwP57kkO6+Ksl/ZNgvbpbk2UleV1W3nln+Pkn+NclhSf4wyV/U4IAMn+sJY2+AH0rysZW2X1XdJckbkzw9w372txmCsZusVFeGfWZp3XfOEMC9MMnvJHnwuC1+uqoekCRV9YgM++2jxvbfP64vNVwy8X+TPHN8D/+R5NgdbKOHZzjYPzTJR5P8XYb/Tzoiw0H8n88s+5qx3jsl+YEkxydZGl+jkjwvw/fge5PcLsmzlq3rp5P8WJLvTnJ0hv17JQ9O8rbu/s4Oah9WXHW3DJ/rR2cmfyjJwTWEbDfKsF1ft9Lrl3lShu35piR3rap7zczbkuSmGT6TWyb5PzPzbpXk5hl6RJyU4fO7b5LvT3LPJMdk+FySVX6Hxu/OLya597ifPSTJOavU+aEkf1ZVj6mq268w/10Z9qdbJvlIktcvm//TuWY/uTJDj5ePjM/fmmR5uPD4sZ47ZvgePnPZ/IwhzDsyfHePyPC79PSqWvrtfFGSF3X3wWM7f7nKewNgGQEGwN7tvAwHE8t9K8mtkxzZ3d/q7vd3d++grWeNZzkvX2X+lu4+q7u/nqGnwE8vnZXeRY9P8sfd/dmxF8JvJXlMXbv3x7O7+/Lu/niGg4brBCEzB2+/1d2Xjme0/yjD5Rm74r4ZQpXnd/c3u/vvk/xNkseO63x0kv/R3d/o7k9lOPhdrwMzdOGfdXGSdffAGP3JeOb+8iTp7rd093nd/Z3ufnOSf89wYLnk89398nGMk9dk2Fc2jfO+k+QeVbV/d3+5uz+5yjp/Jsk7u/s93f2tJC/IcOnCD61W1+j3u/uK7j4tQ/Dzxu6+oLu/lCGkWDoj/nNJntfdnx5Dmf+V5PvHXhgPTfKp7n7ruO4XJtnRwKLv7+6/G9t6S4aD6uePr39TkqOq6pCq2pQhrHn6+H24IMMB/GPGbfuZ8T1f2d3bMhwAP2DZuv5k3P5fy3Cg+/2r1HSLJF/eQd1J8pGqumhs6xVJXrVs/lIvjB9NcnaSL63V2BgEPDDJG7r7/CSnZwzexqDrhCQnd/dF4+/HP8y8/DsZ9vcrx8/18UmeM36G2zIEZkvfudV+h76d5LuS3K2qbtzd54zh3Up+KsN+8btJPlfD+B73XprZ3a8cv+9XZgiS7lljD67Rqd19RndfkeTUJFd092vHff/NuWZ/W/LicZ/9WpL/mZnLc2bcO8nh3f2c8Tfhs0lennEfGd/3narqsO6+rLs/tMp7A2AZAQbA3u2IJF9bYfr/zjD2wmk1XA7wm+to69ydmP/5JDfOcBZzV91mbG+27X1zzQF1cu2D029kOPBf7rAkN1mhrSN2Q33nLjtLvtTu4WOts9tmR9tx1mVJDl427eAkl1bVD9c1A3WuFiKsuM6qeuJ4oLe9qrZnuORg9rO6ent29zfGhweO4dTPJDk5yZer6p1VdddV1nmtz23cPufm2tt7pW0xO37B5Ss8X/psj0zyopn38LUMvR+OGNd9ddvjQfGOtvvy9VzY1wxSuxSwHDiu98YZ3v/Suv88wxn+VNUtq+pN4+UPl2To7bD8e7Ce/TVJvprhAH9H7tXdh3b3Hbv7mSv02NiSYWyMJ2d9l4+cmOTT3f2x8fnrkzxu7G10uyRf6+6LVnnttjEMWLLS93dpUNAVf4e6+zMZeu48K8kF4/ZccSDRMUT5ze6+e4bfhI8l+auxx9CNqur5NVwudUmu6cUx+3msd39bsvx3bqW6jsxwOdX2mX3kt3PNb9ZTM/TeOLuGy8oettJ7A+C6BBgAe6nxLOQRST6wfN54RvIZ3X2HDF3nf3XmGvPVemLsqIfG7WYe3z7DWcYLM5xFv+lMXTfKcGC/3nbPy3BAMNv2Vdn5gfouHGta3taaZ6PX4bwkt6uZO77MtLstQ623nZk3u5125JOZ6U1SVXfIcGb638az1UsDdd591RYGV2/jsYfCyzN00b9Fdx+S5KysMUbEtRoaein8aIYD67PHtlZyrc+tqirDe5/d3jv67NdybpKf6+5DZv727+4PZui1cPV2nln37nBuhksNDptZ78Ezn8HzMryvo8dLBJ6QdW7bFfy/JI9ctm/ttO7+fIbBPB+a5G3reMkTk9xhHDviKxl6kRyWoefFuUluXjPjeixf3bLnK31/zxvrWvV3qLvf0N33G1/bSf5gHe/zwgw9fW6ToefZ45L8ZIZLcW6W5Khx0ev7eSTX/Z07b4Vlzk3yuWX75kHd/dCxzn/v7sdmCL3+IMlbx8uzANgBAQbAXqaqDh7P6L0pyeu6+8wVlnlYDYMmVpJLMnTZXjrbfH6GMSJ21hOq6m5VddMMYwa8dTyD/W9J9quqHx/P4D4zw0H4kvMzdM9f7b9Jb0zy36rqu6tq9i4bV62y/IrGWv4yyf+sqoPGA/lfzfrGA7haVe03+5fkXzKENL9eVTeuqh/JcDD2pnGdb0vyrKq66dhb4YnL2rvx2M4+SfYd21269Ob1SR4+9rY4IMN2fVt3rzTg53odkOGAcNu4/qfk2oM+rvXeN9UwMOcBGQ7iL8s1+81yf5nkx6vquPFzf8b4mg/uQu2zTknyWzUOyljDwJo/Nc57Z5K7V9WjxkuNfjnD2Ay7rLu/nGG8kz8av2v71DBw59JlIgdl2C7bq+qIJP99F1b3xxl63Lxm3F9TVUfUcKvUo3eyracmedDYi2ZVVfWDGcZlOCbDpS3fn2H/eEOSJ43v/11JXlJVh4777/3XaPKNSZ5ZVYfXMDbJ72X8zq32O1RV31NVD6phsM8rMvSEWHE/q6o/qKp71DA47UFJ/muSz3T3VzN8Fldm6Mly0wy/HbvqF2q4devNM/SqePMKy/xLkktqGIh0/7EnyD3GUDlV9YSqOnzsKbN9fM0u35Ia4IZAgAGw93hHVV2a4ezf72Q4+HnKKsveOcPZ3csyDFr3ku5+7zjveRkOOLZX1a/txPq3JHl1hu7x+2U4aEx3X5zk5zNcm/+lDAf7s3clecv471er6iMrtPvKse33ZTiLfEWGgTKvj18a1//ZDD1T3jC2v15HZDiYmv27XZKfyHB2+sIMtzl9Yg93hEiGng43y7BdtmQ4oLtyps2Xj+08NsPndnnGMQLG8SVOzhBkXJDhgOznd6Le6+hhHI4/yvC5n5/k+5L84zpfvk+GIOK8DJdsPGC1err7XzP0PvjTDNvl4RkGh/3mrtQ/0/6pGc5ev2m8POCsDJ/B0pn4n8pwO+GvZtjf1/se1+OJGS5H+lSSizIM9rh0qcezk9wrw1gl78z6ejysaBxn4Ycy9Bz65/H7ffrY9md2sq3/6O6t61j0SUne3t1ndvdXlv4yDDz5sPHA/cSxprMz7JdPX6O95ybZmuQTSc7MMEDmc8d5q/0OfVeGz+7CDN+bW2YIC1Zy0wxjV2zP8L0+MsP3MRkul/l8ht+dT2UY8HNXvSFDgPXZ8e+5yxcYg8uHZwh/Pje+j1dk+B1IhgFcP1lVl2XYro9ZdtkNAKuo3uGYbQDA7lJVf5DkVt29M3cjARasqs5J8rQe7roDwALogQEAc1RVd62qo8dBBY/J0JX/1EXXBQCwp9l3x4sAALvgoAyXjdwmQ3f7P0ry9oVWBACwB3IJCQAAADB5LiEBAAAAJm+Pu4TksMMO66OOOmrRZQAAAABzcMYZZ1zY3Ycvnz7XAGMcrfnSDPe2vqq7Ny+b/5NJfj/Jd5JcleTp3f2Btdo86qijsnXreu4CBgAAAOxpqurzK03fiB4YDxzvx76S05P8dXd3VR2d5C+T3HUDagIAAAD2IAu9hKS7L5t5ekASI4oCAAAA1zHvQTw7yWlVdUZVnbTSAlX1yKo6O8k7k/zsKsucVFVbq2rrtm3b5lguAAAAMEXzDjCO7e57JTkhyS9U1f2XL9Ddp3b3XZM8IsN4GNfR3S/r7s3dvfnww68zjgcAAACwl5trgNHd543/XpDk1CTHrLHs+5LcsaoOm2dNAAAAwJ5nbgFGVR1QVQctPU5yfJKzli1zp6qq8fG9ktwkyVfnVRMAAACwZ5rnIJ6bkpw65hP7JnlDd7+7qk5Oku4+Jcmjkzyxqr6V5PIkP9PdBvIEAAAArqX2tLxg8+bNvXXr1kWXAQAAAMxBVZ3R3ZuXT5/3IJ4AAAAAu0yAAQAAAEzePMfAAGAFp512Ws4///xFlwF7hIsuuihJcuihhy64EtgzbNq0Kccff/yiywCYCwEGADBZ3/zmNxddAgAwEQIMgA3mzBis35YtW5IkJ5544oIrAQAWzRgYAAAAwOQJMAAAAIDJE2AAAAAAkyfAAAAAACZPgAEAAABMngADAAAAmDwBBgAAADB5AgwAAABg8gQYAAAAwOQJMAAAAIDJE2AAAAAAkyfAAAAAACZPgAEAAABMngADAAAAmDwBBgAAADB5AgwAAABg8gQYAAAAwOQJMAAAAIDJE2AAAAAAkyfAAAAAACZPgAEAAABMngADAAAAmDwBBgAAADB5AgwAAABg8gQYAAAAwOQJMAAAAIDJE2AAAAAAkyfAAAAAACZPgAEAAABMngADAAAAmDwBBgAAADB5AgwAAABg8gQYAAAAwOQJMAAAAIDJE2AAAAAAkyfAAAAAACZPgAEAAABMngADAAAAmDwBBgAAADB5AgwAAABg8gQYAAAAwOQJMAAAAIDJ23eejVfVOUkuTfLtJFd19+Zl8x+f5DfGp5cl+a/d/fF51gQAAADseeYaYIwe2N0XrjLvc0ke0N0XVdUJSV6W5D4bUBMAAACwB9mIAGNV3f3BmacfSnLbRdUCAAAATNe8x8DoJKdV1RlVddIOln1qknetNKOqTqqqrVW1ddu2bbu9SAAAAGDa5t0D49juPq+qbpnkPVV1dne/b/lCVfXADAHG/VZqpLtfluHykmzevLnnWTAAAAAwPXPtgdHd543/XpDk1CTHLF+mqo5O8ookP9ndX51nPQAAAMCeaW4BRlUdUFUHLT1OcnySs5Ytc/skb0tyYnf/27xqAQAAAPZs87yEZFOSU6tqaT1v6O53V9XJSdLdpyT5vSS3SPKScbnr3GoVAAAAYG4BRnd/Nsk9V5h+yszjpyV52rxqAAAAAPYO874LCQAAAMAuE2AAAAAAkyfAAAAAACZPgAEAAABMngADAAAAmDwBBgAAADB5AgwAAABg8gQYAAAAwOQJMAAAAIDJE2AAAAAAkyfAAAAAACZPgAEAAABMngADAAAAmDwBBgAAADB5AgwAAABg8gQYAAAAwOQJMAAAAIDJE2AAAAAAkyfAAAAAACZPgAEAAABMngADAAAAmDwBBgAAADB5AgwAAABg8gQYAAAAwOQJMAAAAIDJE2AAAAAAkyfAAAAAACZPgAEAAABMngADAAAAmDwBBgAAADB5AgwAAABg8gQYAAAAwOQJMAAAAIDJE2AAAAAAkyfAAAAAACZPgAEAAABM3r5rzayqH0zyhCQ/nOTWSS5PclaSdyZ5XXdfPPcKAQAAgBu8VXtgVNW7kjwtyd8l+bEMAcbdkjwzyX5J3l5VP7ERRQIAAAA3bGv1wDixuy9cNu2yJB8Z//6oqg6bW2UAAAAAo1V7YKwQXlyvZQAAAAB21fUaxLOqztzdhQAAAACsZtVLSKrqUavNSnKr+ZQDAAAAcF1rjYHx5iSvT9IrzNtvPuUAAAAAXNdaAcYnkrygu89aPqOqHjy/kgAAAACuba0xMJ6e5JJV5j1y95cCAAAAsLJVe2B09/vXmLd1PuUAAAAAXNf1ugsJAAAAwEaaa4BRVedU1ZlV9bGquk6vjaq6a1X9U1VdWVW/Ns9aAAAAgD3XWoN47i4P7O4LV5n3tSS/nOQRG1AHAAAAsIdadw+MqvrR3b3y7r6guz+c5Fu7u20AAABg77Ezl5D8wfVov5OcVlVnVNVJ1+P1SZKqOqmqtlbV1m3btl3fZgAAAIA91LwH8Ty2u++V5IQkv1BV978+jXT3y7p7c3dvPvzww3dvhQAAAMDkrTkGRlW9KkMvikpy+6p65dK87v7ZHTXe3eeN/15QVacmOSbJ+3apYgAAAOAGZ0eDeL565vH9krxmvQ1X1QFJ9unuS8fHxyd5zk5XCAAAANzgrRlgdPc/LD2uqktnn6/DpiSnVtXSet7Q3e+uqpPHtk+pqlsl2Zrk4CTfqaqnJ7lbd1+yc28DAAAA2JvtzG1Uv7kzDXf3Z5Pcc4Xpp8w8/kqS2+5MuwAAAMANz7oH8ezu+86zEAAAAIDVzPsuJAAAAAC7TIABAAAATJ4AAwAAAJi8nRnEE1Z02mmn5fzzz190GQDshZb++7Jly5YFVwLA3mbTpk05/vjjF10GO2HVAKOqLk3Ss5PG55Wku/vgOdfGHuL888/Pl7/85Rx6yCGLLgWAvcw++wydRa+4/PIFVwLA3uSi7dsXXQLXw1o9ME5Pcqskb0vypu7+wsaUxJ7o0EMOyXHHHbfoMgAAAHbo9NNPX3QJXA+rjoHR3Y9I8pAk25K8vKr+oap+vqpuvlHFAQAAACQ7GMSzuy/u7lclOSHJKUmek+TJG1AXAAAAwNXWHMSzqn4oyWOT/HCSDyR5ZHe/fyMKAwAAAFiy1iCe5yTZnuRNSU5KctU4/V5J0t0fmX95AAAAAGv3wDgnw11HHjL+zeokD5pTTQAAAADXsmqA0d0/soF1AAAAAKxq1UE8q+rOVfVXVXVWVb2xqo7YyMIAAAAAlqx1F5JXJnlnkkcn+UiSP92QigAAAACWWWsMjIO6++Xj4/9dVQbtBAAAABZirQBjv6r6gSQ1Pt9/9rm7kAAAAAAbZa0A4ytJ/niV5+5CAgAAAGwYdyEBAAAAJm+tu5A8oapOXGH6f6mqx823LAAAAIBrrHUXkmck+asVpr95nAcAAACwIdYKMG7U3Zcun9jdlyS58fxKAgAAALi2tQKMG1fVAcsnVtVBSW4yv5IAAAAArm2tAOMvkry1qo5amjA+ftM4DwAAAGBDrHUXkhdU1WVJ/qGqDsxw69SvJ3l+d790owoEAAAAWDXASJLuPiXJKWOAUSuNiQEAAAAwbzu6jeo+SdLdly0PL6rqjlV1v3kXCAAAALBWD4xbJPloVZ2R5Iwk25Lsl+ROSR6Q5MIkvzn3CgEAAIAbvLXGwHhRVb04yYOSHJvk6CSXJ/l0khO7+wsbUyIAAABwQ7ejMTC+neQ94x8AAADAQqwZYFTVQ5I8IskRGe5Ccl6St3f3u+dfGgAAAMBg1QCjql6Y5C5JXpvki+Pk2yb55ao6obt/Zf7lAQAAAKzdA+Oh3X2X5ROr6s1J/i2JAAMAAADYEKveRjXJFVV1zArT753kijnVAwAAAHAda/XAeHKSl1bVQbnmEpLbJblknAcAAACwIda6jepHktynqm6VYRDPSvLF7v7KRhUHAAAAkOzgLiRJMgYW1wotququ3X323KoCAAAAmLHWGBhrOW23VgEAAACwhrVuo/onq81KcshcqgEAAABYwVqXkDwlyTOSXLnCvMfOpxwAAACA61orwPhwkrO6+4PLZ1TVs+ZWEQAAAMAyawUY/znJFSvN6O7vnk85AAAAANe11m1Uv7aRhQAAAACsZoe3Ua2qM5P0sskXJ9ma5Lnd/dV5FAYAAACwZIcBRpJ3Jfl2kjeMzx8z/ntJklcnefjuLwsAAADgGusJMI7t7mNnnp9ZVf/Y3cdW1RPmVRgAAADAkn3WscyBVXWfpSdVdUySA8enV82lKgAAAIAZ6+mB8bQkr6yqA5NUhktHnlpVByR53lovrKpzklya4RKUq7p787L5leRFSR6a5BtJntzdH9nZN8FiXXTRRbniiity+umnL7oUAACAHbpo+/bsd8WKN91kwnYYYHT3h5N8X1XdLEl19/aZ2X+5jnU8sLsvXGXeCUnuPP7dJ8lLx38BAAAArraeu5DcLMn/SHL/8fk/JHlOd1+8G9b/k0le292d5ENVdUhV3bq7v7wb2maDHHroobni8stz3HHHLboUAACAHTr99NOz3/77L7oMdtJ6xsB4ZYbLQH56/LskyavW2X4nOa2qzqiqk1aYf0SSc2eef3Gcdi1VdVJVba2qrdu2bVvnqgEAAIC9xXrGwLhjdz965vmzq+pj62z/2O4+r6pumeQ9VXV2d79vZn6t8Jq+zoTulyV5WZJs3rz5OvMBAACAvdt6emBcXlX3W3pSVccmuXw9jXf3eeO/FyQ5Nckxyxb5YpLbzTy/bZLz1tM2AAAAcMOxngDj5CR/VlXnjHcVeXGSn9vRi6rqgKo6aOlxkuOTnLVssb9O8sQa3DfJxca/AAAAAJZbz11IPp7knlV18Pj8kqp6epJP7OClm5KcOtwpNfsmeUN3v7uqTh7bOSXJ32a4hepnMtxG9SnX830AAAAAe7H1jIGRZAguZp7+apIX7mD5zya55wrTT5l53El+Yb01AAAAADdM67mEZCUrDb4JAAAAMBfXN8BwJxAAAABgw6x6CUlVXZqVg4pKsv/cKgIAAABYZtUAo7sP2shCAAAAAFZzfS8hAQAAANgwAgwAAABg8gQYAAAAwOQJMAAAAIDJE2AAAAAAkyfAAAAAACZPgAEAAABMngADAAAAmDwBBgAAADB5AgwAAABg8gQYAAAAwOQJMAAAAIDJE2AAAAAAkyfAAAAAACZPgAEAAABMngADAAAAmDwBBgAAADB5AgwAAABg8gQYAAAAwOQJMAAAAIDJE2AAAAAAkyfAAAAAACZPgAEAAABMngADAAAAmDwBBgAAADB5AgwAAABg8gQYAAAAwOQJMAAAAIDJ23fRBbB3uGj79px++umLLgOAvcyll12WJDnowAMXXAkAe5OLtm/Prffff9FlsJMEGOyyTZs2LboEAPZSF19ySZJkP/+TCcBudOv993ccswcSYLDLjj/++EWXAMBeasuWLUmSE088ccGVAACLZgwMAAAAYPIEGAAAAMDkCTAAAACAyRNgAAAAAJMnwAAAAAAmT4ABAAAATJ4AAwAAAJg8AQYAAAAweQIMAAAAYPIEGAAAAMDkCTAAAACAyRNgAAAAAJM39wCjqm5UVR+tqr9ZYd6hVXVqVX2iqv6lqu4x73oAAACAPc9G9MD4lSSfXmXebyf5WHcfneSJSV60AfUAAAAAe5i5BhhVddskP57kFasscrckpydJd5+d5Kiq2jTPmgAAAIA9z7x7YLwwya8n+c4q8z+e5FFJUlXHJDkyyW2XL1RVJ1XV1qraum3btjmVCgAAAEzV3AKMqnpYkgu6+4w1Fnt+kkOr6mNJfinJR5NctXyh7n5Zd2/u7s2HH374XOoFAAAApmvfObZ9bJKfqKqHJtkvycFV9brufsLSAt19SZKnJElVVZLPjX8AAAAAV5tbD4zu/q3uvm13H5XkMUn+fja8SJKqOqSqbjI+fVqS942hBgAAAMDV5tkDY0VVdXKSdPcpSb43yWur6ttJPpXkqRtdDwAAADB9GxJgdPd7k7x3fHzKzPR/SnLnjagBAAAA2HPN+y4kAAAAALtMgAEAAABMngADAAAAmDwBBgAAADB5AgwAAABg8gQYAAAAwOQJMAAAAIDJE2AAAAAAkyfAAAAAACZPgAEAAABMngADAAAAmDwBBgAAADB5AgwAAABg8gQYAAAAwOQJMAAAAIDJE2AAAAAAkyfAAAAAACZPgAEAAABMngADAAAAmDwBBgAAADB5AgwAAABg8gQYAAAAwOQJMAAAAIDJE2AAAAAAkyfAAAAAACZPgAEAAABMngADAAAAmDwBBgAAADB5AgwAAABg8gQYAAAAwOQJMAAAAIDJE2AAAAAAk7fvogsAuKE57bTTcv755y+6DNgjLH1XtmzZsuBKYM+wadOmHH/88YsuA2AuBBgAwGTd5CY3WXQJAMBECDAANpgzYwAAsPOMgQEAAABMngADAAAAmDwBBgAAADB5AgwAAABg8gQYAAAAwOQJMAAAAIDJE2AAAAAAkyfAAAAAACZPgAEAAABMngADAAAAmDwBBgAAADB51d2LrmGnVNW2JJ9fdB0AwIY5LMmFiy4CANgwR3b34csn7nEBBgBww1JVW7t786LrAAAWyyUkAAAAwOQJMAAAAIDJE2AAAFP3skUXAAAsnjEwAAAAgMnTAwMAAACYPAEGAAAAMHkCDAAAAGDyBBgAAADA5AkwAIDJqar9Vph22CJqAQCmQYABAEzRh6vqvktPqurRST64wHoAgAXbd9EFAACs4HFJXllV701ymyS3SPKghVYEACxUdfeiawAAuI6qekSSLUkuTXL/7v7MYisCABZJDwwAYHKq6i+S3DHJ0UnukuQdVfXi7v6zxVYGACyKMTAAgCk6K8kDu/tz3f13Se6b5F4LrgkAWCCXkAAAAACT5xISAGByqurOSZ6X5G5Jrr6lanffYWFFAQAL5RISAGCKXpXkpUmuSvLAJK/NMKAnAHADJcAAAKZo/+4+PcPlrp/v7mfFbVQB4AbNJSQAwBRdUVX7JPn3qvrFJF9KcssF1wQALJBBPAGAyamqeyf5dJJDkvx+koOT/GF3//Mi6wIAFkeAAQBMTlVtTvI7SY5McuNxcnf30YurCgBYJAEGADA5VfWvSf57kjOTfGdpend/fmFFAQALZQwMAGCKtnX3Xy+6CABgOvTAAAAmp6qOS/LYJKcnuXJpene/bWFFAQALpQcGADBFT0ly1wzjXyxdQtJJBBgAcAMlwAAApuie3f19iy4CAJiOfRZdAADACj5UVXdbdBEAwHQYAwMAmJyq+nSSOyb5XIYxMCpuowoAN2gCDABgcqrqyJWmu40qANxwCTAAAACAyTMGBgAAADB5AgwAAABg8gQYAMBuVVW/U1WfrKpPVNXHquo+c1zXe6tq87zaBwCmY99FFwAA7D2q6geTPCzJvbr7yqo6LMlNFlwWALAX0AMDANidbp3kwu6+Mkm6+8LuPq+qfq+qPlxVZ1XVy6qqkqt7UPyfqnpfVX26qu5dVW+rqn+vqueOyxxVVWdX1WvGXh1vraqbLl9xVR1fVf9UVR+pqrdU1YHj9OdX1afG175gA7cFALAbCTAAgN3ptCS3q6p/q6qXVNUDxukv7u57d/c9kuyfoZfGkm929/2TnJLk7Ul+Ick9kjy5qm4xLvM9SV7W3UcnuSTJz8+udOzp8cwkD+7ueyXZmuRXq+rmSR6Z5O7ja587h/cMAGwAAQYAsNt092VJ/lOSk5JsS/LmqnpykgdW1T9X1ZlJHpTk7jMv++vx3zOTfLK7vzz24PhsktuN887t7n8cH78uyf2Wrfq+Se6W5B+r6mNJnpTkyAxhxxVJXlFVj0ryjd31XgGAjWUMDABgt+rubyd5b5L3joHFzyU5Osnm7j63qp6VZL+Zl1w5/vudmcdLz5f+X6WXr2bZ80rynu5+7PJ6quqYJMcleUySX8wQoAAAexg9MACA3aaqvqeq7jwz6fuT/Ov4+MJxXIr/fD2avv04QGiSPDbJB5bN/1CSY6vqTmMdN62qu4zru1l3/22Sp4/1AAB7ID0wAIDd6cAkf1pVhyS5KslnMlxOsj3DJSLnJPnw9Wj300meVFV/nuTfk7x0dmZ3bxsvVXljVX3XOPmZSS5N8vaq2i9DL43/dj3WDQBMQHUv74EJADAdVXVUkr8ZBwAFAG6gXEICAAAATJ4eGAAAAMDk6YEBAAAATJ4AAwAAAJg8AQYAAAAweQIMAAAAYPIEGAAAAMDk/X+kUIWDVfhtqgAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# Boxplot of CPM-normalized max counts\n", "plt.figure(figsize=(15, 5))\n", "sns.boxplot(data=np.log10(raw_summary.T[['max']] + 1), palette=\"coolwarm\")\n", "plt.title(\"Distribution of Log10-Transformed mean CPM Across Samples\")\n", "plt.xlabel(\"Samples\")\n", "plt.ylabel(\"Log10(CPM + 1)\")\n", "plt.xticks(rotation=90)\n", "plt.tight_layout()\n", "plt.show()\n" ] }, { "cell_type": "code", "execution_count": null, "id": "necessary-wonder", "metadata": {}, "outputs": [], "source": [ "# Boxplot of CPM-normalized max counts\n", "plt.figure(figsize=(15, 5))\n", "sns.boxplot(data=np.log10(raw_summary.T[['max']] + 1), palette=\"coolwarm\")\n", "plt.title(\"Distribution of Log10-Transformed mean CPM Across Samples\")\n", "plt.xlabel(\"Samples\")\n", "plt.ylabel(\"Log10(CPM + 1)\")\n", "plt.xticks(rotation=90)\n", "plt.tight_layout()\n", "plt.show()" ] }, { "cell_type": "code", "execution_count": 38, "id": "consecutive-death", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
countmeanstdmin25%50%75%max
SRX984553432707.0165.4400531075.3328470.02.000023.895103.660091474.000
SRX984554132707.0171.200147957.3850640.02.000025.000110.792085682.000
SRX984552632707.0192.5777871048.7473220.04.000033.037131.232587475.000
SRX984553032707.0211.8516581345.6283960.02.489028.248126.0000137044.000
SRX984554332707.0136.141405686.8966270.02.000020.40696.000069707.000
SRX984552732707.0179.9949301453.5557550.02.000022.00098.0000147753.000
SRX984550632707.0183.0990911201.9312380.02.000025.318112.5605112650.000
SRX984551832707.0180.849091974.8529110.05.581532.227125.000086589.000
SRX984550132707.0159.277506911.7367450.02.403525.640111.000086268.000
SRX984553632707.0236.0523461271.6522150.05.000043.000172.0000100197.000
SRX984552932707.0211.1035461354.0492510.02.000027.000129.0000120295.000
SRX984552432707.0251.5698221655.7177460.05.000041.000164.9465164959.000
SRX984550732707.0166.6715541003.0541310.01.016519.48996.7280102536.000
SRX984554032707.0154.656215845.1411900.01.000018.12694.179584558.000
SRX984551732707.0153.901693878.5965500.02.537524.088103.999582089.000
SRX984553232707.0238.1128361340.4376560.03.000029.000145.0000124108.000
SRX984550332707.0155.634818847.0855870.02.000026.000110.383583942.000
SRX984550832707.0164.501283940.2548340.03.000025.654108.0000101989.000
SRX984553732707.0242.0174521390.8399930.05.000042.000170.6250112791.000
SRX984552032707.0210.4081551202.5193900.01.000016.035104.0000108452.000
SRX984552332707.0138.826165729.0335000.02.000023.000101.000071338.000
SRX984552832707.0200.7261271217.1333250.03.000028.208126.0380116776.000
SRX984551632707.0209.7039701261.8417960.03.000029.000129.0225120823.000
SRX984551432707.0246.9316431886.7954100.04.000036.667148.7890192241.000
SRX984552232707.0276.1270391382.4036120.06.000050.302201.2240112666.000
SRX984550232707.0138.307777821.7174050.01.272518.18185.089080409.000
SRX984550932707.0132.864739728.7802360.02.000022.00092.000074524.000
SRX984553932707.0168.8281541117.8061140.02.000023.000103.1225103742.000
SRX984554232707.0143.742371902.3377290.01.307020.01592.832097866.000
SRX984550432707.0144.648181773.4704710.02.000023.27599.000073309.000
SRX984553832707.0236.8938081763.4671410.04.000034.332143.0000168449.000
SRX984551532707.0247.2775371320.9288920.04.000039.000168.0000121237.000
SRX984551932707.0178.8650501036.7228150.03.999529.826125.524577518.988
SRX984553332707.0247.8087221637.1091460.02.000022.000124.0000152243.000
SRX984550532707.0171.201700986.9122110.03.000027.522113.682097256.000
SRX984553132707.0206.4631991131.4805460.02.000023.721125.000098612.000
SRX984551032707.0165.7417861038.3593130.03.000028.000119.000093615.000
SRX984553532707.0199.5034171311.9556810.03.000028.601118.9225130487.000
SRX984551332707.0191.5383721232.9820140.04.000030.000121.0000113471.000
SRX984551232707.0131.001250721.9620740.02.000020.70093.943577585.000
SRX984552532707.0230.7787461613.9789830.03.000030.670139.0000145177.000
SRX984552132707.0210.0076631196.1413460.03.000031.306134.0000116651.000
SRX984550032707.0148.738319959.9657340.01.000018.32088.000095149.000
SRX984551132707.0129.222067688.4283950.02.000020.00090.000067648.000
\n", "
" ], "text/plain": [ " count mean std min 25% 50% 75% \\\n", "SRX9845534 32707.0 165.440053 1075.332847 0.0 2.0000 23.895 103.6600 \n", "SRX9845541 32707.0 171.200147 957.385064 0.0 2.0000 25.000 110.7920 \n", "SRX9845526 32707.0 192.577787 1048.747322 0.0 4.0000 33.037 131.2325 \n", "SRX9845530 32707.0 211.851658 1345.628396 0.0 2.4890 28.248 126.0000 \n", "SRX9845543 32707.0 136.141405 686.896627 0.0 2.0000 20.406 96.0000 \n", "SRX9845527 32707.0 179.994930 1453.555755 0.0 2.0000 22.000 98.0000 \n", "SRX9845506 32707.0 183.099091 1201.931238 0.0 2.0000 25.318 112.5605 \n", "SRX9845518 32707.0 180.849091 974.852911 0.0 5.5815 32.227 125.0000 \n", "SRX9845501 32707.0 159.277506 911.736745 0.0 2.4035 25.640 111.0000 \n", "SRX9845536 32707.0 236.052346 1271.652215 0.0 5.0000 43.000 172.0000 \n", "SRX9845529 32707.0 211.103546 1354.049251 0.0 2.0000 27.000 129.0000 \n", "SRX9845524 32707.0 251.569822 1655.717746 0.0 5.0000 41.000 164.9465 \n", "SRX9845507 32707.0 166.671554 1003.054131 0.0 1.0165 19.489 96.7280 \n", "SRX9845540 32707.0 154.656215 845.141190 0.0 1.0000 18.126 94.1795 \n", "SRX9845517 32707.0 153.901693 878.596550 0.0 2.5375 24.088 103.9995 \n", "SRX9845532 32707.0 238.112836 1340.437656 0.0 3.0000 29.000 145.0000 \n", "SRX9845503 32707.0 155.634818 847.085587 0.0 2.0000 26.000 110.3835 \n", "SRX9845508 32707.0 164.501283 940.254834 0.0 3.0000 25.654 108.0000 \n", "SRX9845537 32707.0 242.017452 1390.839993 0.0 5.0000 42.000 170.6250 \n", "SRX9845520 32707.0 210.408155 1202.519390 0.0 1.0000 16.035 104.0000 \n", "SRX9845523 32707.0 138.826165 729.033500 0.0 2.0000 23.000 101.0000 \n", "SRX9845528 32707.0 200.726127 1217.133325 0.0 3.0000 28.208 126.0380 \n", "SRX9845516 32707.0 209.703970 1261.841796 0.0 3.0000 29.000 129.0225 \n", "SRX9845514 32707.0 246.931643 1886.795410 0.0 4.0000 36.667 148.7890 \n", "SRX9845522 32707.0 276.127039 1382.403612 0.0 6.0000 50.302 201.2240 \n", "SRX9845502 32707.0 138.307777 821.717405 0.0 1.2725 18.181 85.0890 \n", "SRX9845509 32707.0 132.864739 728.780236 0.0 2.0000 22.000 92.0000 \n", "SRX9845539 32707.0 168.828154 1117.806114 0.0 2.0000 23.000 103.1225 \n", "SRX9845542 32707.0 143.742371 902.337729 0.0 1.3070 20.015 92.8320 \n", "SRX9845504 32707.0 144.648181 773.470471 0.0 2.0000 23.275 99.0000 \n", "SRX9845538 32707.0 236.893808 1763.467141 0.0 4.0000 34.332 143.0000 \n", "SRX9845515 32707.0 247.277537 1320.928892 0.0 4.0000 39.000 168.0000 \n", "SRX9845519 32707.0 178.865050 1036.722815 0.0 3.9995 29.826 125.5245 \n", "SRX9845533 32707.0 247.808722 1637.109146 0.0 2.0000 22.000 124.0000 \n", "SRX9845505 32707.0 171.201700 986.912211 0.0 3.0000 27.522 113.6820 \n", "SRX9845531 32707.0 206.463199 1131.480546 0.0 2.0000 23.721 125.0000 \n", "SRX9845510 32707.0 165.741786 1038.359313 0.0 3.0000 28.000 119.0000 \n", "SRX9845535 32707.0 199.503417 1311.955681 0.0 3.0000 28.601 118.9225 \n", "SRX9845513 32707.0 191.538372 1232.982014 0.0 4.0000 30.000 121.0000 \n", "SRX9845512 32707.0 131.001250 721.962074 0.0 2.0000 20.700 93.9435 \n", "SRX9845525 32707.0 230.778746 1613.978983 0.0 3.0000 30.670 139.0000 \n", "SRX9845521 32707.0 210.007663 1196.141346 0.0 3.0000 31.306 134.0000 \n", "SRX9845500 32707.0 148.738319 959.965734 0.0 1.0000 18.320 88.0000 \n", "SRX9845511 32707.0 129.222067 688.428395 0.0 2.0000 20.000 90.0000 \n", "\n", " max \n", "SRX9845534 91474.000 \n", "SRX9845541 85682.000 \n", "SRX9845526 87475.000 \n", "SRX9845530 137044.000 \n", "SRX9845543 69707.000 \n", "SRX9845527 147753.000 \n", "SRX9845506 112650.000 \n", "SRX9845518 86589.000 \n", "SRX9845501 86268.000 \n", "SRX9845536 100197.000 \n", "SRX9845529 120295.000 \n", "SRX9845524 164959.000 \n", "SRX9845507 102536.000 \n", "SRX9845540 84558.000 \n", "SRX9845517 82089.000 \n", "SRX9845532 124108.000 \n", "SRX9845503 83942.000 \n", "SRX9845508 101989.000 \n", "SRX9845537 112791.000 \n", "SRX9845520 108452.000 \n", "SRX9845523 71338.000 \n", "SRX9845528 116776.000 \n", "SRX9845516 120823.000 \n", "SRX9845514 192241.000 \n", "SRX9845522 112666.000 \n", "SRX9845502 80409.000 \n", "SRX9845509 74524.000 \n", "SRX9845539 103742.000 \n", "SRX9845542 97866.000 \n", "SRX9845504 73309.000 \n", "SRX9845538 168449.000 \n", "SRX9845515 121237.000 \n", "SRX9845519 77518.988 \n", "SRX9845533 152243.000 \n", "SRX9845505 97256.000 \n", "SRX9845531 98612.000 \n", "SRX9845510 93615.000 \n", "SRX9845535 130487.000 \n", "SRX9845513 113471.000 \n", "SRX9845512 77585.000 \n", "SRX9845525 145177.000 \n", "SRX9845521 116651.000 \n", "SRX9845500 95149.000 \n", "SRX9845511 67648.000 " ] }, "execution_count": 38, "metadata": {}, "output_type": "execute_result" } ], "source": [ "raw_summary.T" ] }, { "cell_type": "code", "execution_count": 40, "id": "blind-restaurant", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
mean
SRX9845534165.440053
SRX9845541171.200147
SRX9845526192.577787
SRX9845530211.851658
SRX9845543136.141405
SRX9845527179.994930
SRX9845506183.099091
SRX9845518180.849091
SRX9845501159.277506
SRX9845536236.052346
SRX9845529211.103546
SRX9845524251.569822
SRX9845507166.671554
SRX9845540154.656215
SRX9845517153.901693
SRX9845532238.112836
SRX9845503155.634818
SRX9845508164.501283
SRX9845537242.017452
SRX9845520210.408155
SRX9845523138.826165
SRX9845528200.726127
SRX9845516209.703970
SRX9845514246.931643
SRX9845522276.127039
SRX9845502138.307777
SRX9845509132.864739
SRX9845539168.828154
SRX9845542143.742371
SRX9845504144.648181
SRX9845538236.893808
SRX9845515247.277537
SRX9845519178.865050
SRX9845533247.808722
SRX9845505171.201700
SRX9845531206.463199
SRX9845510165.741786
SRX9845535199.503417
SRX9845513191.538372
SRX9845512131.001250
SRX9845525230.778746
SRX9845521210.007663
SRX9845500148.738319
SRX9845511129.222067
\n", "
" ], "text/plain": [ " mean\n", "SRX9845534 165.440053\n", "SRX9845541 171.200147\n", "SRX9845526 192.577787\n", "SRX9845530 211.851658\n", "SRX9845543 136.141405\n", "SRX9845527 179.994930\n", "SRX9845506 183.099091\n", "SRX9845518 180.849091\n", "SRX9845501 159.277506\n", "SRX9845536 236.052346\n", "SRX9845529 211.103546\n", "SRX9845524 251.569822\n", "SRX9845507 166.671554\n", "SRX9845540 154.656215\n", "SRX9845517 153.901693\n", "SRX9845532 238.112836\n", "SRX9845503 155.634818\n", "SRX9845508 164.501283\n", "SRX9845537 242.017452\n", "SRX9845520 210.408155\n", "SRX9845523 138.826165\n", "SRX9845528 200.726127\n", "SRX9845516 209.703970\n", "SRX9845514 246.931643\n", "SRX9845522 276.127039\n", "SRX9845502 138.307777\n", "SRX9845509 132.864739\n", "SRX9845539 168.828154\n", "SRX9845542 143.742371\n", "SRX9845504 144.648181\n", "SRX9845538 236.893808\n", "SRX9845515 247.277537\n", "SRX9845519 178.865050\n", "SRX9845533 247.808722\n", "SRX9845505 171.201700\n", "SRX9845531 206.463199\n", "SRX9845510 165.741786\n", "SRX9845535 199.503417\n", "SRX9845513 191.538372\n", "SRX9845512 131.001250\n", "SRX9845525 230.778746\n", "SRX9845521 210.007663\n", "SRX9845500 148.738319\n", "SRX9845511 129.222067" ] }, "execution_count": 40, "metadata": {}, "output_type": "execute_result" } ], "source": [ "raw_summary.T[['mean']]" ] }, { "cell_type": "code", "execution_count": 62, "id": "ranging-teddy", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
SRX9845534SRX9845541SRX9845526SRX9845530SRX9845543SRX9845527SRX9845506SRX9845518SRX9845501SRX9845536...SRX9845505SRX9845531SRX9845510SRX9845535SRX9845513SRX9845512SRX9845525SRX9845521SRX9845500SRX9845511
00.1365980.0769200.0688470.0000000.2236870.1788590.1763660.2660680.0762610.280275...0.4449620.2406580.0735240.4290710.4143480.1663650.5295220.1573850.0811880.338823
11.7764730.4161911.0614600.3615170.5932340.8071780.9215440.2243370.2922090.000000...0.1863241.3369820.4813250.7480011.0763060.4574470.8114910.1573851.2467520.634725
40.9761841.6150990.8809310.7362561.8613101.4611231.6704011.6308800.5004121.236466...1.5119331.5483890.8301421.2881071.5332381.1931311.2579321.2452440.9048440.868448
50.1229060.4378910.4141290.4521160.4466350.2731110.4103520.4812610.5328310.668347...0.3013210.3328360.5817910.4959570.4617560.4863700.4936210.4458080.4918500.397706
60.9889171.0671171.2533531.2692431.2145001.1733120.9454271.1955961.0605921.361088...1.0289701.3000611.0398371.2921611.4261751.0843101.1302401.0329471.1602860.942124
\n", "

5 rows × 44 columns

\n", "
" ], "text/plain": [ " SRX9845534 SRX9845541 SRX9845526 SRX9845530 SRX9845543 SRX9845527 \\\n", "0 0.136598 0.076920 0.068847 0.000000 0.223687 0.178859 \n", "1 1.776473 0.416191 1.061460 0.361517 0.593234 0.807178 \n", "4 0.976184 1.615099 0.880931 0.736256 1.861310 1.461123 \n", "5 0.122906 0.437891 0.414129 0.452116 0.446635 0.273111 \n", "6 0.988917 1.067117 1.253353 1.269243 1.214500 1.173312 \n", "\n", " SRX9845506 SRX9845518 SRX9845501 SRX9845536 ... SRX9845505 \\\n", "0 0.176366 0.266068 0.076261 0.280275 ... 0.444962 \n", "1 0.921544 0.224337 0.292209 0.000000 ... 0.186324 \n", "4 1.670401 1.630880 0.500412 1.236466 ... 1.511933 \n", "5 0.410352 0.481261 0.532831 0.668347 ... 0.301321 \n", "6 0.945427 1.195596 1.060592 1.361088 ... 1.028970 \n", "\n", " SRX9845531 SRX9845510 SRX9845535 SRX9845513 SRX9845512 SRX9845525 \\\n", "0 0.240658 0.073524 0.429071 0.414348 0.166365 0.529522 \n", "1 1.336982 0.481325 0.748001 1.076306 0.457447 0.811491 \n", "4 1.548389 0.830142 1.288107 1.533238 1.193131 1.257932 \n", "5 0.332836 0.581791 0.495957 0.461756 0.486370 0.493621 \n", "6 1.300061 1.039837 1.292161 1.426175 1.084310 1.130240 \n", "\n", " SRX9845521 SRX9845500 SRX9845511 \n", "0 0.157385 0.081188 0.338823 \n", "1 0.157385 1.246752 0.634725 \n", "4 1.245244 0.904844 0.868448 \n", "5 0.445808 0.491850 0.397706 \n", "6 1.032947 1.160286 0.942124 \n", "\n", "[5 rows x 44 columns]" ] }, "execution_count": 62, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Log-transform CPM data for Combat\n", "log_cpm_data = np.log10(cpm_counts.iloc[:, 1:] + 1)\n", "log_cpm_data.head()\n", "\n", "\n" ] }, { "cell_type": "code", "execution_count": 63, "id": "aquatic-species", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "SRX9845534 2017\n", "SRX9845541 2017\n", "SRX9845526 2015\n", "SRX9845530 2015\n", "SRX9845543 2017\n", "Name: batch, dtype: int64" ] }, "execution_count": 63, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Extract batch information\n", "batch = metadata_subset.set_index('Sample').loc[log_cpm_data.columns, 'batch']\n", "batch.head()\n" ] }, { "cell_type": "code", "execution_count": 72, "id": "opponent-tomorrow", "metadata": {}, "outputs": [ { "ename": "TypeError", "evalue": "'(slice(None, 2, None), slice(None, None, None))' is an invalid key", "output_type": "error", "traceback": [ "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", "\u001b[0;31mTypeError\u001b[0m Traceback (most recent call last)", "\u001b[0;32m\u001b[0m in \u001b[0;36m\u001b[0;34m\u001b[0m\n\u001b[1;32m 3\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 4\u001b[0m \u001b[0;31m# Apply Combat's fit_transform method\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 5\u001b[0;31m combat_adjusted_data = combat.fit_transform(\n\u001b[0m\u001b[1;32m 6\u001b[0m \u001b[0mY\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mlog_cpm_data\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mT\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;31m# Transpose: genes as rows, samples as columns\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 7\u001b[0m \u001b[0mb\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mbatch\u001b[0m \u001b[0;31m# Batch information\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m~/bin/anaconda3/lib/python3.8/site-packages/pycombat/pycombat.py\u001b[0m in \u001b[0;36mfit_transform\u001b[0;34m(self, Y, b, X, C)\u001b[0m\n\u001b[1;32m 311\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 312\u001b[0m \"\"\"\n\u001b[0;32m--> 313\u001b[0;31m \u001b[0;32mreturn\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mfit\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mY\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mb\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mX\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mC\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mtransform\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mY\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mb\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mX\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mC\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 314\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 315\u001b[0m \u001b[0;32mdef\u001b[0m \u001b[0m_validate_for_transform\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mY\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mb\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mX\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mC\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m~/bin/anaconda3/lib/python3.8/site-packages/pycombat/pycombat.py\u001b[0m in \u001b[0;36mfit\u001b[0;34m(self, Y, b, X, C)\u001b[0m\n\u001b[1;32m 173\u001b[0m \u001b[0;31m# Find intercepts\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 174\u001b[0m alpha_hat = np.matmul(sample_per_batch/float(n_samples),\n\u001b[0;32m--> 175\u001b[0;31m beta_hat[:n_batch, :])\n\u001b[0m\u001b[1;32m 176\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mintercept_\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0malpha_hat\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 177\u001b[0m \u001b[0;31m# Find slopes for covariates/effects\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m~/bin/anaconda3/lib/python3.8/site-packages/pandas/core/frame.py\u001b[0m in \u001b[0;36m__getitem__\u001b[0;34m(self, key)\u001b[0m\n\u001b[1;32m 3022\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mcolumns\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mnlevels\u001b[0m \u001b[0;34m>\u001b[0m \u001b[0;36m1\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 3023\u001b[0m \u001b[0;32mreturn\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_getitem_multilevel\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mkey\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 3024\u001b[0;31m \u001b[0mindexer\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mcolumns\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mget_loc\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mkey\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 3025\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mis_integer\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mindexer\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 3026\u001b[0m \u001b[0mindexer\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;34m[\u001b[0m\u001b[0mindexer\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32m~/bin/anaconda3/lib/python3.8/site-packages/pandas/core/indexes/base.py\u001b[0m in \u001b[0;36mget_loc\u001b[0;34m(self, key, method, tolerance)\u001b[0m\n\u001b[1;32m 3078\u001b[0m \u001b[0mcasted_key\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_maybe_cast_indexer\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mkey\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 3079\u001b[0m \u001b[0;32mtry\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m-> 3080\u001b[0;31m \u001b[0;32mreturn\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_engine\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mget_loc\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mcasted_key\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 3081\u001b[0m \u001b[0;32mexcept\u001b[0m \u001b[0mKeyError\u001b[0m \u001b[0;32mas\u001b[0m \u001b[0merr\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 3082\u001b[0m \u001b[0;32mraise\u001b[0m \u001b[0mKeyError\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mkey\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;32mfrom\u001b[0m \u001b[0merr\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", "\u001b[0;32mpandas/_libs/index.pyx\u001b[0m in \u001b[0;36mpandas._libs.index.IndexEngine.get_loc\u001b[0;34m()\u001b[0m\n", "\u001b[0;32mpandas/_libs/index.pyx\u001b[0m in \u001b[0;36mpandas._libs.index.IndexEngine.get_loc\u001b[0;34m()\u001b[0m\n", "\u001b[0;31mTypeError\u001b[0m: '(slice(None, 2, None), slice(None, None, None))' is an invalid key" ] } ], "source": [ "# Initialize Combat with desired mode ('p' for parametric or 'np' for non-parametric)\n", "combat = Combat(mode='p')\n", "\n", "# Apply Combat's fit_transform method\n", "combat_adjusted_data = combat.fit_transform(\n", " Y=log_cpm_data.T, # Transpose: genes as rows, samples as columns\n", " b=batch # Batch information\n", ")\n" ] }, { "cell_type": "code", "execution_count": 73, "id": "individual-tourism", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Shape of log_cpm_data.T: (44, 32707)\n", "Length of batch: 44\n", "Batch index: Index(['SRX9845534', 'SRX9845541', 'SRX9845526', 'SRX9845530', 'SRX9845543',\n", " 'SRX9845527', 'SRX9845506', 'SRX9845518', 'SRX9845501', 'SRX9845536',\n", " 'SRX9845529', 'SRX9845524', 'SRX9845507', 'SRX9845540', 'SRX9845517',\n", " 'SRX9845532', 'SRX9845503', 'SRX9845508', 'SRX9845537', 'SRX9845520',\n", " 'SRX9845523', 'SRX9845528', 'SRX9845516', 'SRX9845514', 'SRX9845522',\n", " 'SRX9845502', 'SRX9845509', 'SRX9845539', 'SRX9845542', 'SRX9845504',\n", " 'SRX9845538', 'SRX9845515', 'SRX9845519', 'SRX9845533', 'SRX9845505',\n", " 'SRX9845531', 'SRX9845510', 'SRX9845535', 'SRX9845513', 'SRX9845512',\n", " 'SRX9845525', 'SRX9845521', 'SRX9845500', 'SRX9845511'],\n", " dtype='object')\n", "log_cpm_data.T columns: Int64Index([ 0, 1, 4, 5, 6, 7, 9, 10, 11,\n", " 12,\n", " ...\n", " 38253, 38254, 38256, 38257, 38258, 38259, 38260, 38261, 38262,\n", " 38263],\n", " dtype='int64', length=32707)\n" ] } ], "source": [ "print(\"Shape of log_cpm_data.T:\", log_cpm_data.T.shape)\n", "print(\"Length of batch:\", len(batch))\n", "print(\"Batch index:\", batch.index)\n", "print(\"log_cpm_data.T columns:\", log_cpm_data.T.columns)\n" ] }, { "cell_type": "code", "execution_count": 74, "id": "threaded-station", "metadata": {}, "outputs": [], "source": [ "# Ensure `log_cpm_data.T` uses sample IDs as indices\n", "log_cpm_data.T.index = batch.index # Align sample IDs\n", "\n", "# Confirm alignment\n", "assert all(log_cpm_data.T.index == batch.index), \"Indices of log_cpm_data and batch do not match!\"\n" ] }, { "cell_type": "code", "execution_count": 76, "id": "initial-contest", "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "/home/shellytrigg/bin/anaconda3/lib/python3.8/site-packages/pycombat/pycombat.py:79: RuntimeWarning: divide by zero encountered in true_divide\n", " (abs(del_sq_post - del_sq_prior) / del_sq_prior).max())\n" ] } ], "source": [ "# Ensure `log_cpm_data.T` and `batch` are NumPy arrays\n", "log_cpm_data_array = log_cpm_data.T.values\n", "batch_array = batch.values\n", "\n", "# Apply Combat with parametric mode\n", "combat = Combat(mode='p')\n", "combat_adjusted_data = combat.fit_transform(\n", " Y=log_cpm_data_array, # Use NumPy array\n", " b=batch_array # Use NumPy array\n", ")\n", "\n", "# Convert adjusted data back to a DataFrame\n", "combat_adjusted_cpm = pd.DataFrame(\n", " combat_adjusted_data.T, # Transpose back to original format\n", " index=cpm_counts['gene_id'], \n", " columns=cpm_counts.columns[1:]\n", ")" ] }, { "cell_type": "code", "execution_count": 77, "id": "accessible-security", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
SRX9845534SRX9845541SRX9845526SRX9845530SRX9845543SRX9845527SRX9845506SRX9845518SRX9845501SRX9845536...SRX9845505SRX9845531SRX9845510SRX9845535SRX9845513SRX9845512SRX9845525SRX9845521SRX9845500SRX9845511
gene_id
LOC1110990290.1678130.1010850.057031-0.0059560.2651900.1576780.2122790.2374640.1003470.250462...0.5126050.2142170.0972880.3865920.3731230.2010970.4784930.1380320.1058570.393927
LOC1110990301.6412510.3915871.1208330.4028640.5542320.8600020.8558450.2621510.2776870.032037...0.1804131.4034500.4514240.7993011.1360610.4294880.8644260.1934751.1546060.592350
LOC1110990330.9815601.6343380.8872600.7488171.8858911.4424601.6908401.6049050.4954651.227480...1.5289341.5259670.8323501.2768971.5114691.2032141.2480221.2358800.9086720.871487
LOC1110990340.1022170.4179780.4302140.4669410.4267440.2938740.3903710.4951190.5131520.675999...0.2810710.3516170.5622330.5093280.4762610.4665770.5070690.4608420.4720700.377694
LOC1110990350.9996551.0865121.2274651.2420551.2502101.1539730.9513511.1744331.0792641.326385...1.0441421.2703511.0562121.2630981.3861461.1056081.1144251.0250931.1899940.947682
\n", "

5 rows × 44 columns

\n", "
" ], "text/plain": [ " SRX9845534 SRX9845541 SRX9845526 SRX9845530 SRX9845543 \\\n", "gene_id \n", "LOC111099029 0.167813 0.101085 0.057031 -0.005956 0.265190 \n", "LOC111099030 1.641251 0.391587 1.120833 0.402864 0.554232 \n", "LOC111099033 0.981560 1.634338 0.887260 0.748817 1.885891 \n", "LOC111099034 0.102217 0.417978 0.430214 0.466941 0.426744 \n", "LOC111099035 0.999655 1.086512 1.227465 1.242055 1.250210 \n", "\n", " SRX9845527 SRX9845506 SRX9845518 SRX9845501 SRX9845536 ... \\\n", "gene_id ... \n", "LOC111099029 0.157678 0.212279 0.237464 0.100347 0.250462 ... \n", "LOC111099030 0.860002 0.855845 0.262151 0.277687 0.032037 ... \n", "LOC111099033 1.442460 1.690840 1.604905 0.495465 1.227480 ... \n", "LOC111099034 0.293874 0.390371 0.495119 0.513152 0.675999 ... \n", "LOC111099035 1.153973 0.951351 1.174433 1.079264 1.326385 ... \n", "\n", " SRX9845505 SRX9845531 SRX9845510 SRX9845535 SRX9845513 \\\n", "gene_id \n", "LOC111099029 0.512605 0.214217 0.097288 0.386592 0.373123 \n", "LOC111099030 0.180413 1.403450 0.451424 0.799301 1.136061 \n", "LOC111099033 1.528934 1.525967 0.832350 1.276897 1.511469 \n", "LOC111099034 0.281071 0.351617 0.562233 0.509328 0.476261 \n", "LOC111099035 1.044142 1.270351 1.056212 1.263098 1.386146 \n", "\n", " SRX9845512 SRX9845525 SRX9845521 SRX9845500 SRX9845511 \n", "gene_id \n", "LOC111099029 0.201097 0.478493 0.138032 0.105857 0.393927 \n", "LOC111099030 0.429488 0.864426 0.193475 1.154606 0.592350 \n", "LOC111099033 1.203214 1.248022 1.235880 0.908672 0.871487 \n", "LOC111099034 0.466577 0.507069 0.460842 0.472070 0.377694 \n", "LOC111099035 1.105608 1.114425 1.025093 1.189994 0.947682 \n", "\n", "[5 rows x 44 columns]" ] }, "execution_count": 77, "metadata": {}, "output_type": "execute_result" } ], "source": [ "combat_adjusted_cpm.head()" ] }, { "cell_type": "code", "execution_count": 78, "id": "individual-therapy", "metadata": {}, "outputs": [], "source": [ "# Remove zero-variance genes\n", "non_zero_variance_genes = log_cpm_data.var(axis=0) > 0\n", "log_cpm_data_filtered = log_cpm_data.loc[:, non_zero_variance_genes]" ] }, { "cell_type": "code", "execution_count": 79, "id": "hairy-blues", "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "/home/shellytrigg/bin/anaconda3/lib/python3.8/site-packages/pycombat/pycombat.py:79: RuntimeWarning: divide by zero encountered in true_divide\n", " (abs(del_sq_post - del_sq_prior) / del_sq_prior).max())\n" ] } ], "source": [ "combat_adjusted_data = combat.fit_transform(\n", " Y=log_cpm_data_filtered.T.values,\n", " b=batch.values\n", ")" ] }, { "cell_type": "code", "execution_count": 84, "id": "sharp-republic", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "True\n" ] } ], "source": [ "# Check for NaN or infinite values in the adjusted data\n", "print(np.isfinite(combat_adjusted_data).all()) # Should return True" ] }, { "cell_type": "code", "execution_count": 81, "id": "primary-telescope", "metadata": {}, "outputs": [], "source": [ "# Remove low-variance genes\n", "low_variance_threshold = 1e-6 # Define a small threshold\n", "non_low_variance_genes = log_cpm_data.var(axis=0) > low_variance_threshold\n", "log_cpm_data_filtered = log_cpm_data.loc[:, non_low_variance_genes]\n" ] }, { "cell_type": "code", "execution_count": 82, "id": "thick-skirt", "metadata": {}, "outputs": [], "source": [ "# Clip values to avoid numerical instability\n", "log_cpm_data_clipped = log_cpm_data_filtered.clip(lower=1e-6)\n" ] }, { "cell_type": "code", "execution_count": 83, "id": "popular-healthcare", "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "/home/shellytrigg/bin/anaconda3/lib/python3.8/site-packages/pycombat/pycombat.py:79: RuntimeWarning: divide by zero encountered in true_divide\n", " (abs(del_sq_post - del_sq_prior) / del_sq_prior).max())\n" ] } ], "source": [ "combat_adjusted_data = combat.fit_transform(\n", " Y=log_cpm_data_clipped.T.values,\n", " b=batch.values\n", ")" ] }, { "cell_type": "code", "execution_count": 85, "id": "molecular-health", "metadata": { "scrolled": true }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAA1kAAALICAYAAACNTDphAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAABt9UlEQVR4nO3dd5xcdb3/8fdnZna2ZNN7QkjfbDYJxYREA5gEkdCliGABLCBe9XKx159SFEEFuVyviohXEJQiGikBpJgAAoaEkp6QAOm9bLbvzs7n98fMhk2yLdkzZXdfz8djH5n5nu8585n5Jrv7zvme7zF3FwAAAAAgGKFMFwAAAAAAnQkhCwAAAAACRMgCAAAAgAARsgAAAAAgQIQsAAAAAAgQIQsAAAAAAkTIAoDDZGY/MrOdZrY107UEwczczMZkuo4Gjesxs9+Y2f/LdE0AABwOQhaArGdm75pZlZmVm9k2M/s/MytstH22mT1vZmVmtsPM5pvZuQcdY2byl/dvtvE1R5pZ3Mx+dVD7MElfk1Ti7oPM7NNm9mIQ7/Og1xlsZneZ2Zbk+1ppZteZWbegX6s9kmNzahv6Nfl5tsbdv+DuNxx5hfvHfmN7jtHO159nZldk6vUBAOlHyALQUZzj7oWS3ifpBEnflyQz+6ikhyTdI+koSQMl/UDSOQftf7mk3ck/2+IySXskXWJmuY3ah0va5e7bj/B9HMDMIk209ZH0sqR8SR9w9+6SPiypl6TRQbxuBjT3eQIA0OkQsgB0KO6+SdITkiaamUm6VdIN7v47dy9197i7z3f3Kxv2MbMCSR+V9CVJY81sShte6jIlglydkoEtecbmaUlDkmfVHpD0G0kfSD7fm+yXa2Y/N7P1yTNvvzGz/OS2mWa20cy+lZxu+H9NvPZXJZVJ+pS7v5t83xvc/b/cfXHyONPN7FUzK03+Ob3R+52XnNL4UrKuR82sr5ndZ2b7kv1HHPSaZ5rZ28lpkD8zs1DyWKPN7Dkz25Xcdp+Z9Upu+6OkoyU9mnydls4SHvJ5Nqr3G8kzdpvN7LMHbfuDmf0o+fiQs4YHTS0808yWJ8/8bTKzryfP/D3RaMzKzWyImV1rZg+Z2b3J/kvMrMjMvmNm281sg5md1uh1ejY6s7gp+fmGG9eVHPM9ZvaOmZ2R3PZjSSdL+mXytX/ZwmcEAOgkCFkAOpTkdL0zJb0uaZykYZL+0spuF0oqV+KM11NK/MLf0mucrMRZsfslPdjQ392fkXSGpM3uXujuF0v6gqSXk897JQ9xs6QiScdJGiNpqBJn1xoMktRHibNin2+ihFMl/dXd483U10fS45Jul9RXiaD5uJn1bdTtEkmXJl97tBJnxv4v+borJP3woMOeL2mKEmcKPyKpIeyYpJ9IGiJpvBKf97XJz+NSSeuVPMvo7j9tpt4mP8/kttMlfV2JM3Vjk+/9SN0l6arkmb+Jkp5z9wodOGaF7r452f8cSX+U1FuJv09PKfFzcaik6yXd0ejYd0uKKTGex0s6TVLjKYDTJK2S1E/STyXdZWbm7t+T9IKkLydf+8vteH8AgA6CkAWgo5iTPFP0oqT5km5UImBI0pZW9r1c0gPuXi/pT5I+bmY5rfR/wt33JPufYWYD2lJk8uzalZK+4u673b0sWesljbrFJf3Q3WvcvaqJw/Rt5T2dJektd/+ju8fc/c+SVurAM0T/5+5r3b1UiTM5a939GXePKRE2jz/omDcn610v6TZJH5ckd1/j7k8na92hRKCb0ZbPopGWPs+PJWtdmgxE1x7msRurk1RiZj3cfY+7v9ZK/xfc/alGn0l/STe5e50SgXCEmfUys4FKBLVr3L0iOVX0FzpwTNe5+53Jv2N3SxqsxNRVAEAXRMgC0FGc5+693H24u38xGU52JbcNbm6n5JmvWZLuSzb9XVKeEkGlqf75ki5q6O/uLytxtuYTbayzv6QCSYvMbG8yGD6ZbG+ww92rWzjGLrXwnpQ4q7TuoLZ1SpyBabCt0eOqJp4X6kAbDjrWEEkyswFmdn9yitw+SfcqcbamTdrweQ5p4rWP1IVKnOVcZ4nFTz7QSv+DP5OdyZDU8FxKfE7DJeVI2tJoTO+Q1Dh4719p0t0rG+0LAOiCCFkAOrJVSvyCfmELfS5V4nvdo8lroN5WImQ1N2XwfEk9JP3KzLYm9xnaQn8/6PlOJX5Bn5AMhb3cvWdy0Y7m9jnYM5LOb7guqgmblfjFv7GjJW1q5bgtGXbQsRqm1P1EiXqPcfcekj6lxBTCBq29l9Y+zy1NvHZzKpQIsJIkMxvUeKO7v+ruH1Ei/MxRYmpiW2pszQZJNZL6NRrTHu4+oY37t/f1AQAdDCELQIfl7q7EIhH/z8w+Y2Y9zCxkZieZ2W+T3S6TdJ0S10c1fF0o6ayDrmFqcLmk30ua1Kj/iZKOM7NJTfTfJukoM4sma4pLulPSLxqmxJnZUDObfRhv7VYlgsndZja80TFuNbNjJM2VVGRmnzCziJldLKlE0mOH8RoH+4aZ9U6e+fsvSQ8k27srcT3bXjMbKukbB+23TdKoFo7b2uf5oKRPm1mJJRYoOfhascbelDTBzI4zszw1mlpoZlEz+6SZ9UxO99snqeGs1DZJfc2sZ0sfQHPcfYukf0i6pdHfsdFm1tZpk619RgCAToaQBaBDc/e/SLpYiYUaNivxC+2PJP3dzN4vaYSk/3X3rY2+HpG0RsnrjhokQ8SHJN12UP9FSkz5a2r59+ckLZO01cx2Jtu+lTz+K8kpds8osUhHW9/TbknTlbjG6N9mVibpWUmlkta4+y5JZytxv65dkr4p6Wx339nMIdvi75IWSXpDiUU17kq2X6fEYhilyfa/HrTfTyR9PzmN7uuNN7Tl83T3J5S4Buw5JT6z55or0N1XK7EgxTOS3lLi+rzGLpX0bvIz/4ISZ93k7isl/VnS28k6h7TpEznQZZKikpYrsRT9X9TylM7G/lvSR5MrD95+BK8NAOhgLPEfwQAAZB8zu0eJYHl9pmsBAKCtOJMFAMhKlrhR8zhJ72S6FgAADgchCwCQrbZK2ivp4QzXAQDAYWG6IAAAAAAEiDNZAAAAABCgSKYLCEK/fv18xIgRrfarqKhQt27dUl8QWsQ4ZB5jkB0Yh8xjDLID45B5jEH7LVq0aKe792+9J7qCThGyRowYoYULF7bab968eZo5c2bqC0KLGIfMYwyyA+OQeYxBdmAcMo8xaD8zW5fpGpA9mC4IAAAAAAEiZAEAAABAgAhZAAAAABCgTnFNFgAAAICmLVq0aEAkEvmdpIniJEsQ4pKWxmKxKyZPnry9qQ6ELAAAAKATi0Qivxs0aND4/v377wmFQtwkt53i8bjt2LGjZOvWrb+TdG5TfUiyAAAAQOc2sX///vsIWMEIhULev3//UiXODDbdJ431AAAAAEi/EAErWMnPs9ksRcgCAAAAgAARsgAAAAAELhwOTy4uLi4ZN25cSUlJyfinn366W0v9d+7cGb7pppv6t3bcqVOnjnv++ecLgqs0eIQsAAAAAIHLzc2Nr1y5cvmqVauW33DDDZu++93vHtVS/127doXvuuuuAemqL5UIWQAAAABSqrS0NNyzZ89Y8nHoAx/4QFFJScn4oqKiknvvvbeXJH3ta187asOGDbnFxcUlV1111VGS9P3vf39gUVFRybhx40q++MUvDm043p///OfekyZNGj9ixIiJTz75ZGFG3lQLWMIdAAAAQOBqampCxcXFJTU1NbZz586cuXPnrpakgoKC+OOPP76mT58+8S1btkSmTZtW/IlPfGLvLbfcsvHss8/OX7ly5XJJevDBB3s8/vjjvRctWrSye/fu8W3btoUbjh2LxWzJkiUrHnjggZ7XX3/9kNNPP311pt5nUziTBQAAACBwDdMF33nnnWV/+9vf3vrMZz4zMh6PKx6P2zXXXHNUUVFRyaxZs4q2b98e3bhx4yEnf55++uken/rUp3Z27949LkkDBw6sb9h20UUX7ZGk6dOnV2zcuDGavnfVNpzJAgAAAJBSp556asWePXsiW7ZsiTz88MM9d+3aFVmyZMmK3NxcHzp06KSqqqpDTv64u8ysyePl5eW5JEUiEdXX1zfdKYM4kwUAAAAgpV5//fW8eDyugQMHxkpLS8P9+vWry83N9UcffbT75s2bo5LUs2fP+oqKiv355PTTT9/3xz/+sV9ZWVlIkhpPF8x2nMkCAAAAELiGa7KkxFmpX//61+9GIhFdccUVu88444wxEydOHD9hwoTKkSNHVkvSoEGD6idPnlw+duzYCaecckrpHXfcsfG1114rOO6448bn5OT4qaeeWvrLX/5yU2bfVdsQsgAAAAAErr6+flFT7YMHD4698cYbK5va9uijj77T+PmNN9649cYbb9zauG3BggWrGh9r06ZNS4KoN0hMFwQAAACAABGyAAAAACBAhCwAAAAACBAhCwAAAAACRMgCAAAAgAARsgAAAAAgQIQsAAAAACm1atWq6NixYye0tf/tt9/e9913381prc9ll112dPurCx73yQIAAMgilbFa1dbHmtwWDUdUEImmuSIg/e69995+xx13XNWIESPqMl3LkSBkAQAAZJHa+pi+tWBOk9tunnoeIQsdViwW0wUXXDBi6dKlBaNGjap+6KGH3r3uuusGPvnkk71qampCU6ZMKb/vvvvW3X333b2XLl1acNlll43Ky8uLL1y4cMWiRYvyr7nmmqMrKytD0WjUn3/++VWStHXr1pyTTz557Pr163PPOOOMvb/5zW82Zvp9SkwXBAAAAJAG7777bt4XvvCFHatXr17evXv3+M9+9rP+3/jGN7YvXbp0xVtvvbWsqqoqdP/99/f8zGc+s2fixImV99xzz9srV65cHolE9MlPfnL0bbfdtn7VqlXL58+fv6qwsDAuScuXLy+YM2fO2ytWrFj2yCOP9F6zZk2LUwzThZAFAAAAIOUGDRpUe9ppp1VI0qWXXrrrpZdeKnziiSe6H3PMMcVFRUUlL730UvelS5fmH7zf4sWL8wYMGFA3Y8aMSknq06dPPCcnkaVOOumkfX379q0vKCjwMWPGVK9duzY3rW+qGUwXBAAAAJByZnbI86997WvD//3vfy8fM2ZM3Ve/+tUh1dXVh5wEcneZmTd1zGg0ur89HA57XV2dNdUv3TiTBQAAACDltmzZEn3mmWe6SdKf/vSnPtOnTy+XpEGDBsVKS0tDjz76aO+GvoWFhfWlpaVhSTr22GOrt23bFp0/f36BJO3ZsydUV5fd62EQsgAAAACk3KhRo6p///vf9y0qKirZs2dP5Otf//qOT37ykztKSkomnHHGGWOOPfbYioa+l1122c7//M//HF5cXFwSi8V03333rb366quPHjduXMnMmTOLKisrszrHMF0QAAAAQEqNGzeudu3atcsObr/99ts333777ZsPbv/0pz+999Of/vTehuczZsyofPPNN1c27nP11VfvkrSr4fk///nPNcFWfeQIWQAAAFkkGo7o5qnnNbsNQPbjXyoAAEAWKYhEuRcW0MFl9VxGAAAAAOhoCFkAAAAAECBCFgAAAAAEiJAFAAAAAAEiZAEAAABIqTVr1uRMmzataNSoURPGjBkz4YYbbhggSdu2bQtPnz597PDhwydOnz597I4dO8KStHXr1vC0adOKCgoKjr/sssuObnysqVOnjhsxYsTE4uLikuLi4pJNmzZl3WJ+hCwAAAAAKZWTk6Nbbrll49tvv73s1VdfXXHXXXcNWLRoUd4Pf/jDwTNnzixbt27d0pkzZ5b94Ac/GCRJBQUFfv3112++9tprNzZ1vHvuueftlStXLl+5cuXyoUOHxtL7blqX8tRnZsMk3SNpkKS4pN+6+3+bWR9JD0gaIeldSR9z9z3Jfb4j6XOS6iVd7e5PpbpOAAAAANL8zW/1eXzDkqGltdXRntG82rOGTdo0Y8jY3e055vDhw+uGDx9eJ0m9e/eOjx49umr9+vXRJ598stf8+fNXSdJVV121a8aMGeMkberRo0d89uzZ5atWrcoN4C2lXTrOZMUkfc3dx0t6v6QvmVmJpG9Letbdx0p6NvlcyW2XSJog6XRJvzKzcBrqBAAAALq0+Zvf6vPQO4uGl9ZWRyWptLY6+tA7i4bP3/xWn6BeY9WqVdHly5cXzJgxo3zXrl2RhvA1fPjwut27d7fpJNAVV1wxori4uOQb3/jG4Hg8HlRpgUl5yHL3Le7+WvJxmaQVkoZK+oiku5Pd7pZ0XvLxRyTd7+417v6OpDWSpqa6TgAAAKCre3zDkqF18fgBGaEuHg89vmHJ0CCOX1paGrrgggtG33TTTRv69OlzROnogQceeHv16tXLX3755ZUvvfRS4a9+9au+QdQWpLReJGZmIyQdL+nfkga6+xYpEcTMbECy21BJrzTabWOy7eBjfV7S5yVp4MCBmjdvXquvX15e3qZ+SC3GIfMYg+zAOGQeY5AdGIfMYwzQoOEMVlvbD0dNTY2dddZZoy+66KLdl19++V5J6tu3b2zdunU5w4cPr1u3bl1Onz59Wr2+auTIkfunHV588cW7FyxY0E3SrvbWF6S0hSwzK5T0sKRr3H2fmTXbtYk2P6TB/beSfitJU6ZM8ZkzZ7Zaw7x589SWfkgtxiHzGIPswDhkHmOQHRiHzGMM0KBnNK+2qUDVM5pX257jxuNxXXLJJcOLioqqr7322m0N7bNnz957xx139L3xxhu33nHHHX1PP/30vS0dp66uTjt37owMHjw4VlNTY3Pnzu15yimnlLWntlRIS8gysxwlAtZ97v7XZPM2MxucPIs1WNL2ZPtGScMa7X6UpM3pqBMAAADoys4aNmnTQ+8sGt54ymBOKBQ/a9ikTe057tNPP104Z86cvmPHjq0qLi4ukaTrrrtu03XXXbfl/PPPHz18+PB+Q4YMqZ0zZ87ahn2GDh06qby8PFxXV2dPPfVUr7lz564eO3Zs7amnnjq2rq7O4vG4nXzyyfu++tWv7mhPbamQjtUFTdJdkla4+62NNj0i6XJJNyX//Huj9j+Z2a2ShkgaK2lBqusEAAAAurqGVQSDXl1w9uzZ5e6+qKltL7/88uqm2jdt2rSkqfZly5ataE8t6ZCOM1knSrpU0hIzeyPZ9l0lwtWDZvY5SeslXSRJ7r7MzB6UtFyJlQm/5O71aagTAAAA6PJmDBm7u72hqqtLechy9xfV9HVWkvShZvb5saQfp6woAAAAAEiRdNwnCwAAAAC6DEIWAAAAAASIkAUAAAAAASJkAQAAAECACFkAAAAAUmrNmjU506ZNKxo1atSEMWPGTLjhhhsGSNK2bdvC06dPHzt8+PCJ06dPH7tjx46wJG3dujU8bdq0ooKCguMvu+yyoxuOs2fPnlBxcXFJw1fv3r2P/exnPzusudfNFEIWAAAAgJTKycnRLbfcsvHtt99e9uqrr6646667BixatCjvhz/84eCZM2eWrVu3bunMmTPLfvCDHwySpIKCAr/++us3X3vttRsbH6d3797xlStXLm/4GjJkSO1FF120JzPvqnnpuE8WAAAAgA4i/uY/+/grjw5VRWlU3XrW2vvP2RQ6dla77ps1fPjwuuHDh9dJiaA0evToqvXr10effPLJXvPnz18lSVddddWuGTNmjJO0qUePHvHZs2eXr1q1Kre5Yy5ZsiR3165dObNnzy5vT22pQMgCAAAAICkZsOY9MFz1dYkZbxWlUZ/3wPC4pPYGrQarVq2KLl++vGDGjBnlu3btijSEr+HDh9ft3r27zfnk7rvv7nPuuefuDoWyb3Je9lUEAAAAICP8lUeH7g9YDerrQv7Ko0ODOH5paWnoggsuGH3TTTdt6NOnT7w9x/rb3/7W59JLLw0k+AWNkAUAAAAgoaI0eljth6GmpsbOOuus0RdddNHuyy+/fK8k9e3bN7Zu3bocSVq3bl1Onz59Ym051ssvv5xfX19vJ598cmV760oFQhYAAACAhG49aw+rvY3i8bguueSS4UVFRdXXXnvttob22bNn773jjjv6StIdd9zR9/TTT9/bluP98Y9/7HP++edn5VksiWuyAAAAACTZ+8/ZdMA1WZIUzonb+8/Z1J7jPv3004Vz5szpO3bs2Kri4uISSbruuus2XXfddVvOP//80cOHD+83ZMiQ2jlz5qxt2Gfo0KGTysvLw3V1dfbUU0/1mjt37urJkydXS9IjjzzS59FHH32rPTWlEiELAAAAgKTE4hZxJa/NCnB1wdmzZ5e7+6Kmtr388surm2rftGnTkuaOt3Hjxma3ZQNCFgAAAID9QsfO2q2AVhLsqrgmCwAAAAACRMgCAAAAgAARsgAAAAAgQIQsAAAAAAgQIQsAAAAAAkTIAgAAAJAyO3fuDN900039s+1YqXwdQhYAAACAlNm1a1f4rrvuGnBweywWC+xYQWvv6xCyAAAAAKTM1772taM2bNiQW1xcXDJx4sTx06ZNKzrnnHNGjhs3bkIsFtNVV1111MSJE8cXFRWV/OxnP+snSaWlpaEPfOADRSUlJeOLiopK7r333l4HH+uqq6466rHHHut+wgknjDvzzDNHjRgxYuIXv/jFob/+9a/7TJo0aXxRUVHJsmXLciVp8+bNkdmzZ4+eOHHi+IkTJ47/xz/+0U2SvvrVrw656KKLRkydOnXcUUcdNelHP/rRgKZe53DfMzcjBgAAAJAyt9xyy8azzz47f+XKlcsfe+yx7hdddNGY119/fVlxcXHtz3/+8349e/asX7p06Yqqqio74YQTis8555x9o0ePrn388cfX9OnTJ75ly5bItGnTij/xiU/sbXwsSXrssce6r1y5Mv8vf/nL2wMGDIgNHz58Um5u7s4lS5asuOGGGwbccsstA37/+99vuOqqq4Z99atf3TZ79uzyt956Kzp79uyxb7/99jJJWrNmTd5LL720au/eveHx48dP/MY3vrHj4Nc5XIQsAAAAAGlzzDHHVBQXF9dK0jPPPNNj5cqVBY888khvSSorKwsvX748b+TIkXXXXHPNUa+88kphKBTS9u3boxs3bmwyu0yaNKli+PDhdZJ09NFH15xxxhmlknTsscdWzZ8/v7sk/etf/+rx1ltv5TfsU15eHt6zZ09Ikk477bS9+fn5np+fH+vTp09dc69zOAhZAAAAANKmoKAg3vDY3e2WW25Zf+GFF+5r3Of222/vu2vXrsiSJUtW5Obm+tChQydVVVU1ealTbm6uNzwOhULKy8vzhsf19fWWfB0tXLhwRWFhobe0fzgcViwWs/a+R67JAgAAAJAyPXv2rK+oqGgyd3z4wx8u/fWvf92/pqbGJGnx4sW5+/btC5WWlob79etXl5ub648++mj3zZs3R1s7VktOOumkfTfffPP+hSxeeuml/Jb6H+nrNOBMFgAAAICUGTRoUP3kyZPLx44dOyE3Nzfev3//uoZtX/nKV3a+++67uZMmTRrv7tanT5+6uXPnrr3iiit2n3HGGWMmTpw4fsKECZUjR46sPvhYp5xySuk555xT2pYafvvb32644oorji4qKiqpr6+3adOmlU2fPn19W2o+5ZRTSu+4446Nh/OeCVkAAAAAUurRRx99p6n2cDisX/7yl5skbTp42xtvvLGyLcc6++yzyxoeL1iwYFXj9oZtgwcPjj3++ONvH3ysW2+9dXPj52+99day1mpuC6YLAgAAAECACFkAAAAAECBCFgAAAAAEiJAFAAAAAAEiZAEAAABAgAhZAAAAABAgQhYAAACAlFu/fn3k7LPPHjVs2LCJo0ePnjBjxowxixcvzj3c41x//fUDysrKDjvHFBQUHH+4+xwpQhYAAACAlIrH4zr33HPHfPCDHyzbsGHD0rVr1y77yU9+smnz5s05h3usO+64Y2B5eXmTOSYWi7W/2AAQsgAAAADs9/DDD/c5/fTTJ02ZMmXy6aefPunhhx/u095jPvbYY90jkYh/85vf3NHQNn369KrTTjut/Kqrrjpq7NixE4qKikruvPPO3g39p06dOu70008fNXLkyAnnnnvuyHg8rh/96EcDtm/fnjNjxoyiadOmFUmJM1TXXHPNkGOOOab42WefLbz22msHjh07dsLYsWMnXH/99QPaW/uRiGTiRQEAAABkn4cffrjPLbfcMry2tjYkSTt37ozecsstwyXpwgsv3H2kx128eHH+scceW3lw+z333NNryZIl+StWrFi2ZcuWyNSpU8efdtpp5ZK0YsWK/DfeeOPtESNG1E2ePLn46aefLvz+97+//de//vXA+fPnrx48eHBMkqqqqkITJ06suu222za/8MILBX/605/6Llq0aIW7a/LkyeM/9KEPlZ144olVR1r7keBMFgAAAABJ0p133jm0IWA1qK2tDd15551DU/F6L7zwQvePfexjuyORiIYNGxabNm1a+YsvvlggSZMmTaoYPXp0XTgc1oQJEyrXrl0bbeoY4XBYn/70p/dI0rx58wrPPPPMvT169Ij37NkzftZZZ+355z//2T0VtbeEkAUAAABAUuLM1eG0t9WkSZOq3nzzzYKD29292X1yc3P3bwyHw4rFYtZUv2g0Go9EIq0eL50IWQAAAAAkSf369as9nPa2Ouecc8pqa2vtlltu6dfQNn/+/ILevXvH/vKXv/SJxWLavHlzZMGCBYUnn3xyRUvH6tatW31paWmTOeaUU04pnzt3bq+ysrLQvn37QnPnzu09a9assvbUfiQIWQAAAAAkSVdeeeWmaDQab9wWjUbjV1555ab2HDcUCumRRx5Z++yzz/YYNmzYxDFjxkz44Q9/OOTTn/707gkTJlSNHz9+wsyZM4uuu+66jUcffXSLSwRefvnlO88444yxDQtfNHbSSSdVfuITn9j1vve9b/zkyZPHX3rppTvSfT2WxMIXAAAAAJIaFre48847h+7cuTPar1+/2iuvvHJTexa9aDBixIi6uXPnvn1w+x133LFR0sbGbWeffXbZ2Wefvf8M1D333LO+4fH3vve97d/73ve2NzyvrKx8vfG+11577bZrr71228Gvc3C/VCJkAQAAANjvwgsv3B1EqOrKmC4IAAAAAAEiZAEAAABAgAhZAAAAABAgQhYAAAAABIiQBQAAAAABImQBAAAASKmCgoLjW+tz8cUXD1+0aFHe4R77pZdeyn/ggQd6Hu5+U6dOHff8888XHO5+bUHIAgAAAJBxDzzwwLrJkydXH+5+CxcuLHj88ccPO2SlEiELAAAAgCSptrbWrrjiirFXXHHF2PLy8lDD49raWgvi+I899lj3qVOnjjv99NNHjRw5csK55547Mh6PSzrwzNJf//rXHscdd1xxSUnJ+DPOOGNUaWlpSJLmz59fcPzxxxePGzeuZNKkSeN37doV/slPfjLk0Ucf7V1cXFxy55139t63b1/ooosuGjFx4sTx48ePL7n33nt7SVJ5ebmdffbZo4qKikrOOuusUdXV1YG8p6ZwM2IAAAAAkqQvfvGLY5YvX14oSWeeeeYxsVjMGtp/97vfvRXEa6xYsSL/jTfeeHvEiBF1kydPLn766acLZ8+eXd6wfcuWLZEbb7xx8PPPP7+6R48e8e9973uDbrjhhoE/+tGPtn7yk58cfd99962dMWNG5e7du0Pdu3ePf+c739m8cOHCbvfcc896Sfryl788dNasWfseeuihd3fu3BmeMmXK+HPPPXffrbfe2j8/Pz++evXq5f/+97/zTzzxxJIg3k9TCFkAAAAADlBbWxuqra2VJEWj0XiQx540aVLF6NGj6yRpwoQJlWvXro023j5v3rxua9euzZs6dWqxJNXV1dnkyZPLFy9enDdgwIC6GTNmVEpSnz59mqxr3rx5PZ566qlet99++yBJqqmpsTVr1kRffPHFwquvvnq7JE2bNq2qqKioMsj31RghCwAAAIAk6bbbblt75plnHtMQsCQpEon4f//3f68N6jVyc3O94XE4HFbD2bIG7q6TTjpp36OPPvpO4/Z///vf+WbmaoW76y9/+cuaY489tubgbWYpmyF4AK7JAgAAACBJuuaaa0YfHHpisZj913/91+h01TBz5syKhQsXFi5dujRXksrKykKLFy/OPfbYY6u3bdsWnT9/foEk7dmzJ1RXV6cePXrUl5eX7881s2bN2nfLLbcMbLjW61//+le+JJ100knl9957bx9JevXVV/NWr16dkpUFJUIWAAAAgINEo9F4QUFBfdBTBVtjZhoyZEjsjjvuePeSSy4ZVVRUVDJ58uTiJUuW5OXl5fl999239uqrrz563LhxJTNnziyqrKwMnXHGGWWrV6/Ob1j44qabbtoci8WsuLi4ZOzYsRO+//3vD5Wkr3/969srKirCRUVFJTfeeOOgSZMmVaTqfTBdEAAAAIAk6Ve/+tWaL37xi2OkxNTBa665ZnRDe3uOW1lZ+boknX322WVnn312WUN7w2IVkrR3795w//79Y5J07rnnlp177rkrDj7OjBkzKt98882Vjdt69uyppUuXHtD3T3/607qD9y0sLPTHHnvs7fa8j7YiZAEAAACQJEWjUW+8imBQKwq2Zvr06WPHjRtXVVxcXNt67+xHyAIAAACQUS+99FJawly6cE0WAAAAAASIkAUAAAAAASJkAQAAAECACFkAAAAAECBCFgAAAICUCofDk4uLi0savlatWhVty36rVq2Kjh07dkKq6wsaqwsCAAAASKnc3Nz4ypUrl2e6jnThTBYAAACAtHvhhRcKTjjhhHETJkwYf9JJJ41dt25dTkP7uHHjSo477rjiW2+9dUCm6zwSnMkCAAAAsN9JJ510fHV19f6TMXl5efEXX3zx9fYcs6amJlRcXFwiScOGDat57LHH3r766quPfvzxx9cMGTIkduedd/b++te/PvShhx5693Of+9yIX/ziF+vPOuus8quuuuqo9r6fTCBkAQAAANivccBq6vmROHi64Kuvvpr31ltv5Z9yyilFkhSPx9W/f/+6Xbt2hcvKysJnnXVWuSR99rOf3fXcc8/1bO/rpxshCwAAAEBaubuNGTOm6o033ljZuH3nzp1hM8tUWYHhmiwAAAAAaXXMMcdU7969O/LMM890k6SamhpbuHBhXr9+/eoLCwvrn3rqqUJJ+sMf/tAns5UeGUIWAAAAgP3y8vLiLT0P6DX8/vvvX/vtb3/7qHHjxpVMmDChZP78+YWSdNddd7179dVXH33ccccV5+fne9CvnQ5MFwQAAACwX3sXuWhKZWXlIcecPn161cKFC1cd3H7yySdXrlq1av/1W7feeuvmoOtJNc5kAQAAAECACFkAAAAAECBCFgAAANC5xePxeMdfsi+LJD/PZq9VI2QBAAAAndvSHTt29CRoBSMej9uOHTt6SlraXB8WvgAAAAA6sVgsdsXWrVt/t3Xr1oniJEsQ4pKWxmKxK5rrQMgCAAAAOrHJkydvl3RupuvoSkiyAAAAABAgQhYAAAAABIiQBQAAAAABImQBAAAAQIAIWQAAAAAQIEIWAAAAAASIkAUAAAAAASJkAQAAAECACFkAAAAAECBCFgAAAAAEiJAFAAAAAAEiZAEAAABAgAhZAAAAABAgQhYAAAAABIiQBQAAAAABImQBAAAAQIAIWQAAAAAQoJSHLDP7vZltN7OljdquNbNNZvZG8uvMRtu+Y2ZrzGyVmc1OdX0AAAAAEKR0nMn6g6TTm2j/hbsfl/yaK0lmViLpEkkTkvv8yszCaagRAAAAAAKR8pDl7s9L2t3G7h+RdL+717j7O5LWSJqasuIAAAAAIGCZvCbry2a2ODmdsHeybaikDY36bEy2AQAAAECHYO6e+hcxGyHpMXefmHw+UNJOSS7pBkmD3f2zZva/kl5293uT/e6SNNfdH27imJ+X9HlJGjhw4OT777+/1TrKy8tVWFgYzJvCEWMcMo8xyA6MQ+YxBtmBccg8xqD9Zs2atcjdp2S6DmSHSCZe1N23NTw2szslPZZ8ulHSsEZdj5K0uZlj/FbSbyVpypQpPnPmzFZfd968eWpLP6QW45B5jEF2YBwyjzHIDoxD5jEGQLAyMl3QzAY3enq+pIaVBx+RdImZ5ZrZSEljJS1Id30AAAAAcKRSfibLzP4saaakfma2UdIPJc00s+OUmC74rqSrJMndl5nZg5KWS4pJ+pK716e6RgAAAAAISspDlrt/vInmu1ro/2NJP05dRQAAAACQOplcXRAAAAAAOh1CFgAAAAAEiJAFAAAAAAEiZAEAAABAgAhZAAAAABAgQhYAAAAABIiQBQAAAAABImQBAAAAQIAIWQAAAAAQIEIWAAAAAASIkAUAAAAAASJkAQAAAECACFkAAAAAECBCFgAAAAAEiJAFAAAAAAEiZAEAAABAgAhZAAAAABCgSKYLAACgM6uM1aq2Ptbktmg4ooJINM0VAQBSjZAFAEAK1dbH9K0Fc5rcdvPU8whZANAJMV0QAAAAAAJEyAIAAACAABGyAAAAACBAhCwAAAAACBAhCwAAAAACxOqCAACkUDQc0c1Tz2t2GwCg8+G7OwAAKVQQibJMOwB0MUwXBAAAAIAAEbIAAAAAIECELAAAAAAIECELAAAAAAJEyAIAAACAABGyAAAAACBAhCwAAAAACBAhCwAAAAACRMgCAAAAgAARsgAAAAAgQIQsAAAAAAgQIQsAAAAAAkTIAgAAAIAAEbIAAAAAIECELAAAAAAIECELAAAAAAJEyAIAAACAABGyAAAAACBAhCwAAAAACBAhCwAAAAACRMgCAAAAgAARsgAAAAAgQIQsAAAAAAgQIQsAAAAAAkTIAgAAAIAAEbIAAAAAIECELAAAAAAIECELAAAAAAJEyAIAAACAABGyAAAAACBAhCwAAAAACBAhCwAAAAACRMgCAAAAgAARsgAAAAAgQIQsAAAAAAgQIQsAAAAdXn08rupYnerj8UyXAiiS6QIAAACAI1UZq1VlrFbzN7+l7dVl6pdXqJmDx6pbJKqCnNxMl4cuipAFAACADqmyrlbPbl6px9YvPaD9mU0r9eGhxTrj6AnqFiFoIf2YLggAAIAO6a192w8JWA2e3rRSy3ZvkbunuSqAkAUAAIAOqKyuWo+sW9xin8fWL1F5XU2aKgLeQ8gCAABAh+MubazY22KfbVVlqudMFjKAkAUAAIAOydrSpy2dgIARsgAAANDhhM1U3GtQi31G9+iviPHrLtKPv3UAAADocLrl5OrCkcfJmjmfZZIuHHmcurGMOzKAkAUAAIAOqX9ed31pwgdVEIke0J4XztGVxSdpSEGvzBSGLo/7ZAEAAKBDyovkaHyvQbp+8tlaV75bO6rL1C+3UCO691VeOKKcML/qIjP4mwcAAIAOKxIKq3s0rIl9hmS6FGA/pgsCAAAAQIAIWQAAAAAQIEIWAAAAAASIkAUAAAAAASJkAQAAAECACFkAAAAAECBCFgAAAAAEiJAFAAAAAAEiZAEAAABAgAhZAAAAABAgQhYAAAAABIiQBQAAAAABImQBAAAAQIAIWQAAAAAQIEIWAAAAAASIkAUAAAAAASJkAQAAAECACFkAAAAAECBCFgAAAAAEiJAFAAAAAAEiZAEAAABAgAhZAAAAABAgQhYAAAAABIiQBQAAAAABSnnIMrPfm9l2M1vaqK2PmT1tZm8l/+zdaNt3zGyNma0ys9mprg8AAAAAgpSOM1l/kHT6QW3flvSsu4+V9GzyucysRNIlkiYk9/mVmYXTUCMAAAAABCLlIcvdn5e0+6Dmj0i6O/n4bknnNWq/391r3P0dSWskTU11jQAAAAAQFHP31L+I2QhJj7n7xOTzve7eq9H2Pe7e28x+KekVd7832X6XpCfc/S9NHPPzkj4vSQMHDpx8//33t1pHeXm5CgsLA3hHaA/GIfMYg+zAOGQeY5AdGIfMYwzab9asWYvcfUqm60B2iGS6gINYE21NpkB3/62k30rSlClTfObMma0efN68eWpLP6QW45B5jEF2YBwyjzHIDoxD5jEGQLAytbrgNjMbLEnJP7cn2zdKGtao31GSNqe5NgAAAAA4YpkKWY9Iujz5+HJJf2/UfomZ5ZrZSEljJS3IQH0AAAAAcERSPl3QzP4saaakfma2UdIPJd0k6UEz+5yk9ZIukiR3X2ZmD0paLikm6UvuXp/qGgEAAAAgKCkPWe7+8WY2faiZ/j+W9OPUVQQAAAAAqZOp6YIAAAAA0CkRsgAAAAAgQIQsAAAAAAgQIQsAAAAAAkTIAgAAAIAAEbIAAAAAIECELAAAAAAIECELAAAAAAJEyAIAAACAAEVa62BmH5D0KUknSxosqUrSUkmPS7rX3UtTWiEAAAAAdCAtnskysyckXSHpKUmnKxGySiR9X1KepL+b2bmpLhIAAAAAOorWzmRd6u47D2orl/Ra8usWM+uXksoAAAAAoANq8UxW44BlZsPN7NTk43wz635wHwAAAADo6tq08IWZXSnpL5LuSDYdJWlOimoCAAAAgA6rrasLfknSiZL2SZK7vyVpQKqKAgAAAICOqq0hq8bdaxuemFlEkqemJAAAAADouNoasuab2Xcl5ZvZhyU9JOnR1JUFAGjM6+vb1AYAADKvrSHr25J2SFoi6SpJc5VYxh0AkGJeUymtXy6vKm/UVnVIGwAAyA6t3ow4KV/S7939Tkkys3CyrTJVhQEAEgHLVy+SP/0HafgEhc78vBQOy1cvPKDN8gszXSoAAEhq65msZ5UIVQ3yJT0TfDkAgIP5iw8nHqxbpvjc38rfnJ8IWMk2bXtXHmfqIAAA2aKtISvP3ffPSUk+LkhNSQCA/XLyFLrku1LDmap1y+QvPLR/s83+rDR4tCwUzlCBAADgYG0NWRVm9r6GJ2Y2WVJVakoCADSwUEjq2U+hT/w/yezAbR/8mGzs+2S5+c3sDQAAMqGt12RdI+khM9ucfD5Y0sUpqQgAcKC6avn6FZIfeOcMX7dMNuHEDBUFAACa06aQ5e6vmlmxpHGSTNJKd69LaWUAgAMXvmhgIcnj+6/RYuELAACyS1unC0rSCZKOkXS8pI+b2WWpKQkA0Nj+hS+UuAYr9NmfHHCNFgtfAACQXdoUsszsj5J+LukkJcLWCZKmpLAuAIDUaOGL7rLZn5WNeZ/UvY9Cl3xvfxsLXwAAkF3aek3WFEkl7gddEAAASCkLheQ9+yXOXkn7F7loqg0AAGSHtoaspZIGSdqSwloAAE2wUEg6KEg11QYAALJDW0NWP0nLzWyBpJqGRnc/NyVVAQAAAEAH1daQdW0qiwAAAACAzqKtS7jPT3UhAAAAANAZtHV1wfeb2atmVm5mtWZWb2b7Ul0cAAAAAHQ0bZ0u+EtJl0h6SImVBi+TNDZVRQEAADRWGatVbX2syW3RcEQFkWiaKwKA5rU1ZMnd15hZ2N3rJf2fmb2UwroAAAD2q62P6VsL5jS57eap5xGyAGSVtoasSjOLSnrDzH6qxFLu3VJXFrKZV5VL8bgUDsvy+GsAAAAANNbWkHWppLCkL0v6iqRhki5MVVHITl5VLt/0lnzRU1L5XqlHX9kJZ8gGjpDlF2a6PAAAACArtHV1wXXJh1WSrktdOchWXlWu+KO/kjaueq+xdId8w0r56OMVOu3TBC0AAABArYQsM3vQ3T9mZksk+cHb3f2YlFWGrOGxmPzNfx4YsBpb+7p8zTHSxJNk1qYFKwEAAIBOq7UzWf+V/PPsVBeCLFZbKX/9mRa7+KtPyEYfLxV0T1NRAAAAQHZqMWS5+xYzC0u6y91PTVNNyDbxuFRV3nKfvdvVxMlOAAACEQ1HdPPU85rdBgDZpNXvSu5eb2aVZtbT3UvTURSyjIUkM8lbCFGRHEmWtpIAAF1LQSTKMu0AOoy2/tdPtaQlZva0pIqGRne/OiVVIbuEw9KISdI7i5vtYsXTkkELAAAA6NraGrIeT36hC7K8bgrNuFjxDSukWN2hHaL5svefK4vmpb84AAAAIMu0dQn3u1NdCLJcj74KXfI9xZ+6S9qx4b32waMVmv1ZqVuvjJUGAAAAZJM2hSwzGyvpJ5JKJO0/XeHuo1JUF7KMRXKkAcMU+ujXpNoaqapMKugh5URl+awoCAAAADRo63TB/5P0Q0m/kDRL0mfEKgddkuV3l/K7Sz37ZboUAAAAICu19c6x+e7+rCRz93Xufq2kU1JXFgAAAAB0TG1eXdDMQpLeMrMvS9okaUDqygIAAACAjqnFM1lmNjD58BpJBZKuljRZ0qckXZ7SygAAAACgA2rtTNabZrZE0p8lrXb3jUpcjwUAAAAAaEJr12QNlfRzSSdLWm1mc8zsYjPLT31pAAAAANDxtBiy3L3e3Z9y989IGqbEKoPnSXrHzO5LQ30AAAAA0KG0dXVBuXutpOWSVkjap8Q9swAAAAAAjbQasszsaDP7hpm9JukxSWFJH3H341NeHQAAAAB0MC0ufGFmLylxXdZDkj7v7gvTUhUAAAAAdFCtrS74HUnPu7unoxgAAAAA6OhaDFnuPj9dhQAAAABAZ9DmhS8AAAAAAK0jZAEAAABAgA4rZJnZ+83sOTP7l5mdl6KaAAAAAKDDam11wUHuvrVR01clnSvJJL0kaU7qSgMAAACAjqe11QV/Y2aLJP3M3asl7ZX0CUlxJW5IDAAAAABopMXpgu5+nqQ3JD1mZpdKukaJgFUg6bzUlgYAAAAAHU+r12S5+6OSZkvqJemvkla5++3uviPFtQEAAABAh9NiyDKzc83sRUnPSVoq6RJJ55vZn81sdDoKBAAAAICOpLVrsn4k6QOS8iXNdfepkr5qZmMl/ViJ0AUAAAAASGotZJUqEaTyJW1vaHT3t0TAAgAAAIBDtHZN1vlKLHIRU2JVQQAAAABAC1o8k+XuOyX9T+M2M+vj7rtTWhUAAAAAdFCtLXzx/UaPS8xstaRFZvaumU1LeXUAAAAA0MG0Nl3wgkaPfybpv9x9pKSPSfpFyqoCAAAAgA6q1ftkNTLE3Z+QJHdfoMRiGAAAAACARlpbXXCUmT0iySQdZWYF7l6Z3JaT2tIAAAAAoONpLWR95KDnIUkys4GSfp2SigAAAACgA2ttdcH5zbRvk/S/KakIAAAAADqw1lYXDJvZVWZ2g5mdeNC27ze3HwAAAAB0Va0tfHGHpBmSdkm63cxubbTtgqZ3AQAAAICuq7WQNdXdP+Hut0maJqnQzP5qZrlKLIYBAAAAAGiktZAVbXjg7jF3/7ykNyQ9J6kwhXUBOExeU9mmNgAAAKRWayFroZmd3rjB3a+X9H+SRqSqqGwXi9errLZaWyv3aUtlqcrqqlVXX5/pstCFefle+dN3yyv3tdgGAACA1GttdcFPNdP+O0m/S0lFWa6irkYvb3tbT2xcrvK6GklSXjhHpwwp0qlDi9UtJzfDFaKr8fK9ij/0U2nPNvme7Qpd+BUpHj+kzQp6ZLpUAACALqG1+2Q1y8wGufvWIIvJdpWxWj21cYWe2rj8gPbq+jrN3bBMO6vLdcnoKQQtpJlL8Xji4Y71ij/4U6m+TirdmWiL10vumSsPAACgi2ltumBL7gqsig6itj6mf2xc0ez2BTvWqTxWk8aKAMkKeyt08beknv0TDbu3vBew+g5V6KNfk3XrmbkCAQAAupgjDlnuflaQhXQEb+7aKFfLZwRe2LJGzlkDpJkV9lborKsOaQ995EsELAAAgDQ74pBlZl1udcGyuupW+5TX1ShOyEKaeflexZ+485D2+KO/YeELAACANGvPdMHlrXfpXEZ079tqn5Hd+ykcas/HChyexgtfSJJ6DZC690k83rFe8Yd/QdACAABIoxYXvjCzrza3SV3wPlnDC/uqWyRXFc1cdxWxkI7vNyzNVQGeWNxC2n8NluJxxR+8KXFtVrxercxyBQAAQIBaO+Vyo6Tekrof9FXYhn07nfxIjv5zwgzlhMKHbAuZ6QslJysvfMQLNgJHxAp7K/Sxb0ujjt2/yIV1b9z2dVk3lm8HAABIl9YSwWuS5rj7ooM3mNkVqSkpe0VCYR1V2FvXTz5bz2xaqTd2bVTcXSW9B+mMYRPUI5qvKCELGWDdeyt0xpWy3PwW2wAAAJB6rSWCz0ja1cy2KQHX0iHkhMLqk9dN5404VqcPK5EkRUMR5UVyMlwZurqmwhQBCwAAIP1aDFnuvqqFbduCL6fjiIYjnLUCAAAAcIj2LOH+2yALAQAAAIDOoLXVBfs0t0nSmcGXAwAAAAAdW2vz3XZIWqdEqGrgyecDUlUUAAAAAHRUrYWstyV9yN3XH7zBzDakpiQAAAAA6LhaC1m3KXGfrENClqSfBl4NAldXH1NVfUxxjytsIUVCYeWzEiIAAACQMq2tLvi/LWz7n/a+uJm9K6lMUr2kmLtPSV4H9oCkEZLelfQxd9/T3tfqispqq/XkhuV6cdsaVdfHFLaQJvcbpvNHHKeeufkKW5e7nzQAAACQci3+lm1mJ7WyvYeZTWxnDbPc/Th3b7jv1rclPevuYyU9m3yOw1RWV61bljyjZzavVHV9TJJU73Et2LFOP3r9Ce2tqcxwhQAAAEDn1NqpjAvN7CUz+4GZnWVmU83sg2b2WTP7o6THJAV9t9OPSLo7+fhuSecFfPxOrz4e10tb39aWyn1Nbq+I1erBt19TVaw2zZUBAAAAnZ+5e8sdzHpL+qikEyUNllQlaYWkx939xXa9uNk7kvYosWLhHe7+WzPb6+69GvXZ4+69m9j385I+L0kDBw6cfP/997f6euXl5SosLGxPyR1C3F1bqkpVH4+32O+obr0VMmuxTyp0lXHIZoxBdmAcMo8xyA6MQ+YxBu03a9asRY1mZqGLazVkpfTFzYa4+2YzGyDpaUn/KemRtoSsxqZMmeILFy5s9fXmzZunmTNntq/oDqCsrlrfeOVvcrU8tj854SPqk9ctTVW9p6uMQzZjDLID45B5jEF2YBwyjzFoPzMjZGG/jK584O6bk39ul/Q3SVMlbTOzwZKU/HN75irsoFzq20p4alhpEAAAAECwMhayzKybmXVveCzpNElLJT0i6fJkt8sl/T0zFXZc3XKi+tCQcS32mdxvmHIIWQAAAEDgWrtPVioNlPQ3S1wTFJH0J3d/0sxelfSgmX1OiftzXZTBGjukkIU0dcAI/Xv7u3q3fNch23tF83XByOO5XxYAAACQAq2GLDPrIam/u689qP0Yd198pC/s7m9LOraJ9l2SPnSkx0VCYU6u/nPiDP1r61o9s2mV9tVVKy8c0UkDx+i0YePVIycv0yUCAAAAnVKLIcvMPibpNknbzSxH0qfd/dXk5j9Iel9Kq0O7FObk6dSh4zV94Gi5JJOUF44oJ5zJE5gAAABA59bab9vflTTZ3beY2VRJfzSz77r7X5X4nR1ZLhwKqXuUs1adndfVSrVVUiRHllvQbBsAAABSr7WFL8LuvkWS3H2BpFmSvmdmV0utrA8OIC28rlba9q7iv/umfMW/5TWVibbt65Jtr8hrKjNdJgAAQJfR2pmsMjMb3XA9VvKM1kxJcyRNSG1pANokVqv4w7dI9TH5c/dKHpf1G6r4X3+RbLtP1n+YfMhomWX0rg0AAABdQmsh6z900LRAdy8zs9MlfSxlVQFou1BENvuz8rl3SnL5P/90wGlme99pUt8hBCwAAIA0ae23rgolllo/2PslvRJ8OQAOl+XmyUYeKzvjikO3HX+q7P1ny1q5OTUAAACC01rIuk1SWRPtVcltALJBKCQr7H1oe/c+krFGDQAAQDq1Nl1wRFP3wnL3hWY2IjUldX6VsVrV1sea3BYNR1QQiaa5InRk+xe++Outh257/kEpEpXGT2OFQQAAgDRpLWS1tPZ3fpCFdCW19TF9a8GcJrfdPPU8QhYOT6OFL6TEFEH1P0r+j7slufy5e2X9j2LhCwAAgDRp7TeuV83syoMbzexzkhalpiQAhyW58IVksvedJvvAubKiE2RnXrm/jYUvAAAA0qe1M1nXSPqbmX1S74WqKZKiks5PYV0A2shy86SRx8o+9f+kHv3eW+SiqTYAAACkXIshy923SZpuZrMkTUw2P+7uz6W8MgBtZrl58v5HyxotctFUGwAAAFKvxZBlZnmSviBpjKQlku5y96ZXbACQUU2FKQIWAABA+rV2kcbdSkwPXCLpDEk/T3lFAAAAANCBtXZNVom7T5IkM7tL0oLUl9T5RcMR3Tz1vGa3AQAAAOi4WvuNvq7hgbvHmHoUjIJIlGXaAQAAgE6qtZB1rJntSz42SfnJ5ybJ3b1HSqsDAAAAgA6mtdUFw+kqBAAAAAA6A+5OCgAAAAABImQBAAAAQIAIWQAAAAAQIEIWAAAAAASIkAUAAAAAASJkAQAAAECACFkAAAAAECBCFgAAAAAEiJAFAAAAAAEiZAEAAABAgAhZ6FQ8Hm9TGwAAAJAqhCx0Gl5dIW1YIa+ubLENAAAASCVCFjoFr66Qv/hXxR++Vb70eXlNZbLt4ffaCFoAAABIg0imCwDay92lfbvlS+Ynnj//kBSrk/btki99IdH20hxZ8fslFWSwUgAAAHQFnMlCh2dmUq/+snO+KJlJSoSqhoClSI5CF31TyuuWwSoBAADQVRCy0ClYNE92dIls9ucO2Ra68GtS/2GySE4GKgMAAEBXQ8hC5xGvlzatPqTZN72VmD4IAAAApAEhC51CwyIXvuT5Q7e9+DALXwAAACBtCFno8BILX+x6L2BFchT6+Pdk53zpgGu0FKvNXJEAAADoMghZ6PASC18MSISqnNzEIhf9h8mGlyQWw2hoY+ELAAAApAFLuKNTsGiedPR42ed/LoVzkotc5EhHlxzUBgAAAKQWIQudhkXz2tQGAAAApBLTBQEAAAAgQIQsAAAAAAgQIQsAAAAAAkTIAgAAAIAAEbIAAAAAIECELAAAAAAIECGri4u7qzpWp9r6WKZLAQAAADoF7pPVRcXi9aqM1enNXRu1dM9mRSykEweN1rDC3uqew72lAAAAgCNFyOqCYvF6barYq1uXPKvqRmewFu5cr8EFPfSVSR9Sz2h+BisEAAAAOi6mC3ZBlbG6QwJWgy2V+3THihdVXleTgcoAAACAjo+Q1cXEPa7Xdq5vMmA1WLtvhypjhCwAAADgSBCyupja+not3bO51X5r9+1MQzUAAABA50PI6oLC1vqwR9rQBwAAAMCh+E26i8mL5Gj6wFEt9jGZxvYcmKaKAAAAgM6FkNUFje7RT/3zCpvdPn3gSEXD4TRWBAAAAHQehKwuqDAnT1875lQNL+xzQLvJNH3AKF0w8ngVRKIZqg4AAADo2LhPVhfVO7dAV0+YqbJYjdaU7lBOKKziXgMVDUcIWAAAAEA7ELK6sMJongqjeRpc0DPTpQAAAACdBtMFAXRJXl0hryprtQ0AAOBwEbIAdDleXSF/aY7iD/5UXlHabBsAAMCRIGQB6FK8tlq+eJ78jeekXZsV/8vP9wes/W1/u01eW53pUgEAQAfFNVkAuhSL5kkTTpIve0naszURqn79X5J7okM4otApn5JC/B8UAAA4MvwWAaDLsW49FfrYN6XeyZtuNw5YF31TGjBMxiqbAADgCBGyAHRN4RxpwNEHtuV1k3r1z2jA8vr6NrUBAIDsRcgC0OUkrsH6m7Tq1QM3VJQq/tDPMrbwhVdXyNe+Lq+qOLStuqKFPQEAQDYhZAHoUg5Y+EKSwhHZ7M9KvQclnmdo4QuvrpC/+oT8sV/L5z8gr6o4sG3eAwQtAAA6CBa+ANClHLDwxb6diWuwBh4tGzFR8Qd/mmjLxMIXdbX7g58v/5fk9VJ+d/lrTyfaVr4iO2F2YkojAADIaoQsAF3O/oUvynZL/YbKwjnSwW3pvi4rv7tCF39b8Qdukupq5CteeW9bKKzQhV+VevRPb00AAOCIMF0QQJdk3XpK/Y8+IEw11Za2eiIRqc8QhT769UO3nf0FadAoWQ4rHgIA0BEQsgB0WRYOt6ktbWI18jWvHdq+5nWpvi799QAAgCNCyAKALLB/kYtXn0i2mGSW2Lb8JRa+AACgAyFkAUA2aLTwhUJhhT76NYU++QMpJ1dSYuELVezNXH0AAKDNCFkAsppXV8gr9iW+6mOZLid18rsrdPG3pNwChS74ijR4dOIarYY2Fr4AAKDDYHVBAFnJqyuknZsUf+VRafs6KZonK5kuHTtLKugus871f0QWicj7DFXoyp9LZvsXuWiqDQAAZDdCFoCs49UV8oVPyRc8/l5jdYX8lUflb85T6BPfk3p2vrM6Fono4G/LTbUBAIDs1rn+KxhA51C2+8CA1VhVmeJz75RXlae3JgAAgDYiZAHIKl5bLV8wt+VOW9ZKtdXpKQgAAOAwEbIAZJdYrXzX5la7+b6daSgGAADg8BGyAGQXC0n53Vrvltd6HwAAgEwgZAHIKpZfKDv+1JY7FfaWuvVKSz0AAACHi5AFIOvY0CJp4IjmtspOvVTKLUhnSQAAAG1GyAKQdSy/UKHzr5EdM0OK5Ly3oc9ghS78qmxokSwczlyBAAAALeDmKwCykhV0lz54sWz6eVJNVSJshXOk/EKZWabLAwAAaBYhC0DWsmiuFM2VCnpkuhQAAIA2Y7ogAAAAAASIkAUAAAAAASJkAQAAAECACFkAAAAAECBCFgAAAAAEiJAFAAAAAAEiZAEAAABAgLhPFpBCHq+Xqiul6gqprloq7CVForLcgkyXBgAAgBQhZAEp4jWV8o2r5f/8k7RvV6LRTBp1nEIfulRW2DOzBQIAACAlmC4IpIDH6+UbVsn//j/vBSxJcpfWvq74gzfLq8oyVyAAAABShpAFpEJ1ReIMVnP2bpOvfVPu8fTVBAAAgLQgZAGpUFUule1usYsvnidVVaSnHgAAAKRN1oYsMzvdzFaZ2Roz+3am6wEOS21V632qKyV5yksBAABAemVlyDKzsKT/lXSGpBJJHzezksxWBRyG7n0kWct9BgyTwjlpKQcAAADpk5UhS9JUSWvc/W13r5V0v6SPZLgmoO0iUWnEhBa7hKaeJcvNT1NBAAAASBdzz77pSmb2UUmnu/sVyeeXSprm7l9u1Ofzkj4vSQMHDpx8//33t3rc8vJyFRYWpqZotFmXGYd4vbR7q1QfO3RbYS+poLtkmfl/ji4zBlmOccg8xiA7MA6Zxxi036xZsxa5+5RM14HskK33yWpqntUBadDdfyvpt5I0ZcoUnzlzZqsHnTdvntrSD6nVlcbBK8vkby2UL34+cZ1W/2EKTTtb6tlflpe5GxJ3pTHIZoxD5jEG2YFxyDzGAAhWtoasjZKGNXp+lKTNGaoFOGJW0F2aNEM2dkriHlmRHKYIAgAAdHLZGrJelTTWzEZK2iTpEkmfyGxJwJGxUCgxNRAAAABdQlaGLHePmdmXJT0lKSzp9+6+LMNlAQAAAECrsjJkSZK7z5U0N9N1AAAAAMDhyNYl3AEAAACgQyJkAQAAAECACFkAAAAAECBCFgAAAAAEiJAFdABeXy+vq5V7vMU2AAAAZB4hC8hyXl8v7d2u+L3XSaW75B5v1Hbt/jYA6Iq8qlxeVd5qGwCkEyELyHZluxS//0Zpz9bEn6W7EgHr/hulPdsSf5aXZrpKAEg7rypX/Ln7FH/yrv2hqqk2AEg3QhaQ7XILpAFHJx5X7lP8vusU/9MNUk1lom3waCmSk7n6ACADvLpS/q+/SasWSO8sToSqulrFn7vvvbbn7pNXV2a6VABdECELyHKWX6jQ2f8hDStONNRUSXU1icejj1fotE/L8gszV2AGeU2lvD7WahuATiiaK5t8mpTXLfH8ncWK/+o/EwFLkqJ5Ck09S8qJZq5GAF0WIQvoCHLzFZp+3iHNoRPPl/IK9j+vqY+prLZau6rLVVpbpfK66jQWmV5eVa74U3+Qdm/dH6qaagPQOVkoLPXsp9DHv/de0Gr4d5+Tp9DF35H6DJKFI5krEkCXxXceIMvtX+Rizu2HbIv/5ecKXfJdec++Kqur1Zx339CCHetUF6+XJI3q3k8fHzNFg/J7KNqJftHwqnLFH79DWr9c8fXLFLrkO/LC3oe29eYXLKAzs1BYXtBd6jtE2vTWext6DZC69+HfP4CM4UwWkO0aFr5ouAZrxERp6NjE48p9it9/o7x8r+5bs0D/2vb2/oAlSW+X7dRNb/xDmys72cIYZlJ+8n+ua6sVv/8nij/4U2n98kRbOEeKRCXjWxzQmXlVueLP/PHAgCVJO9Yr/sSdLHwBIGP4DQTIdo0Xvhh9vEJnXKnQuV9+7xqtwaNVJWnV3m1N7l7vcd371gKVdaKpg5bXTaEPXSqNOyHRUFst7dyYeJzfXaGPf1fq0U8W4lsc0FkdsPCFJEXzZB+69MBrtFj4AkCGcB4dyHINC1/4m/+UHTtr/yIXDW1+zEz97+qXVVVf1+wxNlTsUV19vdSJFiG0vG4KnfYZxd96TWp09s5OujAxTYiABXRu0VzZ5A/LV78qxevfuwbr6GLF/3xjoo2FLwBkCCEL6AAsv1CafJosJ/eQtnKPa0d161NiWgphHdH+67IaBSxJ8vn3ywaP5HosoJOzUFjes3/izHUstn+Ri6baACDd+K9eIE0qY7XaW1PZ5FdlrLbV/RsHrMZtYQtpYEGPlveVqbCJ/TuqxgtfSJLyu0tDixKPk9doaQ8rDAKdXWKFwf4HhKmm2gAg3fjuA6RJbX1M31owp8ltN089TwWRI5vS0i0nV2cMm6C3Src322dSn8HKCYWP6PhZyUxquDdYwzVY+YWKP3OPtOpVFr4AuhBr4ntbU20AkE6ELKATGFHYR9P6j9C/d7x7yLZe0Xx9YszUIw5x2Six8MWn5HndEjcjTS5yEfrQpfK8wgPaAAAA0o2QBXQC3XJydfHoyXpfv2F6YsNybaksVUEkqpMGjdYHB49V9040VbCB5XWTTrwgsaJYMkw11QYAAJBuhCygk+iWk6vj+g3TmB79FU+2FURyFOnE02Ysr6BNbQAAAOlEyAI6mcJoXqZLAAAA6NKYTwMAAAAAAeJMFtBOdfWx/fegygmFld/MAhPRcEQ3Tz2v2W0AAADoHPjNDjhCsXi9yupq9PTGFVq0c73q3VXUc4DOOnqi+uZ2U14k54D+BZFop1rhDwAAAE0jZAFHoN7j2lK5Tz9782nVxN+74e2inev12s71urzo/Tq+77BDghYAAAA6P67JAo5AZV2t/nf5/AMCVgOXdPfqf++fQggAAICuhZAFHIHt1WXaU1PZ7HaX619b1yoejzfbBwAAAJ0TIQs4Ahsr9rbaZ335btXG61NfDAAAALIKIQs4Aj1yWr8XVWFOnsIh/okBAAB0NSx8ARyBMT36KycUVl0LZ6pmDSlSTiicxqqQ7TxWK9VWJ560IagDAICOiZAFHIFoOKzzhh+rh955rcntx/QZqt65BWmuCtnKY3VSVZl80T/kb78pSbKRx0g+QB6rlbG0PwAAnQohCzgCueEcTR80St1yovr7u4u1pzaxCEZeOEczB4/Vh48ar8Kc3AxXiWzgsTpp2zrFH/65FHtvxUl//RnpqBnS1nflg0bKWO4fAIBOg5AFHKGCSFTT+o/QxN5DVBOPqd7jKohElRfKUU6YaYJIqq1WfM5tBwSs/dwVn/PfCn32JomQBQBAp0HIAtohFAqpezRP3TNdCLKWb1ot1VQ136G2Wr5+uax4WvqKAgAAKcXSZwCQSpveCqYPAADoMAhZAJBKBT2C6QMAADoMQhYApFBbpgFayfQ0VAIAANKFkAUAqRTNk02a0exmKzlRys1PY0EAACDVWPgCAFLI8rpJJ10gFfaSL/qHVJtcBCOaJxX2lM04NdEHAAB0GoQsAGgjr6+XHbQ8f1NtB7P8QumE02XHzpQqShON3XpKCxYltgEAgE6FkAUAbeCxOmnXJnlhb1m3ns22NcciUSkSPWiRC0thxQAAIFO4JgsAWuGxOmn7esXv/4nif7lFXlHaZBsAAIBEyAKA1nlc8Wfukepj0q5NiVD11iLFH/rp/jZf9i95bXWmKwUAAFmAkAUArbCcXIUu+IrUa2CiYdcm+RN3JgKWJDvuFNkxM2TRvAxWCQAAsgUhCwDawAp7KXTxt6X87gduGDdVNv18VggEAAD7EbIAoA08VieV7nhvCfYGOzdK9XWZKQoAAGQlQhYAtGL/IhcN12BJ2r8y4K7NLHwBAAAOQMgCgNY0XvhCiWuwQv9xm9S70TVaLHwBAACSCFkA0IrGC1/YcafIpp8nyy9U6KJvSr0Hyo6dxcIXAABgP25GDABtYIW9FPr4dyQL7V/kwgp7KXTJgW0AAACELABoIzt4ZcFm2gAAQNdGyAKALFEZq1Xt/oU1DhQNR1QQiaa5IgAAcCQIWQCQJWrrY/rWgjlNbrt56nmELAAAOggWvgAAAACAABGyAAAAACBAhCwAAAAACBAhCwAAAAACRMgCAAAAgACxuiAAZIloOKKbp57X7DYAANAx8FMbALJEQSTKMu0AAHQChCwASIOmbjRc73HtrankRsMAAHQyhCwASIOmbjR8YkWeHlwwhxsNAwDQybDwBQAAAAAEiJAFAAAAAAEiZAEAAABAgAhZAAAAABAgQhYAAAAABIiQBQAAAAABYgl3AEiDaDiim6eed0Dbopde0c1T369omG/FAAB0JvxkB4A0KIhED7kXVthC6pVbkKGKAABAqjBdEAAAAAACRMgCAAAAgAARsgAAAAAgQIQsAAAAAAgQIQsAgAB5TVWb2gAAnRchCwCAgPi+3fLXnpZXV7TYBgDo3FjCHQCAAPi+3Yo/8BOpbLdUXSFNP0+qqTqw7QPnyvK6ZbpUAECKEbIAAGgnj9VJFXukyn2J568/Iy/bLW19Ryrfk2jbslbmnskyAQBpwnRBAEAgPFYrryxLfMXrM11OWlkkR+o/TKGPfl0KJ///cs1r+wOWBo1U6Lz/kuUXZq5IAEDacCYLAFKsMlar2vrYIe31HldlrFYFkWgGqgqO19VK1eXy15+Vr1suhUKy4mlS8ful/EJZqGv8f55FovKBI2QzLpY/d997G8IRhT76dVk0L3PFAQDSipAFAClWWx/TtxbMOaT9xIo81dbHOnTI8littHmN4nP+W2oUJH3bu/J/Py675Duq7dFHuVn8Ht3jUlVFIhwmr5dqqq1NKsvkrz5xYFt9TP6vv3E9FgB0IV3jvxcBAKlRU6X4328/IGDtV10u/9ttKivbpSW7N6m8rib99bVF6U7F7/5/8oVPyasrEgHroLa2OGDhC0nq0Xf/1EF//Rn5y4+wwqCU+IxrKg9qjB/aBgAdGCELAHBE3F2+eqEUq2u+U+kORcv26MG3X9Mti5/Rvtosu1+Uu+IP3ixVlckXPC5f+JS0b5fif77xvbY1r6u+lYB48MIXGjRSoU/98IBrtHzLWqmLL3zhVeWKP/V/8qUv7g9VXlUu7dt1QBsAdHSELADAkYnVyTeubrVbZNs6Dcjrrs2Vpfr9qpdVkWVntEKzPyeFwpIkX/C44nd9W6oqSzwfVqzdQ8fo8U2rVFpbJW8mJB2w8MXQosQiF3ndpIHDD2zrwgtfeHWl4k/fLa19XT7/gUSoqqtR/B9/kKorE20rXuHGzQA6BUIWAODIhEJtWswhHs1TXXK1wRV7t6qmqamFmWImDR6t0AXXJB434keN055TL9VNq17S4xuW6sbXn9TeFs7EWSSaCFXnXy0r6N5sW5cVCsmK3y8p8Tn7/AcU/7/vSWtfT2zPLZANGy9l8fV7ANBWhCwAaKSmvk7V9XWKezzTpWQ9C0dkx85spVNI8aPH6+2ynfubNlbsTWldhy0nR+o5QArnHNAc6zNIu2urVBFLnHnbW1ulv77zuqpamB5pkagsmt9qW1dk0TzZiAmys/9DDUFr/xL3oZBCl3xX6jVAFg5nrEYACAqrCwKApPK6ar1btlv/2rZWsXhcE3oP1vv6DVNBJKpIqH2/9EXDEd089bxD2he99Iqi4Q7+bbhnf2lokbSp6WmDsWNn6fXS7fvPZElSJNuWdC/dpfj9N0qx2gOacxbP1/Bonq4umqrb1yxQ3F2Ldm7QRaMmK185zRwMLbFonjRigrzXAGnvtvc25OZL3XsTsAB0Gln2kw4A0q+0tko/ffMZ/c+yeXpt5wYt3r1Jf167UN9/9VFtKN+jWDtvrFsQiapXbsEhX2ELdejl2yXJ8gsVOveLsvEf2H9dkyQpJ0+xaWdp66ST9ODG5e81h8I6qlvvDFTajIaFL5ILVviwYlWf9pn97yVn4VMatX2jzj2qRFLi3mbt/fvQlXlVueJP/O7AgCVJVRXyJc+z8AWAToOQBaBLK6+r0R0rXtS2qn2HbKuJx/SLJc+p4qAzHDiQ5XeXnfJJha78qexj31TlhV/Vvkt/oKcGHa2frnpJsUZTLz88tFi5WXa2IjT7s4lQdfR4lX74Mr3QrVDV53xRCoXlw4q1b9g4Pbt1baKvLPvOxHUQjRe+kJS4BmvyaWp8jRYLXwDoLDr4PBUAaJ/KWK3W7tvR7PaaeEyLdqzXzCFjFTJ+uW6O5eZLufmybr0Ur6nSHcuf1zvlu/Zvj1hIpw4t1qlDi5UbzqKpdg0LX1zyHalHPy3es1l/27hc1YPGavZFX1d5QQ/dvOolldVVS5KO6TtUOe2cPtplJRe+8DWvS7n5iWuwevSRBo+RXl+aCF1Hs/AFgM6BkAWgS3t7385W+yzbs0UfGDhK+RFCVlv0zM3XlyfOVHldjdaX71Y0HNGo7n0VDUeUl00BK8miufIBw2WhkKbkDNf8rWv0xNa3tL3vUVq9Yen+gNUtEtXHRr1P+YSAI9JwPZbO/aKs9+D3FrkYMUF6Z5tCH/6u1JOFLwB0DoQsAF1aW85KREIh1Xtce5u5XiQajnT4a6uCVpiTq8KcXA0q6JHpUtrEklMAu+fk6WuTPqQnNyzXi9vWqLo+prCFNLnfMJ0/4jj1zGWVwPawaJ50dIkUztkfpiyaJ0XzCVgAOhVCFoAubUzP/grJFFfTN5mVpOkDR8nd9a0Fc5rcfvPU8whZnUj3aJ7OG3GMZg8rUdzjCllIOaGw8iPZdxauI2ry3mpmBCwAnQpzXwB0adFQRCcNHtPs9gF53TW6Rz/Vc9+sTsWrKw5ZYKFxW044oh7RPPXKLVCPaB4BCwBwWAhZALq0/EiOzht+jE4eNEahhhukJo3o3ldfO+ZDKsxp4n/e0WF5dYV84VPyd5dI7oe0sbodAKC9mC4IoMvrlpOrC0cep7OPnqS3SrepLl6vUT36Ja8rImB1Jl5dIX/1SfmrcyWZNOECeW31gW1nXSWNmJhYMREAgCOQkZBlZtdKulJSw7rJ33X3uclt35H0OUn1kq5296cyUSOAriU/ElV+RDphwIhMl4JUMpMKeyWfuLR3h+IP3yJteTvRFArJ8gsl7oUFAGiHTJ7J+oW7/7xxg5mVSLpE0gRJQyQ9Y2ZF7l6fiQIBAJ2L5RZIJR+QJPk//5Ro3B+wwgpd8BVp8ChZTm6GKgQAdAbZNl3wI5Lud/caSe+Y2RpJUyW9nNmyAHR10XBEN089r9lt6Dgst0CaeJJ8yfMHthedIA0aScACALRbJudDfNnMFpvZ782sd7JtqKQNjfpsTLYBQEYVRKLqlVvQ5BfLt3csXl0hf+UxaefGA9tX/lv+zmIWvgAAtJu5N39vmHYd2OwZSYOa2PQ9Sa9I2inJJd0gabC7f9bM/lfSy+5+b/IYd0ma6+4PN3H8z0v6vCQNHDhw8v33399qTeXl5SosLDzCd4SgMA6ZxxhkB8YhAzwuVeyTKkolSeXRQhWG6qXqRsGqV//EzXHNmjkIgsa/hcxjDNpv1qxZi9x9SqbrQHZI2RwXdz+1Lf3M7E5JjyWfbpQ0rNHmoyRtbub4v5X0W0maMmWKz5w5s9XXmjdvntrSD6nFOGQeY5AdGIf085pK+fKX5f98VAqF9ULJBZoxc6Z82b8S12iFwgp9gOuy0o1/C5nHGADBytTqgoPdfUvy6fmSliYfPyLpT2Z2qxILX4yVtCADJQIAOqH9C19YSNZnkLR2qyyad2AbAQsA0E6Zulr7p2Z2nBLTBd+VdJUkufsyM3tQ0nJJMUlfYmVBAECQLLdAGj9NCoWlt7cd0kbAAgC0V0ZClrtf2sK2H0v6cRrLAQB0MZZb0KY2AACOBHdbBAAAAIAAEbIAAAAAIECELAAAAAAIUKYWvgAAoMPzumqptlb+zmJp7zap92DZyIlSTi4LaABAF0bIAgDgCHhNpXz1Qvlz90n1sffaIzmyD39aGnWsLDc/cwUCADKG6YIAAByJre/Kn777gIAlSYrVyZ+4U9q1KTN1AQAyjpAFAMBh8soyxV/8S4t94i8+LK8qT1NFAIBsQsgCAOBweVzatq7lPhtXJ/oBALocQhYAAIfLPdMVAACyGCELAIDDFQpL/Y5quc/gUZLxYxYAuiJWFwSANqiM1ar24AUOkqLhiAoi0TRX1LlU1NWo3uMyMxVGcmVmmS6pRdXRXIXff47Cj/262T6h6efL8gvTWBUAIFsQsgCgDWrrY/rWgjlNbrt56nmErCNUUVejd8t368kNy7Sjqlw9onn60NBxmtB7iAqz9D5T5XU1umPFC/pwv+Ea+4GPKPrKIwdOH7SQbOYl0qARmSoRAJBhhCwAQEZU1NXoT2te1cKd6/e37amt1O9XvazhhX30nxNmqns0L4MVNu2dsp1aXbpda0p36GNHlWjy5T9SdPWriu7bqdqe/RUbd4Kiud0UzS3IdKkAgAxhsjgAIO3cXcv2bDkgYDW2rny3ntywXHXNTNHMlIq6Gj2zaaUkKS7X/RuX6fsrntcfe/XRgyMn6O4ePfXtZfO0K16f4UoBAJnEmSwAQNqV19XoyQ3LW+zz4rY1mj2sRDnh7PlRFXfXvtrqA9pq4jEt3HHgcu5lddUarJ7pLA0AkEU4kwUAyIjNlXtb3F5dH1N9lt1nKhIKaUB+91b79c7tloZqAADZipAFAMiI/DYsFhLKslUG8yNRnX5USYt9hhf2UX4kJ00VAQCyUfbMwQCALBYNR3Tz1POa3YbDkxfO0fQBo/TM5pXN9inqOUCRLLzP1ID87jp50Bi9sHXNIdvywzn67LjpWbsyIgAgPfjNAADaoCASZZn2AOWEwzpt2Hj9e8c7KqurOWR7xEK6ZPQUdcvCsNItJ1fnjzhWE3oP1twNS7WhfK9ywxG9f8AInT5sgnpkYc0AgPQiZAEADlsQN2funpOrbx83W/es/rdWlW7b3z6sW29dOnaa+udl7418u+Xk6vh+wzSmR//9bbnhCGc1AQCSCFkAgCMQxM2ZQxZSv7xCXTX+JNXF61VWV6OCSFTRcFjdc7Lv/lhNycb7eAEAMo+QBQDIqIYpgb24eS8AoJPIviuKAQAAAKADI2QBAAAAQIAIWQAAAAAQIEIWAAAAAASIhS8AAIeNmzMDANA8fhICAA4bN2cGAKB5TBcEAAAAgAARsgAAAAAgQIQsAAAAAAgQIQsAAAAAAkTIAgAAAIAAEbIAAAAAIECELAAAAAAIECELAAAAAAJEyAIAAACAABGyAAAAACBAkUwXAADInMpYrWrrY01ui4YjKohE01wRAAAdHyELALqw2vqYvrVgTpPbbp56HiELAIAjQMgCAHQ4cXdV1NWotLZKu2oq1DOar7653VQQiSocYiY8ACCzCFkAgA6lLh7T5opS3bnyX9pRXb6/vVc0X5cXvV+ju/dTbiQngxUCALo6/rsPANCh7K2p0s8WP3NAwJKkvbVV+p+l87Slal+GKgMAIIGQBQDoMKpidfr7u4tVF69vcntcrofefk3ldTVprgwAgPcQsgAAHUYsXq/Xdm1osc+afTtU7/E0VQQAwKG4JgsAurBoOKKbp57X7LZs1JYA5e5pqAQAgKZl509QAEBaFESiHWuZdpOO6tZLGyv2NtulZzRfIWOiBgAgc/gpBADoMLrn5Omsoye12OfDQ8erW0cKjgCAToeQBQDoUMb1HKhZQ4qa3Pa+vsM0feBI7pUFAMgopgsCADqUbjlRnXP0JJ08aIz+sWGFdtSUqVc0Xx8eOl7987urW05upksEAHRxhCwAQIfTLSdX3XJy9YmxU1QXjysSCikvzA2IAQDZgZAFAOiwcsM5yg1nugoAAA7EpHUAAAAACBAhCwAAAAACRMgCAAAAgAARsgAAAAAgQIQsAAAAAAgQIQsAAAAAAkTIAgAAAIAAEbIAAAAAIECELAAAAAAIECELAAAAAAJEyAIAAACAABGyAAAAACBAhCwAAAAACBAhCwAAAAACRMgCAAAAgAARsgAAAAAgQIQsAAAAAAgQIQsAAAAAAkTIAgAAAIAAEbIAAAAAIECELAAAAAAIECELAAAAAAJEyAIAAACAAJm7Z7qGdjOzHZLWtaFrP0k7U1wOWsc4ZB5jkB0Yh8xjDLID45B5jEH7DXf3/pkuAtmhU4SstjKzhe4+JdN1dHWMQ+YxBtmBccg8xiA7MA6ZxxgAwWK6IAAAAAAEiJAFAAAAAAHqaiHrt5kuAJIYh2zAGGQHxiHzGIPswDhkHmMABKhLXZMFAAAAAKnW1c5kAQAAAEBKEbIAAAAAIECdMmSZ2UVmtszM4mY2pVH7CDOrMrM3kl+/abRtspktMbM1Zna7mVlmqu88mhuH5LbvJD/rVWY2u1E745BCZnatmW1q9G/gzEbbmhwTBM/MTk9+zmvM7NuZrqcrMbN3k99j3jCzhcm2Pmb2tJm9lfyzd6br7EzM7Pdmtt3MljZqa/Yz53tRajQzDvxMAFKkU4YsSUslXSDp+Sa2rXX345JfX2jU/mtJn5c0Nvl1eurL7PSaHAczK5F0iaQJSnzOvzKzcHIz45B6v2j0b2Cu1OqYIEDJz/V/JZ0hqUTSx5OfP9JnVvLvf8N//nxb0rPuPlbSs8nnCM4fdOj38iY/c74XpdQf1PTPVH4mACnQKUOWu69w91Vt7W9mgyX1cPeXPbESyD2SzktVfV1FC+PwEUn3u3uNu78jaY2kqYxDRjU5JhmuqbOaKmmNu7/t7rWS7lfi80fmfETS3cnHd4vvO4Fy9+cl7T6oubnPnO9FKdLMODSHcQDaqVOGrFaMNLPXzWy+mZ2cbBsqaWOjPhuTbUiNoZI2NHre8HkzDunxZTNbnJw60jBFp7kxQfD4rDPLJf3DzBaZ2eeTbQPdfYskJf8ckLHquo7mPnP+faQfPxOAFIhkuoAjZWbPSBrUxKbvufvfm9lti6Sj3X2XmU2WNMfMJkhq6rof1rZvgyMch+Y+b8YhAC2NiRLTMW9Q4nO9QdItkj4rPvt04rPOrBPdfbOZDZD0tJmtzHRBOAD/PtKLnwlAinTYkOXupx7BPjWSapKPF5nZWklFSvwPzVGNuh4laXMQdXZ2RzIOSnzewxo9b/i8GYcAtHVMzOxOSY8lnzY3Jggen3UGufvm5J/bzexvSkyB2mZmg919S3La8vaMFtk1NPeZ8+8jjdx9W8NjfiYAwepS0wXNrH/DhZtmNkqJhRXeTk5VKDOz9ydXs7tMUnNnYdB+j0i6xMxyzWykEuOwgHFIveQvMw3OV2JxEqmZMUl3fV3Eq5LGmtlIM4sqcXH5IxmuqUsws25m1r3hsaTTlPg38Iiky5PdLhffd9Khuc+c70VpxM8EIHU67JmslpjZ+ZL+R1J/SY+b2RvuPlvSByVdb2YxSfWSvuDuDReB/ocSK+/kS3oi+YV2aG4c3H2ZmT0oabmkmKQvuXt9cjfGIbV+ambHKTHt411JV0lSK2OCALl7zMy+LOkpSWFJv3f3ZRkuq6sYKOlvyTtDRCT9yd2fNLNXJT1oZp+TtF7SRRmssdMxsz9Lmimpn5ltlPRDSTepic+c70Wp08w4zORnApAalljEDQAAAAAQhC41XRAAAAAAUo2QBQAAAAABImQBAAAAQIAIWQAAAAAQIEIWAAAAAASIkAUAbWBm9Wb2hpktNbOHzKwg2T7IzO43s7VmttzM5ppZUXLbk2a218wea+XYt5nZB5OPv2xma8zMzaxfoz4zzaw0WcMbZvaDVo75P2ZW3tr+yfsHvph8X+c16v93MxvSxHFnmtnLB7VFzGzbQffcaam2IWb2l7b0ba/k+3syHa8FAEADQhYAtE2Vux/n7hMl1Ur6QvKm2X+TNM/dR7t7iaTvKnE/Jkn6maRLWzqomfWR9H53fz7Z9C9Jp0pa10T3F5I1HOfu17dwzCmSerVx/49LulvSByR9I7n/OZJec/fNTRzjeUlHmdmIRm2nSlqavKF4i8ws4u6b3f2jrfUNgrvvkLTFzE5Mx+sBACARsgDgSLwgaYykWZLq3P03DRvc/Q13fyH5+FlJZa0c66OS9p9pcffX3f3dIy3MzMJKhLtvtnGXOiVu/p0rKW5mEUnXJI9xCHePS3pI0sWNmi+R9Gczm2pmL5nZ68k/xyVr+nTy7N+jkv5hZiPMbGly2wgze8HMXkt+TU+2zzSzeWb2FzNbaWb3JUOtzOyE5PHfNLMFZtbdzMJm9jMze9XMFpvZVY3qmyPpk238PAAAaDdCFgAchmQIOUPSEkkTJS1q5yFPPIxjfCAZLJ4wswnN9PmypEeaOavU1P5/kjRbiaB3raQvSrrH3StbqOPPSgQrmVmupDMlPSxppaQPuvvxkn4g6cbGry3pcnc/5aBjbZf0YXd/nxLB7fZG245XIvCVSBol6UQzi0p6QNJ/ufuxSpxFq5L0OUml7n6CpBMkXWlmI5PHWSjp5BbeDwAAgYpkugAA6CDyzeyN5OMXJN0l6QsBHHewpB1t6PeapOHuXm5mZypxdmZs4w7Ja6gukjSzrfu7e6mks5L795b0LUkXmNmdknpLusXdD7gGy91fNbPC5Jmq8ZJecfc9ZjZM0t1mNlaSS8pptNvT7r67ibpyJP3SzI6TVC+pqNG2Be6+MVnbG5JGSCqVtMXdX03Wsi+5/TRJx5hZwzTEnsnP5x0lgtwh15cBAJAqhCwAaJsqdz+ucYOZLVNiul+7jispr7VODWEi+Xiumf3KzPq5+85G3Y5XYhrjmuTMugIzW+PuY9q4/w8k/ViJ67QWKXGW6+9KTIs82P1KnM0ar8SZLUm6QdI/3f385DVb8xr1r2jmrX1F0jZJxyoxu6K60baaRo/rlfiZZUoEuIOZpP9096ea2JanxOcMAEBaMF0QAI7cc5JyzezKhobk9UIzDuMYK5QIRi1KrmLYcE3SVCW+f+9q3MfdH3f3Qe4+wt1HSKp09zFt2T959mmIu8+XVCAprkSYaS4A/lnSpySdIumRZFtPSZuSjz/d2ntqtM+W5LVel0oKt9J/paQhZnZCsu7uySmcT0n6DzPLSbYXmVm35D5Fkpa2sR4AANqNkAUAR8jdXdL5kj5siSXclylxXdNmSTKzF5RYJOJDZrbRzGY3cZjH1Wh6n5ldbWYbJR0labGZ/S656aOSlprZm0pct3RJ8vVliWXjW5sO1+z+ST+W9P3k4z8rEZJekfTzZt77ckmVkp5z94azVD+V9BMz+5daD0sNfiXpcjN7RYkw1NwZr4bXrVXi2q3/Sb6Xp5UIgr+TtFzSa8lFNe7Qe7M1ZinxOQMAkBZ24M9YAEC6mdmLks52972ZrqUzMrPnJX3E3fdkuhYAQNdAyAKADDOzaUpc87U407V0NmbWX9KJ7j4n07UAALoOQhYAAAAABIhrsgAAAAAgQIQsAAAAAAgQIQsAAAAAAkTIAgAAAIAAEbIAAAAAIED/H3g2A5Wxicq5AAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# Perform PCA on Combat-adjusted data\n", "scaler = StandardScaler()\n", "scaled_data = scaler.fit_transform(combat_adjusted_cpm.T)\n", "pca = PCA(n_components=2)\n", "pca_results = pca.fit_transform(scaled_data)\n", "\n", "# Create a PCA DataFrame\n", "pca_df = pd.DataFrame(pca_results, columns=['PC1', 'PC2'])\n", "pca_df['Sample'] = combat_adjusted_cpm.columns\n", "pca_df = pd.merge(pca_df, metadata_subset, left_on='Sample', right_on='Sample')\n", "\n", "# Plot PCA after Combat adjustment\n", "plt.figure(figsize=(12, 10))\n", "sns.scatterplot(\n", " data=pca_df,\n", " x='PC1', y='PC2',\n", " hue='batch', style='treatment',\n", " s=100, palette='Set2'\n", ")\n", "plt.title('PCA After Combat Adjustment')\n", "plt.xlabel(f'PC1 ({pca.explained_variance_ratio_[0]*100:.2f}% Variance)')\n", "plt.ylabel(f'PC2 ({pca.explained_variance_ratio_[1]*100:.2f}% Variance)')\n", "plt.legend(title='Batch', bbox_to_anchor=(1.05, 1), loc='upper left')\n", "plt.grid()\n", "plt.tight_layout()\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "severe-macedonia", "metadata": {}, "source": [ "If the PCA plot looks the same before and after batch effect correction, this could indicate one of the following:\n", "\n", "Possible Explanations\n", "\n", "Insufficient Batch Adjustment: Combat may not have fully corrected the batch effects.\n", "This can happen if batch is highly confounded with other factors (e.g., treatment or time point).\n", "\n", "Low Impact of Batch: Batch effects might not have been a major source of variance in the data to begin with.PCA often highlights dominant sources of variation, which may be unrelated to batch.\n", "\n", "Residual Batch Effects in High Variance Genes:The Combat adjustment could have reduced batch effects in most genes, but residual effects might remain in genes with the highest variance.\n", "\n", "Biological Signals are Confounded with Batch: If batch and biological variables (e.g., treatment or time point) are correlated, Combat may struggle to disentangle them.\n" ] } ], "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.8.8" } }, "nbformat": 4, "nbformat_minor": 5 }