{
"cells": [
{
"cell_type": "markdown",
"id": "f57f260f-e5d3-424e-ae55-b11d4e5fab59",
"metadata": {},
"source": [
"# Stage 4: Intervention design\n",
"\n",
"In this tutorial, we'll walk through how to use the CASCADE model trained in\n",
"[stage 2](training.ipynb) to perform targeted intervention design, using\n",
"K562-to-erythroid differentiation as an example."
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "72c9dd12-2c29-424f-9e55-3f9166928ac2",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T09:32:18.074620Z",
"iopub.status.busy": "2025-03-23T09:32:18.074186Z",
"iopub.status.idle": "2025-03-23T09:32:22.820404Z",
"shell.execute_reply": "2025-03-23T09:32:22.819488Z",
"shell.execute_reply.started": "2025-03-23T09:32:18.074574Z"
}
},
"outputs": [],
"source": [
"import anndata as ad\n",
"import numpy as np\n",
"import pandas as pd\n",
"import scanpy as sc\n",
"import seaborn as sns\n",
"from scipy.sparse import csr_matrix\n",
"\n",
"from cascade.data import (\n",
" configure_dataset,\n",
" encode_regime,\n",
" get_all_targets,\n",
" get_configuration,\n",
")\n",
"from cascade.model import CASCADE, IntervDesign\n",
"from cascade.plot import plot_design_error_curve, plot_design_scores, set_figure_params"
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "548b5c2c-b4ce-42e7-83c5-26faee143247",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T09:32:22.822530Z",
"iopub.status.busy": "2025-03-23T09:32:22.822071Z",
"iopub.status.idle": "2025-03-23T09:32:22.832240Z",
"shell.execute_reply": "2025-03-23T09:32:22.831514Z",
"shell.execute_reply.started": "2025-03-23T09:32:22.822507Z"
},
"scrolled": true
},
"outputs": [],
"source": [
"set_figure_params()"
]
},
{
"cell_type": "markdown",
"id": "eb5caada-c122-4f0c-9962-4537722adebf",
"metadata": {},
"source": [
"## Prepare source data\n",
"\n",
"First, we will extract the source state, i.e., unperturbed K562 cells from our\n",
"preprocessed data:"
]
},
{
"cell_type": "code",
"execution_count": 3,
"id": "98243314-f2e5-4f32-a062-5f9b6c84baec",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T09:32:22.833057Z",
"iopub.status.busy": "2025-03-23T09:32:22.832875Z",
"iopub.status.idle": "2025-03-23T09:32:27.863915Z",
"shell.execute_reply": "2025-03-23T09:32:27.863055Z",
"shell.execute_reply.started": "2025-03-23T09:32:22.833039Z"
}
},
"outputs": [
{
"data": {
"text/plain": [
"AnnData object with n_obs × n_vars = 86744 × 1064\n",
" obs: 'guide_id', 'gemgroup', 'ncounts', 'knockup'\n",
" var: 'perturbed', 'highly_variable', 'highly_variable_rank', 'means', 'variances', 'variances_norm', 'selected'\n",
" uns: '__CASCADE__', 'hvg', 'log1p', 'pca'\n",
" obsm: 'X_pca', 'covariate'\n",
" varm: 'PCs'\n",
" layers: 'counts', 'interv'"
]
},
"execution_count": 3,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"adata = sc.read_h5ad(\"adata.h5ad\")\n",
"adata"
]
},
{
"cell_type": "code",
"execution_count": 4,
"id": "fdc65fb6-3cf8-470b-8ddf-954c31853eb9",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T09:32:27.865249Z",
"iopub.status.busy": "2025-03-23T09:32:27.865026Z",
"iopub.status.idle": "2025-03-23T09:32:28.005043Z",
"shell.execute_reply": "2025-03-23T09:32:28.004292Z",
"shell.execute_reply.started": "2025-03-23T09:32:27.865230Z"
}
},
"outputs": [
{
"data": {
"text/plain": [
"AnnData object with n_obs × n_vars = 11855 × 1064\n",
" obs: 'guide_id', 'gemgroup', 'ncounts', 'knockup'\n",
" var: 'perturbed', 'highly_variable', 'highly_variable_rank', 'means', 'variances', 'variances_norm', 'selected'\n",
" uns: '__CASCADE__', 'hvg', 'log1p', 'pca'\n",
" obsm: 'X_pca', 'covariate'\n",
" varm: 'PCs'\n",
" layers: 'counts', 'interv'"
]
},
"execution_count": 4,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"source = adata[adata.obs[\"knockup\"] == \"\"].copy()\n",
"source"
]
},
{
"cell_type": "markdown",
"id": "2bf889c4-558d-4a3c-a398-1f2903cadf69",
"metadata": {},
"source": [
"## Prepare target data\n",
"\n",
"Next, we use an erythroid scRNA-seq dataset from\n",
"[Xu, et al. (2022)](https://doi.org/10.1038/s41590-022-01245-8) as the\n",
"target state. The data file can be downloaded from here:\n",
"\n",
"- http://ftp.cbi.pku.edu.cn/pub/cascade-download/Xu-2022.h5ad"
]
},
{
"cell_type": "code",
"execution_count": 5,
"id": "256d5848-e696-497e-a016-77869ee27d69",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T09:32:28.006029Z",
"iopub.status.busy": "2025-03-23T09:32:28.005804Z",
"iopub.status.idle": "2025-03-23T09:32:33.886187Z",
"shell.execute_reply": "2025-03-23T09:32:33.885668Z",
"shell.execute_reply.started": "2025-03-23T09:32:28.006010Z"
}
},
"outputs": [
{
"data": {
"text/plain": [
"AnnData object with n_obs × n_vars = 74460 × 27069\n",
" obs: 'Clusters', 'batch', 'stage', 'ncounts'"
]
},
"execution_count": 5,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"target = sc.read_h5ad(\"Xu-2022.h5ad\")\n",
"target"
]
},
{
"cell_type": "markdown",
"id": "abb7eaae-b400-45d9-b14e-7600c9398d11",
"metadata": {},
"source": [
"Some genes in our source dataset could be missing from the target data, so we\n",
"manually expand the target data with zero-paddings:"
]
},
{
"cell_type": "code",
"execution_count": 6,
"id": "8a66b31e-0cfc-427f-bc4e-da25bb1d2e0c",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T09:32:33.887171Z",
"iopub.status.busy": "2025-03-23T09:32:33.886965Z",
"iopub.status.idle": "2025-03-23T09:32:33.894189Z",
"shell.execute_reply": "2025-03-23T09:32:33.893720Z",
"shell.execute_reply.started": "2025-03-23T09:32:33.887153Z"
}
},
"outputs": [
{
"data": {
"text/plain": [
"(930, 134)"
]
},
"execution_count": 6,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"common_vars = [i for i in adata.var_names if i in target.var_names]\n",
"missing_vars = [i for i in adata.var_names if i not in target.var_names]\n",
"len(common_vars), len(missing_vars)"
]
},
{
"cell_type": "code",
"execution_count": 7,
"id": "a385c2b0-2b74-4652-8045-4bf208f3b908",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T09:32:33.895042Z",
"iopub.status.busy": "2025-03-23T09:32:33.894858Z",
"iopub.status.idle": "2025-03-23T09:32:35.415799Z",
"shell.execute_reply": "2025-03-23T09:32:35.415069Z",
"shell.execute_reply.started": "2025-03-23T09:32:33.895025Z"
}
},
"outputs": [
{
"data": {
"text/plain": [
"AnnData object with n_obs × n_vars = 74460 × 27203\n",
" obs: 'Clusters', 'batch', 'stage', 'ncounts'"
]
},
"execution_count": 7,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"empty = ad.AnnData(\n",
" X=csr_matrix((target.n_obs, len(missing_vars)), dtype=target.X.dtype),\n",
" obs=target.obs,\n",
" var=pd.DataFrame(index=missing_vars),\n",
")\n",
"target = ad.concat([target, empty], axis=1, merge=\"first\")\n",
"target"
]
},
{
"cell_type": "markdown",
"id": "ffd42550-66ad-4284-835c-41559e5e99cf",
"metadata": {},
"source": [
"Similar to the training data, we backup the raw counts and log-normalize the\n",
"dataset."
]
},
{
"cell_type": "code",
"execution_count": 8,
"id": "800aa648-ead5-4d40-9854-5c0295f19bea",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T09:32:35.418310Z",
"iopub.status.busy": "2025-03-23T09:32:35.418092Z",
"iopub.status.idle": "2025-03-23T09:32:36.544245Z",
"shell.execute_reply": "2025-03-23T09:32:36.543449Z",
"shell.execute_reply.started": "2025-03-23T09:32:35.418290Z"
}
},
"outputs": [],
"source": [
"target.layers[\"counts\"] = target.X.copy()\n",
"sc.pp.normalize_total(target, target_sum=1e4)\n",
"sc.pp.log1p(target)"
]
},
{
"cell_type": "markdown",
"id": "814ccbd3-0811-4f8d-88ff-84c41b36bb07",
"metadata": {},
"source": [
"We'd also need to configure the target dataset. Of note, the target data\n",
"configuration should provide the same covariate as the training data. However,\n",
"it obviously does not have batch info compatible with the training data.\n",
"\n",
"Here we are using a rough approach of simply taking the covariate\n",
"of random training cells. Another possibility is to take the covariate of\n",
"training cells most similar to the target cells, based on the expression levels\n",
"of housekeeping genes."
]
},
{
"cell_type": "code",
"execution_count": 9,
"id": "d709a98f-23d6-4989-917d-bdb2635d47d9",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T09:32:36.545162Z",
"iopub.status.busy": "2025-03-23T09:32:36.544968Z",
"iopub.status.idle": "2025-03-23T09:32:36.551924Z",
"shell.execute_reply": "2025-03-23T09:32:36.551201Z",
"shell.execute_reply.started": "2025-03-23T09:32:36.545143Z"
}
},
"outputs": [],
"source": [
"target.obsm[\"covariate\"] = source.obsm[\"covariate\"][\n",
" np.random.choice(source.n_obs, target.n_obs)\n",
"]"
]
},
{
"cell_type": "code",
"execution_count": 10,
"id": "8f197cf9-b2d6-4d34-bae2-a6ed74f8cfc6",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T09:32:36.552785Z",
"iopub.status.busy": "2025-03-23T09:32:36.552604Z",
"iopub.status.idle": "2025-03-23T09:32:36.558087Z",
"shell.execute_reply": "2025-03-23T09:32:36.557400Z",
"shell.execute_reply.started": "2025-03-23T09:32:36.552767Z"
}
},
"outputs": [
{
"data": {
"text/plain": [
"{'regime': None,\n",
" 'covariate': 'covariate',\n",
" 'size': 'ncounts',\n",
" 'weight': None,\n",
" 'layer': 'counts'}"
]
},
"execution_count": 10,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"configure_dataset(\n",
" target, use_covariate=\"covariate\", use_size=\"ncounts\", use_layer=\"counts\"\n",
")\n",
"get_configuration(target)"
]
},
{
"cell_type": "markdown",
"id": "18c6c5a9-1c36-4019-b76e-8f0116973622",
"metadata": {},
"source": [
"## Define gene weights\n",
"\n",
"Given that the source and target data were produced with different experimental\n",
"protocols in different studies, batch effect would be a substantial problem.\n",
"To mitigate this issue, we give higher weight to known erythroid markers when\n",
"comparing counterfactual states with the target state during intervention\n",
"design, which helps avoid learning batch effect between source and target."
]
},
{
"cell_type": "code",
"execution_count": 11,
"id": "846c3207-d75a-4dfc-b3b7-55daba11b294",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T09:32:36.558885Z",
"iopub.status.busy": "2025-03-23T09:32:36.558698Z",
"iopub.status.idle": "2025-03-23T09:32:36.571089Z",
"shell.execute_reply": "2025-03-23T09:32:36.570380Z",
"shell.execute_reply.started": "2025-03-23T09:32:36.558868Z"
}
},
"outputs": [
{
"data": {
"text/plain": [
"25"
]
},
"execution_count": 11,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"markers = [\n",
" \"AHSP\",\n",
" \"ALAS2\",\n",
" \"ALDOA\",\n",
" \"BLVRB\",\n",
" \"BPGM\",\n",
" \"CLIC1\",\n",
" \"ENO1\",\n",
" \"GYPA\",\n",
" \"GYPB\",\n",
" \"HAMP\",\n",
" \"HBA1\",\n",
" \"HBA2\",\n",
" \"HBB\",\n",
" \"HBD\",\n",
" \"HBE1\",\n",
" \"HBG1\",\n",
" \"HBG2\",\n",
" \"HBZ\",\n",
" \"HEMGN\",\n",
" \"LDHA\",\n",
" \"PRDX1\",\n",
" \"PRDX2\",\n",
" \"SLC25A37\",\n",
" \"SLC4A1\",\n",
" \"SMIM1\",\n",
"]\n",
"assert not set(markers) - set(adata.var_names)\n",
"assert not set(markers) - set(target.var_names)\n",
"len(markers)"
]
},
{
"cell_type": "code",
"execution_count": 12,
"id": "7306b5f2-d9da-48d9-b490-289531c9d20c",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T09:32:36.571774Z",
"iopub.status.busy": "2025-03-23T09:32:36.571600Z",
"iopub.status.idle": "2025-03-23T09:32:36.576925Z",
"shell.execute_reply": "2025-03-23T09:32:36.576248Z",
"shell.execute_reply.started": "2025-03-23T09:32:36.571758Z"
}
},
"outputs": [
{
"data": {
"text/plain": [
"905"
]
},
"execution_count": 12,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"non_markers = [i for i in common_vars if i not in markers]\n",
"len(non_markers)"
]
},
{
"cell_type": "code",
"execution_count": 13,
"id": "defa71f8-50c0-423e-9446-59063c9adecc",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T09:32:36.577700Z",
"iopub.status.busy": "2025-03-23T09:32:36.577529Z",
"iopub.status.idle": "2025-03-23T09:32:36.584341Z",
"shell.execute_reply": "2025-03-23T09:32:36.583652Z",
"shell.execute_reply.started": "2025-03-23T09:32:36.577685Z"
}
},
"outputs": [],
"source": [
"target.var[\"weight\"] = 0.0\n",
"target.var.loc[non_markers, \"weight\"] = len(common_vars) / len(non_markers) / 2\n",
"target.var.loc[markers, \"weight\"] = len(common_vars) / len(markers) / 2"
]
},
{
"cell_type": "code",
"execution_count": 14,
"id": "1b05ef02-d9be-4fdd-8124-7a72ee31f461",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T09:32:36.585121Z",
"iopub.status.busy": "2025-03-23T09:32:36.584950Z",
"iopub.status.idle": "2025-03-23T09:32:36.592525Z",
"shell.execute_reply": "2025-03-23T09:32:36.591843Z",
"shell.execute_reply.started": "2025-03-23T09:32:36.585105Z"
}
},
"outputs": [
{
"data": {
"text/plain": [
"SAMD11 0.513812\n",
"HES4 0.513812\n",
"ISG15 0.513812\n",
"RNF223 0.513812\n",
"TNFRSF4 0.513812\n",
"Name: weight, dtype: float64"
]
},
"execution_count": 14,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"target.var.loc[non_markers, \"weight\"].head()"
]
},
{
"cell_type": "code",
"execution_count": 15,
"id": "184f9d06-b1ca-4250-8d67-85af7e5f5786",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T09:32:36.593369Z",
"iopub.status.busy": "2025-03-23T09:32:36.593199Z",
"iopub.status.idle": "2025-03-23T09:32:36.599672Z",
"shell.execute_reply": "2025-03-23T09:32:36.598995Z",
"shell.execute_reply.started": "2025-03-23T09:32:36.593353Z"
}
},
"outputs": [
{
"data": {
"text/plain": [
"AHSP 18.6\n",
"ALAS2 18.6\n",
"ALDOA 18.6\n",
"BLVRB 18.6\n",
"BPGM 18.6\n",
"Name: weight, dtype: float64"
]
},
"execution_count": 15,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"target.var.loc[markers, \"weight\"].head()"
]
},
{
"cell_type": "markdown",
"id": "253d3b75-0754-4ea7-9565-12a12fc1ef5d",
"metadata": {},
"source": [
"In this case we scaled the gene weights so that the total weight of marker genes\n",
"equals and the total weight of non-marker genes, and genes that were zero-padded\n",
"in the target data were given zero weight, so they will not bias the result.\n",
"\n",
"## Specify candidate genes\n",
"\n",
"We can also provide a candidate gene pool, e.g., we'll just use genes perturbed\n",
"in the CRISPRa dataset that show higher expression in the target cells:"
]
},
{
"cell_type": "code",
"execution_count": 16,
"id": "f85e1dec-0e0b-4d29-8336-4738b7081895",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T09:32:36.600604Z",
"iopub.status.busy": "2025-03-23T09:32:36.600425Z",
"iopub.status.idle": "2025-03-23T09:32:37.439052Z",
"shell.execute_reply": "2025-03-23T09:32:37.438288Z",
"shell.execute_reply.started": "2025-03-23T09:32:36.600588Z"
}
},
"outputs": [
{
"data": {
"text/plain": [
"['ARRDC3',\n",
" 'ATL1',\n",
" 'BPGM',\n",
" 'CDKN1A',\n",
" 'CDKN1B',\n",
" 'CDKN1C',\n",
" 'CEBPA',\n",
" 'CEBPB',\n",
" 'COL1A1',\n",
" 'CSRNP1',\n",
" 'EGR1',\n",
" 'ETS2',\n",
" 'FOSB',\n",
" 'FOXF1',\n",
" 'HNF4A',\n",
" 'HOXA13',\n",
" 'IER5L',\n",
" 'IGDCC3',\n",
" 'IRF1',\n",
" 'JUN',\n",
" 'KLF1',\n",
" 'MAML2',\n",
" 'MAP2K3',\n",
" 'MAP2K6',\n",
" 'MAP3K21',\n",
" 'MAP4K5',\n",
" 'MEIS1',\n",
" 'MIDEAS',\n",
" 'PRDM1',\n",
" 'PRTG',\n",
" 'PTPN13',\n",
" 'RUNX1T1',\n",
" 'SGK1',\n",
" 'SLC38A2',\n",
" 'SLC4A1',\n",
" 'SNAI1',\n",
" 'SPI1',\n",
" 'TBX3',\n",
" 'TGFBR2',\n",
" 'TSC22D1',\n",
" 'ZBTB1',\n",
" 'ZBTB10']"
]
},
"execution_count": 16,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"all_targets = sorted(get_all_targets(adata, \"knockup\"))\n",
"target_mask = (\n",
" target[:, all_targets].to_df().mean() > adata[:, all_targets].to_df().mean()\n",
")\n",
"candidates = target_mask.index[target_mask].to_list()\n",
"candidates"
]
},
{
"cell_type": "markdown",
"id": "f57d86ed-7264-4759-aba4-85f703e6ddc8",
"metadata": {},
"source": [
"## Run intervention design\n",
"\n",
"> (Estimated time: 10 min – 20 min, depending on computation device)\n",
"\n",
"Finally, we can load our trained model for intervention design:"
]
},
{
"cell_type": "code",
"execution_count": 17,
"id": "6f93041e-f931-44e6-969d-f88489e7dc80",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T09:32:37.440036Z",
"iopub.status.busy": "2025-03-23T09:32:37.439834Z",
"iopub.status.idle": "2025-03-23T09:32:37.931893Z",
"shell.execute_reply": "2025-03-23T09:32:37.931188Z",
"shell.execute_reply.started": "2025-03-23T09:32:37.440019Z"
}
},
"outputs": [],
"source": [
"cascade = CASCADE.load(\"tune.pt\")"
]
},
{
"cell_type": "markdown",
"id": "14fffd8e-849e-41d3-b6a2-bb3e4ad90ee6",
"metadata": {},
"source": [
"For the sake of speed, we will subsample both the source and target data to\n",
"5,000 cells."
]
},
{
"cell_type": "code",
"execution_count": 18,
"id": "b886d7c0-31b4-45a2-b785-7b069bf3cfb2",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T09:32:37.932958Z",
"iopub.status.busy": "2025-03-23T09:32:37.932694Z",
"iopub.status.idle": "2025-03-23T09:32:38.133815Z",
"shell.execute_reply": "2025-03-23T09:32:38.133298Z",
"shell.execute_reply.started": "2025-03-23T09:32:37.932933Z"
}
},
"outputs": [],
"source": [
"sc.pp.subsample(source, n_obs=5000)\n",
"sc.pp.subsample(target, n_obs=5000)"
]
},
{
"cell_type": "markdown",
"id": "5b72dbcb-0081-44e7-bc06-123092f89cfc",
"metadata": {},
"source": [
"CASCADE model provides a dedicated [design](api/cascade.model.CASCADE.design.rst)\n",
"method for intervention design, which uses differentiable optimization to\n",
"optimize interventions that produces effects more similar to the target state.\n",
"\n",
"We'll need to pass the source and target datasets, along with a candidate gene\n",
"pool, a maximal combination order (`design_size=1` for designing single-gene\n",
"perturbations), as well as the target gene weight we just assigned."
]
},
{
"cell_type": "code",
"execution_count": 19,
"id": "14cc407b-4d29-4b6f-8e3f-c82855234f45",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T09:32:38.134753Z",
"iopub.status.busy": "2025-03-23T09:32:38.134560Z",
"iopub.status.idle": "2025-03-23T09:47:16.615367Z",
"shell.execute_reply": "2025-03-23T09:47:16.614649Z",
"shell.execute_reply.started": "2025-03-23T09:32:38.134735Z"
}
},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"\u001b[32m17:32:38.143\u001b[0m | \u001b[33m\u001b[1mWARNING \u001b[0m | \u001b[33m3810563\u001b[0m:\u001b[36mmodel\u001b[0m:\u001b[36malign_vars\u001b[0m - \u001b[33m\u001b[1m26139 variables are not in the `scaffold` and will thus be ignored.\u001b[0m\n"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"\u001b[32m17:32:38.151\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[33m3810563\u001b[0m:\u001b[36mutils\u001b[0m:\u001b[36mautodevice\u001b[0m - \u001b[1mUsing GPU [1] as computation device.\u001b[0m\n"
]
},
{
"data": {
"text/html": [
"╭────────────────────────────── cascade-reg ──────────────────────────────╮\n",
"│ │\n",
"│ Training on 1064 variables with 12294 scaffold edges and 5000 samples │\n",
"│ │\n",
"╰──────────────────────────────── v0.4.0 ─────────────────────────────────╯\n",
"\n"
],
"text/plain": [
"╭────────────────────────────── cascade-reg ──────────────────────────────╮\n",
"│ │\n",
"│ Training on \u001b[1;35m1064\u001b[0m variables with \u001b[1;35m12294\u001b[0m scaffold edges and \u001b[1;35m5000\u001b[0m samples │\n",
"│ │\n",
"╰──────────────────────────────── v0.4.0 ─────────────────────────────────╯\n"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"\u001b[32m17:32:38.325\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[33m3810563\u001b[0m:\u001b[36mcore\u001b[0m:\u001b[36mfit_stage\u001b[0m - \u001b[1mNumber of topological generations: [68, 88, 70, 101]\u001b[0m\n"
]
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"\n",
" | Name | Type | Params | Mode \n",
"------------------------------------------------------\n",
"0 | scaffold | Edgewise | 49.2 K | eval \n",
"1 | sparse | L1 | 0 | eval \n",
"2 | acyc | SpecNorm | 0 | eval \n",
"3 | kernel | RBF | 0 | eval \n",
"4 | latent | EmbLatent | 6.3 K | eval \n",
"5 | lik | NegBin | 0 | eval \n",
"6 | func | Func | 7.9 M | eval \n",
"7 | design | IntervDesign | 8.6 K | train\n",
" | other params | n/a | 8.5 K | n/a \n",
"------------------------------------------------------\n",
"10.7 K Trainable params\n",
"8.0 M Non-trainable params\n",
"8.0 M Total params\n",
"31.955 Total estimated model params size (MB)\n",
"1 Modules in train mode\n",
"16 Modules in eval mode\n"
]
},
{
"data": {
"application/vnd.jupyter.widget-view+json": {
"model_id": "56b07ca487474564abdd22ec32d2e4fa",
"version_major": 2,
"version_minor": 0
},
"text/plain": [
"Training: | …"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
"Restoring best model: log_dir/design/lightning_logs/version_0/checkpoints/epoch=8-step=1200.ckpt.\n",
"\n"
],
"text/plain": [
"Restoring best model: log_dir/design/lightning_logs/version_0/checkpoints/\u001b[33mepoch\u001b[0m=\u001b[1;36m8\u001b[0m-\u001b[33mstep\u001b[0m=\u001b[1;36m1200\u001b[0m\u001b[1;36m.\u001b[0mckpt.\n"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"scores, design = cascade.design(\n",
" source, target, pool=candidates, design_size=1, target_weight=\"weight\"\n",
")"
]
},
{
"cell_type": "markdown",
"id": "13a2df02-628e-4871-b844-31221d441aef",
"metadata": {},
"source": [
"The `design` method returns two objects:\n",
"\n",
"- The `scores` object is a data frame containing scores of candidate genes\n",
" (or gene combinations if `design_size` was larger than 1).\n",
"- The `design` object is an [IntervDesign](api/cascade.nn.IntervDesign.rst)\n",
" module that contains both the scores and optimized interventional scales and\n",
" biases for the designed interventions, which can also be saved and loaded just\n",
" like the CASCADE model."
]
},
{
"cell_type": "code",
"execution_count": 20,
"id": "f87871fb-601c-4124-aaf3-4311ef9f7c75",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T09:47:16.617144Z",
"iopub.status.busy": "2025-03-23T09:47:16.616547Z",
"iopub.status.idle": "2025-03-23T09:47:16.629375Z",
"shell.execute_reply": "2025-03-23T09:47:16.628825Z",
"shell.execute_reply.started": "2025-03-23T09:47:16.617116Z"
}
},
"outputs": [
{
"data": {
"text/html": [
"\n",
"\n",
"
\n",
" \n",
" \n",
" | \n",
" score | \n",
"
\n",
" \n",
" \n",
" \n",
" | KLF1 | \n",
" 10.006369 | \n",
"
\n",
" \n",
" | BPGM | \n",
" 4.527291 | \n",
"
\n",
" \n",
" | ATL1 | \n",
" 2.883361 | \n",
"
\n",
" \n",
" | PTPN13 | \n",
" 1.928048 | \n",
"
\n",
" \n",
" | ARRDC3 | \n",
" -0.666598 | \n",
"
\n",
" \n",
"
\n",
"
"
]
},
"metadata": {
"image/png": {
"height": 302,
"width": 318
},
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"plot_design_error_curve(curve, cutoff=cutoff)"
]
},
{
"cell_type": "markdown",
"id": "70fbf16f-a787-44c6-b3e8-c8f255b43e05",
"metadata": {},
"source": [
"We can also visualize the scores with a scatter plot:"
]
},
{
"cell_type": "code",
"execution_count": 24,
"id": "8d427685-9732-4fea-801a-a07ba1ad311a",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T10:07:45.751653Z",
"iopub.status.busy": "2025-03-23T10:07:45.751373Z",
"iopub.status.idle": "2025-03-23T10:07:46.026165Z",
"shell.execute_reply": "2025-03-23T10:07:46.025373Z",
"shell.execute_reply.started": "2025-03-23T10:07:45.751631Z"
}
},
"outputs": [
{
"data": {
"text/plain": [
""
]
},
"execution_count": 24,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAngAAAJcCAYAAABqqtkaAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjkuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8hTgPZAAAACXBIWXMAABibAAAYmwFJdYOUAACQL0lEQVR4nOzdd3xb1f3/8ddHw9uOs4eTkAAJIeyEPUqAQthQKDOstpTSFrpbWlooUMpov235FfjSLyuBssouGwoUQtgBElYYSchw9na8Len8/riyYju2LNmSZcnv5+MhdKV77rkf6djh43PvOcecc4iIiIhI7vBlOgARERERSS0leCIiIiI5RgmeiIiISI5RgiciIiKSY5TgiYiIiOQYJXgiIiIiOUYJnoiIiEiOUYInIiIikmOU4ImIiIjkGCV4IiIiIjlGCZ6IiIhIjlGCJyIiIpJjlOCJiIiI5BgleDnGzGab2exMxyEiIiKZE8h0AJJywyoqKioAl0jhjRs3AlBeXp6+iCRl1F7ZRe2VXdRe2aWPtJd19UD14ImIiIjkGCV4IiIiIjlGCZ6IiIhIjlGCJyIiIpJjlOCJiIiI5BgleCIiIiI5xpxLaDYNyRJmVllRUVFRWVmZ6VBERESkezRNioiIiIh4lOCJiIiI5BgleCIiIiI5RgmeiIiISI5RgiciIiKSY5TgiYiIiOQYJXgiIiIiOUYJXh8XiUSIRCKZDkMSpPbKLmqv7KL2yi5qr/iU4PVxVVVVVFVVZToMSZDaK7uovbKL2iu7qL3iU4InIiIikmOU4ImIiIjkGCV4IiIiIjlGCZ6IiIhIjlGCJyIiIpJjlOCJiIiI5BgleCIiIiI5RgmeiIiISI5RgiciIiKSY5TgiYiIiOQYJXgiIiIiOUYJnoiIiEiOUYInIiIikmOU4ImIiIjkGCV4IiIiIjlGCZ6IiIhIjglkOgDJrGAwmOkQJAlqr+yi9souaq/sovaKTwleH1dcXJzpECQJaq/sovbKLmqv7KL2ik+XaEVERERyjBI8ERERkRyjBE9EREQkxyjBExEREckxSvBEREREcowSPBEREZEcowRPREREJMcowevjGhoaaGho6FYdV199NYMGDWLQoEHU1dV1WG7ZsmVMnjyZQYMGMXToUJ544onYvj/84Q+xOpKJZ/Xq1bHjOnscc8wxWx2/atUqnn76aa666iqOP/54Ro8eHSs/d+7c5L6IHpCK9pKeo/bKLmqv7KL2ik8THfdxzQlZfn5+l+uora1l3bp1ccssXLiQr3/963z11Vfk5eXxwAMPcPzxxydVR3uccwkft2nTpq3e23vvvVmyZEm75cPhcNLxpFsq2kt6jtoru6i9sovaKz4leJJ2n3zyCYcffjgrVqygqKiIxx9/nMMPPzzl5/ne977HWWed1eH+srKydt83M3bYYQcmT55MOBzmgQceSHlsIiIiPUkJnqTV7NmzmTp1KuvXr6e8vJynn36a/fffPy3nGjt2LAceeGBSxzz66KOMHz+e0tJSAB588EEleCIikvWU4EnazJw5k2OPPZbNmzczZMgQXnjhBXbbbbdMh9XK5MmTMx2CiIhIymmQhaTFM888w9SpU9m8eTOjRo3itdde63XJnYiISK5Sgicp9+CDD3LCCSdQX1/PuHHjmDVrFuPHj890WCIiIn2GLtFKSt1xxx38+Mc/JhKJsNtuu/H8888zdOjQHjl3ZWUl7733Xrv7fD4fe+yxR4/EISIikmlK8CSlLr74YgD22msvnn/+efr3799j577pppu46aab2t1XVFRETU1Nj8UiIiKSSUrwJC2WLVvG2rVrezTBCwaDHc6HVFRU1GNxiIiIZJruwZOUuvXWW/H5fCxfvpxDDjmE+fPn99i5//CHP7B58+Z2H6tWreqxOERERDJNCZ6k1FlnncVdd92Fz+dj2bJlTJkypUeTPBEREVGC14qZVZjZiWZ2tZk9b2brzMxFH4u6UJ+Z2alm9oSZLTazejNbaWavmdmPzawkDR8j45qTPDNj2bJlHHLIISxYsCDTYYmIiPQZugcvysyuAy5JYX1DgQeAKW12DY0+DgR+bGanOudmp+q8vcVZZ52Fc45zzz2XyspKpkyZwiuvvMJ2222X6dBERERynnrwtiho5731XanIzIqAp9iS3NUB9wCXAzcDzTeEjQWeM7NxXTlPb3f22WfHevIqKys55JBDWLhwYabDEhERyXnqwdtiPfACMBt4N/ocAL7qQl2XA3tGt78CjnDOxW5EM7NLgEeAqcBA4E7goC5HnmYLPljA87c+x7LPK3EORowbwZHfO4rtJ2/f6bFnn302zjnOO+88li5dGuvJ23bbbXsgchERkb5JCV6Uc+6qtu+Z2Zhk6zGzIcCPWrw1rWVyFz1XjZmdDswDhgEHmtlRzrlnkz1fOn01dyF/PfuvXmIXcYRDYQA+fe1TXr77Za599Tp22GeHTus555xzcM7xrW99i6VLl3LIIYfwyiuvMHbs2HbLf/zxx+Tl5cWtc9ttt6W4uDj5D9XGypUrWbNmTez1kiVLYttffvklwWAw9nrkyJE9Ou2LiIhIVynBS72TgcLo9svOuTfbK+Sc22hmtwBXRt86C+g1Cd5nb87jd4f9lqb6JpxzrfZFwmEi4TANNfUJ13fuuecC8K1vfYslS5bEevLaS/L23HPPrd5r6/XXX2f//fdP+Pwdufnmm7n66qvb3Xf66ae3ej1jxozY5xAREenNlOCl3nEtth/tpOyjbEnwjjUzv3MunJ6w2tfeBMBVazfx+6mX01jXmFAdfr+/0x438JI85xwXXnghK1eu5Mgjj+TFF19k1KhRCdfRzMxi283H+f3+hI9PNnbwljvLNE3YnF3UXtlF7ZVd1F7xWdveGdkieom2+R68xc65MQkcswLvsivAJOfcB3HKGrAB6Bd9a7xz7ssuB+zVWVlRUVFRWVnZ5Toevu4hHrjq/k4TvB/eejFTvzu1y+cRERGRuKzzIu3LfJdEDjGzcrYkdwBxJ39zXnbdcljpjmkIKynhcJgnbvh3Qr13j17/MKsXr+6BqERERCQZukSbWi2Tu2rnXFUCxywD9mjn+IxYs2QNm9dvTqjsigUruGDb86mYMJLCkgJ+9eCvGbLNkDRHKCIiIp1RgpdaLVemqEnwmNoW26WJHGBm8a6/DnPOsXHjxgRP39qa5Wvw+/2EmxK7FTASibD00yXkFeax7Ktl5PVL/D66ZoWFheTn58ctU1NTQ1NTU9J1t1VWVtbpvXRd/e5aCgQClJTEX6iksbGR2trauGUSoe+ve/T9dY++v+7R99c9ffX7S4Qu0aZWYYvtxEYoQEMHx2dEYWkhoVAo6ePCoTBFZRkPX0RERFAPXqrVtdhOtCur5Z8edR2WasE5N7KjfWZWaWYV5eXlCZ6+tbKyMsqH9mf9snVJHVc6oJRxu49P20jTVMx5l6iufnfJysvLS2rkcHfo++sefX/do++ve/T9dU8ufn+JUA9ealW32E70J6rlOO/Ebn5LI5/Px4k/O5H8ovhd3i3lFeZxws9O7BXTiIiIiIgSvFRb1WK7xMwSuaeuooPjM+br3zqcYH6w1Vxz8QTzgxz+nSPSHJWIiIgkSgleCjnnNgArW7y1Xbzy0XnwWi7K+mk64oqnurqa6urqVu+V9C/hDy/9kfzifMzXcZJnZhQUF3DVf66mbGBZukMV2m8v6b3UXtlF7ZVd1F7xKcFLvTkttjtbS2siWyY53kzrOfF6RCgUandQxXZ7bMdfZ9/AjvvviD/oJ5C35XbNQF4Af9DPDvtN4C/v/o1xe47ryZD7tI7aS3ontVd2UXtlF7VXfBpkkXpPAkdGt08C/jdO2W+02H6qp5cp68zIHUZy3Wt/YvmXy3jxzhep/KISHFSMr+Dr3z6civEVnVciIiIiPU4JXuo9AvwFKAAONbN9nHNvty1kZv2AH7R4694eii9pI8ZVcM6152Y6DBEREUmQLtGmmHNuFXBj9KUB95lZq3vxzKwIuA8YHn3rDefc0z0XpYiIiOQy9eBFmdlg4JY2b7ecwmSwmT3c9jjn3Dfbqe4q4DBgEt4gio+ix34BDAFOYcuyZBuA73QvehEREZEtlOBtUQycHGd/USf7Y5xz1WZ2DPAAcDDeChVnt1N0EXCac+6z5EIVERER6Zgu0aaJc24lcAhwBvAUUIm3fNka4HXgZ8Cuzrl3MhakiIiI5CT14EU55xbh3TOXyjodXi/eA6msV0RERCQe9eCJiIiI5BgleCIiIiI5RgmeiIiISI5RgiciIiKSY5TgiYiIiOQYJXgiIiIiOUYJnoiIiEiOUYInIiIikmOU4ImIiIjkGPMWW5BcYWaVFRUVFZWVlZkORURERLqnyytsqQdPREREJMcowRMRERHJMUrwRERERHKMEjwRERGRHKMET0RERCTHKMETERERyTFK8ERERERyjBK8Pi4SiRCJRDIdhiRI7ZVd1F7ZRe2VXdRe8SnB6+OqqqqoqqrKdBiSILVXdlF7ZRe1V3ZRe8WnBE9EREQkxyjBExEREckxSvBEREREcowSPBEREZEcowRPREREJMcowRMRERHJMUrwRERERHKMEjwRERGRHKMET0RERCTHKMETERERyTFK8ERERERyjBI8ERERkRyjBE9EREQkxyjBExEREckxSvBEREREcowSPBEREZEcE8h0AJJZwWAw0yFIEtRe2UXtlV3UXtlF7RWfErw+rri4ONMhSBLUXtlF7ZVd1F7ZRe0Vny7RioiIiOQYJXgiIiIiOUYJnoiIiEiOUYInIiIikmOU4ImIiIjkGCV4IiIiIjlGCZ6IiIhIjlGC18c1NDTQ0NCQ6TAkQWqv7KL2yi5qr+yi9opPEx33cXV1dQDk5+dnOBJJhNoru6i9sovaK7uoveJTD56IiIhIjlGCJyIiIpJjlOCJiIiI5BgleCIiIiI5RgmeiIiISI5RgiciIiKSY5TgiYiIiOQYJXgiIiIiOUYJnoiIiEiOUYInIiIikmOU4ImIiIjkGCV4IiIiIjlGCZ6IiIhIjlGCJyIiIpJjlOCJiIiI5BgleCIiIiI5JpDpACSzioqKMh2CJEHtlV3UXtlF7ZVd1F7xKcHr4/Ly8jIdgiRB7ZVd1F7ZRe2VXdRe8ekSrYiIiEiOUYInIiIikmOU4ImIiIjkGCV4IiIiIjmmRwdZmNloYDBQBCx0zi3ryfOLiIiI9AVp78Ezs13MbIaZrQa+At4BXgFOa6fsX8zsVjP7RbrjEhEREclVae3BM7NLgSsAP2AtdrkODokA5wONZjbdObcunfGJiIiI5KK09eBFk7ur8ZLICPA88OdODrst+hwEjk9XbLJFdXU11dXVmQ5DEqT2yi5qr+yi9souaq/40pLgmdl4vJ47B3wJ7OqcO8o5d0m845xzXwCfR18elo7YpLVQKEQoFMp0GJIgtVd2UXtlF7VXdlF7xZeuHrwf4vXchYDjnHPzkjj2HbzLubumIzARERGRXJeue/AOw+u9eybaK5eMxdHnkakNSaRjjY2NfPnll6xdu5ampibKy8vZfvvtKS8vz3RoIiIiSUtXgtecnH3QhWNros8lKYpFeonzzz+ff/3rX+3uKygoYOTIkRx22GGceOKJHHjggR3Ws+uuu/LVV1+1W0dFRQVHHnkkxx9/PPvttx9m1k4NHuccH3zwATNmzOC+++5j3bqtx/SMHz+ek08+me985ztst912W+3/4osvmDx5cuz1hAkTeOeddzo97x577MGCBQti73344YeMHTu2w2NERESSka5LtM0rADd14dh+0WfdOZlj6uvrYzfFtn2sXbuWOXPm8Je//IWDDjqI7373u4TD4Xbrqamp6bCOuXPncv3113PAAQdw9tlnd3h/xsqVK5kyZQqTJ0/mxhtvZN26deTl5bHNNtswYcIEBg0aBHgJ3LXXXssOO+zA9ddfv1U9zrlWMcyePZv3338/7vfw5ptvMnfu3FbHOdfRwHIREZHkpSvBWx19Ht2FYydFn1ekKBbpZfr3789nn30We3z66afMnDmTv//977FerNtvv52HHnoobj2HHHJIrI558+bx7rvvcv311zNypNeBfO+993LddddtddyqVavYZ599mDlzJgAHH3ww06dPZ+3atSxatIh58+axZs0ali9fzh133MGuu+5KOBzmzTffjBvPgAEDAJgxY0bccs37Bw4cGLeciIhIV6UrwXsfb6DE1GQOiq50cSje/XtvpCEu6QUCgQA77LBD7LHjjjty0EEHcfHFF/POO+9QXFwMwA033BC3ntLS0lgdEyZMYM899+RXv/oV77zzTizZuvXWW1v1BDrnmDZtGkuWLAHg5ptv5pVXXuG8886jtLS0Vf3Dhw/n29/+Nh988AF33HEHJSXx7xqYNm0aAPfddx8NDQ3tlqmtreWBBx5oVV5ERCTV0pXgPR593sbMfpjIAWYWAKaz5b7AR9IQl/RygwYN4phjjgHgk08+6VIdw4cPjyVPS5cu5d13343tmzNnDi+99BIAF154IT/4wQ86rc/n8/Htb3+bW2+9NW65Aw88kG233Zb169fz1FNPtVvmscceY/PmzUycOJE999wz0Y8kIiKSlHQlePcC8/F68f5mZj82sw7PZWZ7AK8CU/B6795zzj2XptgkS3TnEubuu+8e216xYsvV/vvvvz+2/dOf/jSpOouKiuLu9/l8nHvuuUDHl2mnT58OwHnnnRd3IIaIiEh3pCXBc86FgdPxBkr4gb8ClWZ2T4tiR5jZLWY2F5gN7IuXEG4CzkxHXNL7rV27lqeffhqA/fffv8v1tBwR23zJNxKJxEbxjh8/nvHjx3cj0vadc845ADz77LOsXLmy1b7Fixfz8ssv4/P5OOuss1J+bhERkWZpW4vWOfe+mR0J3A+MAoYBZ7BlHdrDWxRv7sqoBE5yzs1PV1ySeaFQiM8++yz2OhKJsGbNGubOncsNN9xATU0NZWVlXH755V2q3zkXS+Ty8vJi05jMmzcvdu/d0Ucf3c1P0b4xY8ZwyCGH8N///pd7772Xn//857F9d999N845jjrqKIYPH56W84uIiEAaEzwA59wbZrYr8GPgu0BFB0XXAv8H/NU5tyGdMUnmbdiwgR133LHdfc33u/3hD39gxIgRSde9du1arrjiCt577z0ATjzxxNil3tWrV8fKjRkzJvnAE/Stb32L//73v8yYMYOf/exnmBnOudhl229961tpO7eIiAikOcEDcM5tAq4CrjKz7YFdgIHRc68DPnfOfZjuOKR3aXt/XX19PTU1NUQiEe69915GjBjBpZdeSmFhYYd1zJo1iylTpgBb5qObO3dubNTsgAEDuPLKK2Plq6qqYtvxRsReddVVHY7gnTVrFhMnToz72U466SR+8IMf8PHHH/P+++8zefJkXnvtNRYuXEj//v057rjj4h4vIiLSXWlJ8MzsR9HNWufc7c3vRy+96vJrHzd48OBWvWnNFi9ezKOPPspvf/tbrr76at59912eeeYZfL72bxVdv349r7766lbvDxw4kGOOOYarrrqKbbbZJvZ+y2lQampqtjquWV1dHRs2tN+R3NHkyy0VFxdz6qmncueddzJjxgwmT54c670788wzyc/P77QOERGR7khXD94NePfaPQDcHr+oiGebbbbhpz/9Kfn5+fzwhz/k+eef59lnn41Nm9LWpEmTWt3jVlhYyIgRI9hzzz3x+/1blR8yZEhse/HixVvtb/bjH/+41Rx1Tz31FL/5zW+S+iznnXced955J/fddx9XXnklDz74YOx9ERGRdEtXgrcZby3ZeWmqX3LYtGnTuOiii3DO8eijj3aY4I0cOZIzz0x8wPXEiRMZOXIklZWVPPPMM/zlL39pt9ywYcMYNmxY7PXcuXOT+wB4c+Jtt912LFiwgPPOO4+amhp23nnnVuvWioiIpEu65sFbHn3Oi1tKpB39+vWL3Xv31Vdfpaxen8/HaaedBsBnn33Gl19+mbK62zKzWG/dk08+CWjuOxER6TnpSvBexpv6ZO801S8pUl5eTnl5eabDaGXx4sXU1tYCpDy2M844I7b997//PaV1t3XOOefEEjq/35+Spcl6Y3tJx9Re2UXtlV3UXvGlK8H7PyAEHGZmuiYlSWk5gvUb3/hGSuueNGlSbOTtTTfdxG233ZbS+lsaPXo0zz//PI899hjPP/98q8u+IiIi6ZSWe/Cccx+a2aXAn4CnzOwU59ysdJxLsk/biY7Bmybl448/5pFHHuHxxx8HYMSIEZx66qkpPbeZcd9997H33ntTWVnJBRdcwIMPPsi5557LSSed1Go5sqamJmbOnMntt3d9nNDhhx/eeSEREZEUS9c0KV8D3sHryfse8KqZzQSeBT7HW44s0lk9zrmZ6YhPMiveRMfNxo4dy2OPPZaWKUWGDx/O22+/zSmnnMIbb7zBiy++yIsvvsgFF1xARUUFpaWlbN68mcWLF9PU1AR4ieGFF17Itttum/J4REREUi1do2hfYcuSZA7vfryvRR+JcvTARMzScwoKClrNRddSMBhkwIABHHTQQRx11FEcc8wxrXrTWiopKaG0tDTuJMidGTFiBLNmzeLdd99lxowZ3H///WzcuJH581tP07jtttty0kkncf7557PDDjtsVY/P54t9pkAguR/XYDAYO1aDL0REJJXMOdd5qWQrNeu0dy4Bzjm39WRmEpeZVVZUVFRUVlZmOpSsUl9fz7x581i3bh1NTU3079+fbbfdttXceSIiIj2sy3/9p6uH7MrOi4j0HgUFBeyxxx6ZDkNERCQl0jXIQgmeiIiISIaka5oUyRKRSIRIJBVX1KUnqL2yi9oru6i9sovaKz4leH1cVVUVVVVVmQ5DEqT2yi5qr+yi9souaq/4lOCJiIiI5JgemYbEzMqBU4D9ge2BAXjToGwA5gNvAA855zb2RDwiIiIiuSytCZ6ZBYGrgR8CHU1atj9wDvA3M7sRuNw515TOuERERERyWdou0ZpZKfA68AugCG8ul3iPIuBXwKzosSKdqq6u5tJLL+WZZ57JdCgiIiK9Rjp78O4H9oxuO+Bp4BHgA2AtXlI3ENgdOBk4NvrensB9wHFpjE1ygHOOb3/72zz00ENMnjyZo48+OtMhiYiI9ArpWot2KnA0XmK3BjjZOfd6O0WXAR8Cd5vZfngJ4DDgaDOb6px7Ph3xSW7461//ykMPPQTAxRdfnOFoREREeo90XaKdFn0OAUd1kNy14px7EzgmekzLOkS28sorr3DJJZcAcOGFF3LuuedmOCIREZHeI10J3v54vXfPOuc+SPSgaNmn8S7VHpCm2CTLLVu2jNNOO41wOMzee+/NDTfckOmQREREepV0JXjDos/vd+HY5oRwaIpikRzS2NjIKaecwurVqxk0aBAPP/ww+fn5mQ5LRESkV0lXgte8doi/C8c2H6P1R2QrP//5z3nzzTfx+Xw88MADjBo1KtMhiYiI9DrpSvCWR5/368KxzcesSFEskiPuuecebrrpJgCuueYaDjvssAxHJCIi0julK8F7Fe8+ukPNbEqiB5nZIcBhePfvvZqWyCQrzZ07lwsuuACAb3zjG/zqV7/KcEQiIiK9V7oSvBnRZwMeN7NvdHZAtMxj0WMApqcnNMk2GzZs4KSTTqKuro7x48czY8YMzKzzA0VERPqotMyD55x708z+CZwNlAIPm9kc4FFgDrAuWrR5ouNvAHvgJXcOuCc6bYr0cZFIhHPOOYeFCxdSXFzMo48+SllZWabDEhER6dXSuZLFd4HBwJHR17tHHx1p7pJ5Djg/bVFJVrnmmmt46qmnALjjjjvYaaedMhyRiIhI75e2tWidc43OuaOBn+ANmOhsLdoVwI+dc8c455rSFZdkj+eff57LL78cgJ/+9KecdtppGY5IREQkO6SzBw8A59zfzex/ga/jTYA8DiiP7t4AzAdeB15yzoXarSTLmNkMIJmlFS52zt2UpnDiCgaDmThtpxYtWsSZZ56Jc46DDjqI66+/PtMh9Qq9tb2kfWqv7KL2yi5qr/jSnuABRBO356IP6UWKi4szHcJW6urqOPnkk1m/fj3Dhw/nwQcf1C9yVG9sL+mY2iu7qL2yi9orvh5J8Pq4a/B6KuOZ1ROBZAPnHD/84Q95//33CQQCPPTQQwwbNqzzA0VERCRGCV763eacW5TpILLF7bffzvTp3gw5f/3rXzngAC1JLCIikqy0JHhmFgCuxVt27C3n3IMJHncasA/QCPzGOefSEZ/0Tu+++y4XXXQRAGeccUZsW0RERJKTrh6844Gf481p93gSx63CG3Xr8C5bPpXqwKR3Wrt2LSeffDKNjY3svPPO3HbbbZrMWEREpIvSNU3KcdHnpc65mYke5Jx7BVgafXlCqoOS3ikcDnPGGWewdOlSysrKePTRR3XzrIiISDekK8HbC68XLuHkroXmdWz3TmlEmXOzmX1lZvVmVmVm883sPjM71czSNg9hNrnssst48cUXAbj77rsZN25chiMSERHJbum6RLtN9Hl+F45d0KaObHd0i+18vKXbtgPOAOaZ2enOuQ8zElkv8Pjjj3PttdcCcOmll3LCCeq4FRER6a50JXh50eeurEjRGH0uSlEsmbQQeAMvaW0AhgAHAZOj+3cE3jSzKc65dxOt1Mwq4+we5pxj48aNXYs4KhAIUFJSErdMY2MjtbW1XT7H/PnzOeeccwA4/PDDueqqqzosW1NTQ1NT9xc4KSsrw+eL33Ha3e8Oeub7a1ZYWEh+fn7cMvr+Oqbvr3v0/XWPvr/uCQaDnd7Sk4vfX0L1dLuG9q0DhgIVXTi2+ZiNKYum5z0D/N059357O81sCnAP3mctAh4zs/HOue7/tCepeaByTw9oqKmp4eyzz2bz5s2MGjWK++67D7/f36MxZKNMtZd0jdoru2jihuzinKOxsVH3bHfA0vEDbWYzgQOB+c658Uke+wXeJczZzrl9Uh5cL2Fm44F3gbLoW790zv1PCuqtrKioqKisjNfJt0XzXxvl5eXdPXXCnHNMmzaN+++/n7y8PF5//XX23HPPHjt/NstEe0nXqb2yi9oru/SR9uryX4fpusn/xejzdmb27UQPMrPvANtHX/4n5VH1Is65L4C/tXjr5EzF0tNuvPFG7r//fgBuvvlmJXciIiIplq4E7w68e87AG0V6emcHRMvcHH3ZBNyepth6k0dbbO9lZjl/jXLWrFn8/Oc/B+A73/kO559/foYjEhERyT1pSfCcc8vwVrIwvAEX95rZq2Z2vpntZmYjzGx4dPt8M3sVuDda1gH/00eW91rQYtsPDMhUID1hxYoVnHLKKYRCISZPnsxNN92U6ZBERERyUtrWonXOXWVm2wLn4CVtB0YfHWm+znyvc+536Yqrl8m5O3qdc6z4cgVVqzfhDwYYtv0wSgeW0tTUxGmnncbKlSsZMGAAjzzyCAUFBZkOV0REJCelLcEDcM6dZ2bvAr8HBnVSfC1whXPuf9MZUy+zbYvtMLA+U4F0V2N9I7Pue40n//xvVi9aQzA/gHOOUGOIScdM5iP/J7z22muYGffffz/bbJMr0xyKiIj0PmlN8ACcczeb2R14gwimABOAgdHd64B5wCvAo865+nTH08uc1GL7PedcOGORdMOm1Zv4w9evYPWi1TTWetMYhptCsf2PPP4Ib4bfBuCqq67iiCOOyEicIiIifUXaEzyAaOJ2b/QhgJltB/y0xVuPZSqW7qivruOqQ3/PygUrCTdtnZ9WuSreicwGoMJfwYTkZs0RERGRLtBaqClmZldGHyPjlPka8F+gPPrWKiArRxw8d/NzrF60ut3kDmC5W0GYMCUUs7ebzGN/fJR1let6OEoREZG+pUd68PqYfsCPgcvM7H3gfWA5UI+3VNmBwF4tytcD33DOVfd0oN0VCYd59v89TVN9x0vAjLUxAIy2UeRZHr6AjxdvfYHTrjqjh6IUERHpezKa4JnZgcAFwA5AHfAa8P+cc2szGVeKGN6as5PjlPkEmOacm9szIaXWvNc+o666Lm6ZfMtngu0Qe91U38RLt72oBE9ERCSN0pLgmdk44KXoy2ucc/9op8wFwP/SehmOg4DvmtmhzrlP0xFbD7gaeBXYH9gTGIE3grgM2AysAN7Cu+fuaZfFix+uX7au04WX27N53Wacc1qfU0REJE3S1YN3JDASiNDO4AEzGwP8nfbvARwCPGJmuzjnQu3s79WivY+PkaWDJpJhPiVoIiIivVG6Bll8Pfr8vnNuVTv7f8yWVSseAfaJPl6O7h8PnJam2CRFho4dSjiU/Mwu/Yf3V++diIhIGqUrwRuDl7y918H+U6L7vwJOd86965x7F29euObBBid1cKz0EtvvM45+Q8uTOiavMI+pFx2VnoBEREQESF+CNyT6vLztDjMbj3dfGsA/W07u65yrAp7Euy9v9zTFJiliZhz/ixPIK8xL+JhIOMIh3zo0jVGJiIhIuhK8AdHn9lam2L/F9vPt7P8q+jw0pRFJu4qKiigqKury8YeefxhjJ21LMD/YadlgfpDv3PxdygaVdfl8fV1320t6ltoru6i9sovaK750JXjNid3AdvYdFH2uA2a3s7/5Em3i3ULSZXl5eeTldf2rDgQD/ObpS9nhgAleT147t9YF8gP4gwHOveFb6r3rpu62l/QstVd2UXtlF7VXfOkaRbsEmEjrCX0xMx9wNN79d291MEq2f/R5c5pikxQrKCnk0ud+x8cvfcSTf3mCj1/+CJ/fRyTsyC/O47DvHs7U709lyFh1yoqIiPSEdCV4bwA7AQeb2UHOudei71+Ad+nVAf/p4NiJ0eclaYpN0sDn87Hr4bux6+G7UV9dR9XazQTzg5QOKiUQ1IIpIiIiPSld/+e9E/gu3gW7F8zseaAQOCy6vxH4Z9uDzMwP7IeXAH6cptgkzQpKCikoKcx0GCIiIn1WWu7Bc869DdyMl+DlA8fhzY3XfL7rnHNbjbDFSwCbB2jMSkdsIiIiIrkuXYMscM5dDPwSqMRL9AxYCfzCOXdlB4dd3GL7mXTFJiIiIpLLrCeWQjWzgYDPObemk3IH4SWdTc65N9IeWA4ys8qKioqKysrKTIciIiIi3dPlZZ965O5359y6BMu91nkpEREREYknbZdoJTtUV1dTXV3deUHpFdRe2UXtlV3UXtlF7RWf5q/o40Kh9qYilN5K7ZVd1F7ZRe2VXdRe8akHT0RERCTHqAdPerWlS5dy8cXe4OpgMMi999671dI0P//5z1mwYEGXz/G3v/2NsWPHMnPmTP76178C8Oc//5lx48YlVU9DQwMfffQRs2fP5r333mPNGm9M0dlnn83JJ5/c5fhERESSpQRPerW77rqLf//737HXp59++lbJ0quvvsp7773X5XNceaU3a8+yZcti5/rtb3+bVB2XXXYZ119/PU1NTVvt23fffbscm4iISFcowZNeKxKJMH36dADMDOcc06dP3yrBu/jii1m5cuVWx2/evJk//vGPABxxxBEceuih7Z5n+PDh3Y515cqVseSuf//+7LLLLsycObPb9YqIiHSFEjzptV577TUWLlwIwO9+9zv+8Ic/8Oyzz7J8+XJGjBgRK3fuuee2e/yqVatiCd7BBx/MJZdckrZYTzjhBKZOncrkyZMZM2YMmzZton///mk7n4iISDwaZCG91p133gnApEmTuPTSS+nXrx+RSIR//nOrZYwz7thjj+Wb3/wmY8eOxazL81KKiIikhBI86ZWqqqp46KGHAK+HrqCggNNOOw3wEr+eWIFFREQkWynBk17pwQcfpK6ujkAgwBlnnAFsuRT7xRdf8Oabb2YyPBERkV5NCZ70Ss2XZ48++mgGDx4MwH777cf222/far+IiIhsTQme9Drz5s2L9dCdc845sffNLPb6X//6FzU1NRmJT0REpLdL+yhaMysADgV2B4YAJSSWWDrn3HfSGJr0UjNmzAC86UaOPfbYVvvOPvtsLr/8cqqrq3n44Yc7HEErIiLSl6UtwTOzAHAp8FOgrIvVKMHrY5qamrjrrrsAb1Lj/Pz8VvvHjBnDwQcfzKuvvsqdd96pBE9ERKQdablEa948EQ8Dvwf6AdaFh/RBzz33HKtWrQJaX55tqTmpmzlzJvPnz++x2ERERLJFunrwzgGOj25HgH8BjwBfANWA5riQdrUcPPGzn/2s3Tnl6urqYtszZszg6quv7pHYREREskW6Erzm62YR4CTn3BNpOo/kkNWrV/PUU0/FXicyFcpdd93FlVdeid/vT2doIiIiWSVdCd6ueL10zyq5k0Tdc889hEIhAG655RZKSko6LLto0SIuu+wyKisrefHFF5k6dWpPhSkiItLrpSvBK44+v52m+iVFysvLMx0CAM652OXZ/fbbjwsvvDBu+XA4zC233MLy5cuZPn16n0nwekt7SWLUXtlF7ZVd1F7xpWsevGXRZ103k4TMnj2bTz75BIBp06Z1Wt7v98dWuHjsscdYv359WuMTERHJJulK8GbhjYTdNU31S45p7r3z+/2ceuqpCR1z5plnAtDY2Mh9992X0nimTJlCWVlZ3Mcdd9wRK//GG2+02jdq1KjYvt///vet9p1//vkpjVVERKStdCV4fwfCwNFmNiZN55AcUVtbG0vQpk6dGluarDN77LEHEyZMAGD69Okpj2nz5s1xH01NTbHy4XC41b7q6urYvsbGxlb76uvrUxqriIhIW2m5B885976Z/Ra4Dvi3mR3lnFuejnNJ9nPO8dJLLwEwYsSIhI8zM5555hnWrVsHQCQSwefb8jfLgAEDePfddwGoqKjotL4jjjgiVj4Ro0ePjm3vvvvuCR87YMCAhM8hIiLSFeZc6qekM7Pm//NdgLeaxSbgNuAFYAmQUBeGc25JyoPLcWZWWVFRUVFZWZnpUERERKR7urzwQ7oSvAitJzM2kp/c2Dnn0r5Wbq5RgiciIpIzupzgpTOBahuUlh+TmNqNNcz65yzeefAtajbUkF+Uz46H7MihFx7O4LGJ3YMnIiIi7UtXgndXmuqVFItEIgCt7l1L9/ke/t2DvHjjC/gCPhprG2P7ln+2jBdv+g8Tv74z37vrQorKi+PU1Df1dHtJ96i9sovaK7uoveJLyyVayZxkL9Fu3LgR6JkJIyORCP931v8y55kPaKpr6rBcID/AgFEDufz1K5TktdGT7SXdp/bKLmqv7NJH2qvLVz+V9kqP+c+Nz3ea3AGEGkKsX7qO/zvnlh6KTEREJLcowZMeEQlHeOZPT3Wa3DULNYT49KVPWL1wdZojExERyT1K8KRHfPjcXBpqG5I6xvzGy/94MU0RiYiI5C4leNIj5r/5ZcK9d81CDSE+ffmTNEUkIiKSu9IyitbM7uzG4Q6oA6qAxcAHwLtOo0GyWv3mOrrShI1J9vqJiIhI+qZJOY/kJzaOZ4WZ3QL8yTmXXDeQ9AolA0vxBXxEQpGkjtMoWhERkeSl8xKttXi0fd3eI16ZEcBVwLtmNjCNMUua7HrUbvj8yf245RXmsfep+6YpIhERkdyVrgRvbPTxV7yePAc8ibc27f7AhOhj/+h7T0bLRIC/AOOAScCZwN1AA16ityvwaJpiljQau+e2DBqT3AoV4VCYA885KE0RiYiI5K60JHjOucXARcBPgSXA3s65E5xztzvn3nLOfRF9vBV97wRgb6AS+BnwfefcHOfcA86584DJePfjARxoZsenI25JHzPjlD+eRiA/sbsCgoVBDr3w65QMKElzZCIiIrknLQmemR0C/Byv522qc+79zo6JlpkKNAI/NbNDW+z7FPgmW+7rOzXlQUva7XHcJL559amdJnnBwjx2PnxXTrv+jB6KTEREJLek6xLt96PPTzrnvkj0IOfc58C/8S7H/qDNvveA16L79klRnNLDjvjRkXz/3osYOm4Ygfxgq/vy8oryKBlYwomXfYMf/uvipO/ZExEREU+6RtHujdfb9lEXjv0Yr4dur3b2vQV8DRja9dAk0/Y4bhK7H7sHC99dyNxnPmDzms0U9itk+33Gsdsxu+MP+DMdYq8VDAYzHYKIiGSBdCV4zQlYQReOzY8+D2ln34boc7rilh5iZmy393Zst/d2mQ4lqxQXa9oYERHpXLoSpfXAMGBKF449JPq8oZ19pdHnNV2oVwSA1atXE4lsPR9fYWEhpaWl+HxbXxru6JhEDR48GL/fT01NDZs3b95qv9/vp6ysjPz8/HaO3mLt2rWEQiEACgoKKC8vj1t+5cqVAJSUlFBSEn/AinOOtWvXEg6HMTOGDk2so7ypqYnPPvuML7/8MvY9DRo0iIMPPjjhOkREJLXSleC9DZwI7GdmZzjn7k/kIDM7E2/qFBeto62dos+rUhGk9E277rorq1a1/yNUVlbG1KlTOeGEEzj99NPx+73LxbvvvjsrVqzo8jlXrVrFkCFD+Mc//sEvfvGLDstNmjSJ4447jnPPPZexY8dutf/QQw/lo4+8Ox+GDx/O/PnzKSoq6rC+4cOHA/DHP/6RSy+9tNW+hoYGnn32WWbPns17773H7NmzWbt2LQClpaVUVVV1+rmOOeYYXn75Zerr69vdv+eee3LllVdy9NFHd1qXiIikTrruYr+1xfZdZnaZmXV4bcnMis3scmBGi7f/r02ZArzePQfMTmGsIjFVVVU89NBDnHXWWRxxxBHU1dX16Pnff/99rrzySnbZZReef/75uGVXrFjBjTfe2OVzLV++nG984xv88Y9/5Lnnnosld8l47bXXqK+vJy8vj4kTJ3LIIYewzz77MHiwN+fh7NmzOeaYY5g+fXqX4xQRkeSlax6854DpeCNe/cAVwCoze9rMrjezX5vZJdHtp/F65H7Plh7Fu6J1tHQ2UBbdfiYdcfdFwWCwz964f/rppxMOhwmHw4RCIdavX8+DDz7IbrvtBsDLL7/MNddcA8C8efPYsGHDVo+ZM2fG6rvpppvaLbNhwwYGDRq01fmXLFkSO399fT1z5szhJz/5CWZGTU0Np59+eqcJ5nXXXceGDe3dzZAYv9/Prrvuyre+9S1uuukmjjnmmKSO/9WvfsWsWbOoqqrik08+4eWXX+att95i5cqVPPbYY7FE75JLLqG2trbLcWarvvz7lY3UXtlF7RVfOgcrnA/U4U2ZYkARcGT00Za12P4H3iTJba1iy/QrL6YuzL6tL9+0b2at7rfr378/p5xyCocffjjjxo1j7dq1TJ8+nSuuuIJ+/fq1W0dpaWlsu6ioqNN74lry+Xyx8+fn57Pbbrvxt7/9jaKiIq655ho2btzII488wllnnbXVsfvvvz/vvfceGzdu5M9//nMsEU3GyJEj2bx5M4WFhbH3Pvnkk6Tq+N3vftfu+z6fjxNPPJH169fzne98hzVr1vDoo4+2+1lyWV/+/cpGaq/sovaKL20TjTnPRXgDLZ4BwnS81mw4WmaKc+4Hzrmt7mZ3zj3hnPu/6KP9G35EUqC8vJwzzvAmWV62bBkffPBBj57/+9//fmz7iSeeaLfMyJEjuegi7++gG264oUv3BwaDwVbJXTq0vPdu+fLlaT2XiIhskfaZZJ1zrznnjgUGAYfj9cL9Brg0un04MMg5d6xzbmbHNYn0nHHjxsW2169f36PnHjlyZGw0bbxz/+Y3v6GsrIy6ujr+8Ic/9FR4SamsrIxtb7edpsQREekpPbZUgHNuk3PupWgP3PXOueui2y855zb1VBwiiVi6dGlsu6ysLE7J1FuzZg0NDQ2dnnvgwIH88pe/BOC2225jwYIFPRJfopqamrjiiisAGDRoEEcddVRmAxIR6UM0YbBIG/X19dx/vzezz8CBA9ljjz169PwzZsyIbR9xxBFxy/7kJz/hxhtvZPXq1Vx++eXce++9aY6ufatWreLll18GvJHIixYt4t5772Xp0qUMHTqUf//733GncxERkdRSgicS5Zzj7bff5rLLLotdWpw2bVqnkw+nypo1a/jnP/8ZG7hQUFAQuxewIyUlJfzud7/jRz/6Effddx+//OUv2X333Xsg2tbmzZvHmWeeudX7J510ErfddhsDBgzo8ZhERPoyJXiSFdYvWc9Xby+koaaBwn6FjPvaeEoGxl+ZoTOPPfYYFRUVgJfcVVdXt1plYqedduLyyy/v1jnimTx5cmwi5aamJtatWxdbLcPn8/GPf/yjw9G7LV1wwQX89a9/ZdGiRfz2t7/l6aefTlvMHenXrx9f+9rXAKirq+PTTz+lpqaGRx99lKVLl3LrrbdmJPEUEemrupXgmdnCFi+dc267dt7vqlh90nd9OetLnrv2ab6Y+QXBgiDOgRmEGkLsdsLuHHnJ0VTsXNGluuvr67ca2Wlm7LPPPhx//PH86Ec/Susw/PZW0xgyZAjHHHMM3/ve99hnn30Sqic/P58rr7ySc889l2eeeYbXXnuNgw46KNXhxrXHHnvw6quvxl43NDTwwgsvcNFFF/Huu+9yyCGHMGfOHLbZZpsejUtEpK/qbg/eGLyVJSz63N77XeU6LyLd1Xwzf09dhkzGyze+yL8ve5xQYwgcNFQ3tNr/wWPv8+GTc/n2P89nt+N2T7r+I444gquvvjr2uqioiGHDhjFw4MDuhp6Qp556iiFDhgDehMP9+vVj7Nix7a6F25lp06bxpz/9iU8++YTf/OY3vPbaa5h159eve/Lz8znuuOMYO3Ysu+yyCxs3buSmm27iz3/+c8ZiyoTe/PslW1N7ZRe1V3zdTfCW0H4i1tH70ss0r5TQ235B3rn/bS+5awh1WMaFHaFwiDvPvp0fPfMTttt/+6TOMXDgQPbaa6/uhtplu+++e+wScXf5/X6uueYaTjjhBF5//XWefvppjj322JTU3R0777wzU6ZM4ZVXXuGpp57qcwleb/39kvapvbKL2iu+biV4zrkxybwvkohwU5iHfv6vuMldS6GGEA/9/F/8+s3fpjmy3u24445jv/3248033+TSSy9tNclwJg0bNgyATZs0G5KISE/psXnwRBL14VNzE07umi3/ZDnLPl4GeN324XA4HaH1ambGtddeC8BHH30Um+olk8LhMLNnzwa8CZxFRKRnKMGTXufNu9+gsbax03IRF+Gdpnd5qvEZ/lXzMJMOmkRpaSkFBQWMHTuWjRs3pj/YXubggw/myCO95Z4vu+yytJ4rkR65f/7zn8yfPx+A8847L63xiIjIFpomRXqdTSs2JlSuylUxO/z+ljdaHFZVVRWbcqSvueaaa3juuef46quvOi374osvsmbNmtjrL7/8EvCmbWnbA3j88ce3GlX8P//zP8ycOZMzzjiDAw44gG233ZaioiKqqqp44403ePzxx7n99tsBbyWLc889NxUfT0REEpDRBM/MtgXOA3YA6oBZwD+dcw3xjpPc5s9L7Meyn/Xj4MBBVLtqCqyAXafsxvE/P4FBgwYxfvz4Hl9irLfYY489OO200/jXv/7Vadkrr7ySWbNmbfV+fX39VhMXL168eKtpY2bOnMnMmfGXkB4+fDhPPPFEWqecERGR1tKS4JnZSOC+6MvbnXN3t1PmOOBfQMvhL2cDPzWzw51zy9seI33D6EnbUDlnKeGm+PfRmRk7BSYCkFeUxzdOOZmvHXlwp/XvuuuurFu3jrFjx3Y71sLCQiZNmgSQ0PQqQ4YMiZUPBoNJn2/HHXckGAx2Gvsf/vCHWG8cbBno0Nb48eOpra1N6Nxt4z3zzDOpr6/nv//9L3PnziUU2nLfZGFhIVOmTGHq1Kmcd955CU3YLCIiqWPOpX42EzM7H7gVb6qU7Zxzi9rsHwzMB0o7qOJt59x+KQ+sDzCzyoqKiormpbY603yfWnl5efqCStLKz1Zwzd5Xd5rgtRTID3Ddkj9TWFaYxsikI6FQiBUrVlBTU0NxcTFDhgzR1AX0zt8v6ZjaK7v0kfbq8oSm6Rpk0bxC+qdtk7uoi/GSOwe8CZwWfXwY3b+3mZ2Qptiklxs2YTijJ43GF0jsxzOQH2Cv0/dRcpdBgUCAUaNGMWHCBEaNGqXkTkQkw9KV4I3DS97e6WB/8wrqq4CvO+cecs49BBwJNA+fPDVNsUkW+Pbd36WgrBDzx//jxZ/nZ8DoAXzzz6f0UGSZtXHjxj45OlhERJKTrgRvcPR5SdsdZrYNsB1eAniPc66ueZ9zbiXwDF6X5OQ0xSZZYMDoAVwy69cMHD2QYGHeVp3U5jMCBQFG7zGaX7xyCQWlBZkJVEREpBdK1yjaQdHn9u7e3r/F9jPt7P88+jwipRFJ1hk0djC//+gqPv3PJ7x0w3/46u2vCDWGCBYGmXj4Thz2o68zdt9tM7rmqoiISG+UrgSvCQgC5e3sO7BFmbfa2b8x+qwuGcHn97Hzkbuw85G7ZDoUERGRrJGuS7TNU5zs1s6+o/Euz852ztW3s795PoXE5m4QERERkVbSleC9g3fX1OFmNqH5TTP7BrBN9OVLHRy7Q/Q5sXk+RERERKSVdF2i/ScwLVr/a2Z2L1CIN5ExQAS4p4Nj98Xr4fskTbGJiIiI5LS0JHjOuRfM7DHgG8AAvHnvYMtYyFudc1+2Pc7MDgCG4yV4r6cjNhEREZFcl861aM8A/gKcz5blyJqA/wN+3sEx32+x/Xz6QpNmRUVFmQ5BkqD2yi5qr+yi9souaq/40rJUWasTmBUDE/Du9/vMObc5Ttkz8UbfNjnn7uuonHTMzCpHDB9eMfc//0n4GH9+Pv23336r9zfMn0+4oSHhegoHD6Z4yJBW74Xq69m4YEHCdQCUb7cdgYLWg6hrVq+mbs2ahOvQZ+qYPlPH9Jnap8/UMX2mjukztS+ZzzRop526PA9YOnvwAHDO1QDvJVhWSV0KNFVXM2PnnRMuP3DiRL71yda3PD5+wgms+/TThOvZ7/e/54Arrmj13sYFC5KKBeC8jz9m0E47tXpvzv/+L29eeWXCdegzdUyfqWP6TO3TZ+qYPlPH9Jnal8xn+kU3OuHSNYpWRERERDIk7T14zcysDNgdbxmzIrx58Ob11PlFRERE+oqeuAfvKOAS4ABa9xj+0jn31zZl7wcGAp86536S1sBylO7B8+gzdUyfqWP6TO3TZ+qYPlPH9Jna11P34KUtwTMzH/AP4DvNb7XY7Wg/wbsMuBJvnrwxzjlNdpwkM6usqKioqKzUVyciIpLlupzgpfMevL/jTZFiQBXe9Cg/7OSYGXjJnwHHpzE2ERERkZyVlgTPzPYEfoCXrL0NjHfOfd85d0u845xzS4EPoy8PSUdsIiIiIrkuXT14F0Sf64BvOOcSvzgNs/F68HbqrKB0X3V1NdXV1ZkOQxKk9souaq/sovbKLmqv+NI1inYKXu/dv51zK5M8dln0eURKI5J2hUKhTIcgSVB7ZRe1V3ZRe2UXtVd86erBGx59TnwWwi3qos9ag0RERESkC9KV4HV51AfeNCkAm1IRiIiIiEhfk64Eb1X0ebsuHLtP9HlZ3FIiIiIi0q50JXhv4fXiHWVmwUQPMrOdgAPx7t+bmabYRERERHJauhK8h6PPQ4DfJ3KAmfUD7mXL5d0H0hCXiIiISM5LS4LnnHsMeA8vWfuNmf0tmsBtxcx8ZvYNvOlRdsXrvfuPc+6NdMQmIiIikuvSuZLFqcBqvCTvR8BKM5vVYv9pZvYssAavx2/b6PvLgHPSGJeIiIhITktbguec+wo4CJiDl+TlA/vh9dAB7AkcAfSP7rdo2YOcc6vTFZeIiIhIrktnDx7OuS+BvYHzgDeAEFuSueZHBO9y7neBvZ1zi9MZk4iIiEiuS9dKFjHOuRBwN3C3mRUCE/DmugsA64AvnXMb0x2HiIiISF+R9gSvJedcHfBBT55TREREpK9J6yVaEREREel5PdaDZ2YlQClQ5Zyr6anzioiIiPQ1aUnwzCwPOA04FG8k7aiW5zKzJmAp3moVLwMPOuea0hGLiIiISF+T8ku0ZvZ9YCEwA28+u7FAkNYjZ/Pw5r07D28AxkIzuyDVsYiIiIj0RSlL8MwsaGZ3ATcBw5vfZsvSY1sd0mJfBXCLmd1hZj068ENEREQk16QymfobcDbeRMYGNAFPAq8Bc4ENQDXefXj9gd3wLt8ey5YevvOAzcBPUhiXxFFeXp7pECQJaq/sovbKLmqv7KL2is+cc52X6qwSs0OB/zS/BO4FfuWcW5HAscOBPwNnRt+KAIc652Z2O7A+yMwqKyoqKiorKzMdioiIiHRPR1dBO5WqS7S/bhHEVc65sxNJ7gCccyucc2cBf4i+ZcClKYpLREREpM/pdoJnZmOBw/Auzb7inLuiK/U4534PvIKX4H3dzLbpbmwiIiIifVEqevCmsqX37o/drOvq6LMBR3azrl7BzPY3s9vN7AszqzGzDWY218z+aGZjMh2fiIiI5J5UDLKYHH1e7Zx7uTsVOedeNrNVwJAW9WYlM/MDNwA/pPU19CKgHNgV+JGZXeScu6vHAxQREZGclYoEbze8y7Nvp6AuovUcD+yeovoy5Sbgwui2A17A+2wFwFHALkAJMN3MapxzD2ckShEREck5qbhEOyj6vCAFdQHMjz4PTlF9Pc7MjmRLctcAHOecO9I593vn3CXOuV1pPajkdjMbmIlYRUREJPekIsHrF33emIK6ADa1qTcbtbwX8Qrn3NNtCzjnLgeej77sB/yqJwJrKxKJEIlEMnFq6QK1V3ZRe2UXtVd2UXvFl4pLtGXR55oU1AVQG30uSVF9PcrMJgCToi83412q7cg1eINUAKaZ2a9dKiYmTEJVVRWgCSOzSV1dHcXFxRk7f1NTE6tXr2bdunXk5+czbNgw+vVL7O+xxsZGKisr2bRpE3l5eYwaNYqysrLOD8xS+v3KLmqv7KL2ii8VCZ4f7x6zVPOnoc6ecHyL7eedc9UdFXTOzTSztXiXuSvwBpbMTnN80o4HHniA3/72twCMGzeO5557rtNjHnnkEX71q607Xn0+H6Wlpey+++4ceeSRTJ06NW4C9Oijj/LLX/4ybj1Tp05l6tSplJeXd5jcPf744/z85z/f6v2CggLKysoYOnQohx12GEceeSTjxo3r9PO1NWfOHO644w7uueceNm7c2Grf+PHjOfHEE7ngggvYbrvtWu1bsGABf/rTn3jjjTeYN28e4XA4ti8QCHDggQfyzW9+kwsuuIBgMJh0XCIi0g7nXLceeCtPhIGfdbeuaH0/b64zFfX19AO4Hy/hdcBPEyj/RIvy303B+SsrKipcojZs2OA2bNiQcPlc9bWvfc21aAf3wQcfdHrMjBkzWh3T0WPw4MHupZde6rCeu+66K6F6Bg0aFLeef/7znwnVA7jtttvO/eY3v3E1NTWdfs6vvvrK7bvvvlvVUV5e7oqKilq9FwgE3HXXXdfq+Mcff3yrYwsKCrZ676CDDnKbN2/uNJ5sot+v7KL2yi59pL26nA+kaiULgP5mNrq7D7x1arPZji2253dYqv0yO3ZYStJm/vz5zJzprYxXVFQEwPTp05Oq484772TmzJnMnDmT//73v8yYMYMTTjgBgDVr1nDCCSds1evVnjvuuGOrek488UQA1q5dywknnMCGDRs6ree2226L1fPyyy/z+OOPc80117DvvvtiZixYsIBrr72WPfbYgyVLlnRYz/z589lnn3146623ADjwwAO5++67qampYcOGDdTU1LB27Vruu+8+DjzwQEKh0Fa9nyUlJUybNo0ZM2bw6aef0tDQQF1dHVVVVbz88stMmTIFgNdee42///3vnX42ERHpXLfXojWzCKm/RGuAc85l3WVaM1sJDI2+3NM5914n5X8O/E/05f3OuTPjlY8eE2+h2WEjRozwf/LJJwnF29z+Zq2XuwsEApSUxL8NsrGxkdra2rhlElFYWEh+fn7cMjU1NTQ1NXX7XGVlZfh8rf+u+d3vfscf//hHRowYwYUXXsjll1/OgAEDWL58Ofn5+Vuduzneu+66i/POOw/wLl/utttuW53viiuu4MorrwTg5ptv5gc/+EFsXzgcxu/3c/fdd3PuuecC8MEHH7D77rtvVc+VV17JFVdcAcCNN97IRRddRCgUorq6mtLSUvx+P/fccw9nn302AO+99x6TJk3aqh6AlStX8utf/5q77vKmX9xpp5145513KCoqoqbGu5W2uLiYhoYG9t13X+bMmYPf72f69Omx+tvjnOPhhx/mnnvu4d///neH5doKhULsvPPOfP755wwaNIiVK1fi93f+q58NP39tf7/a+/lrK5E/BDrTl35/2+rO99fcXsFgUN9fN/TUz19ze/XvH79fKMu/v4yvRdvMUvDIdi1/qhMZeNLyJ7w0xbFIJ8LhMDNmzADgrLPO4rzzzsPMWL9+PU888QQAeXl5CdfX2NhIXV1d7D6ziy++mEDAu9X1scceS7qeUCgUq6f5/rTmehJJgpo1NTXF6ho2bBgzZszgu9/9LgCffPIJjz/+OOB91ubzPPnkk8yZMwfwkuB4yR14Scwpp5zCP//5T8Ab4dbY2Ah4/xA3/4NeXV3N5s2bqampwTlHIBDgrLPOArxeyrVr1yb8uUREpH2pGGQBqU3Msj3JK2yx3ZhA+YYOju2Qc25kR/vMrNLMKhIdVdT810ZXRiHl5eUllfx0R7pGjf7nP/9h2bJlAJxzzjmMGjWKQw89lJdeeonp06dzyimnEAwGE/5+2n4fAwcOZOjQoSxbtox169a12hcvQWtbz4ABAxg6dCiVlZWxesws4bjaG7zwpz/9iRkzZtDU1MQ///lPzjzzzFblHnjggVgsv/jFLxI6DxAbFevz+WKfw8zi/rw0/3WdyF/+zbLh568rv189NSIwG76/rujO99fcXon8DOr7655UfH+J9pbl4veXiG734DnnfGl6ZN3l2ai6FtuJ/PS27Juv67CUpEXzvXaTJ09mp512ArxED+D555+PJX9dVVNTE0vIhg0b1uV6amtrU1JPS+Xl5Rx22GEAvPDCC6xatSq2r6qqiqef9qZvPPTQQxNOurqipqYm1ot65JFHZnQKGBGRXJHqS7QCLadFSeT/VEUttjenOBaJY926dbFLk833wAGcdNJJFBcXE4lEuPvuu7t1jnvvvZf6+noADjjggC7Xc99991FXV9fteto66aSTAO9y6vPPPx97/8UXX4zFfeyxx6bsfOvWrWPJkiUsWrSId955h1tvvZXJkyezZMkSxowZww033JCyc4mI9GWpukQrW6xiyyCLEUDcQRZ489+1PFZ6yH333UdjYyOBQIDTTz899n5JSQnf/OY3ueuuu7jzzjv59a9/vdUglLaqq6tjk26GQiEWL17Mo48+yp///GcAhgwZ0mqARSL1hMNhFi1axGOPPcaf/vQnAAYPHpxQPYkaP358bLvlvW8rVqxot0x3/epXv+LOO+9s9V5BQQGXXnopP/nJTxg8OGtXKBQR6VWU4KXep8Cu0e3t4hVsp8ynqQ9HOtJ8efaYY47ZKrE499xzueuuu5g/fz6zZs3ioIMOilvXgQce2OG+/fffn+nTp3c60guIe5799tuPGTNmMHBg6pYtLi3dMq5n06ZN7W63LNNdPp8Pv99PJBKJjYCrr6/nnnvuYeLEiZx55pmdJtMiItI5XaJNvTkttvdPoPx+HRwrafTBBx/wwQcfAFvuuWvp4IMPZvTo0UDyc+K1NGDAAC677LJu94Klqp62mqdFgdaJXMt77lIxFUSz2267jVAoRDgcZuPGjTz++OMceeSRLFmyhLPOOovf//73KTuXiEhfpgQv9Z5ssT3VzDq8D8/MDgSau46Wo2XKekxz0jZgwACOOeaYrfb7fL7YtCAPPvgg1dUdrjgHwBtvvMHmzZvZvHkzlZWVPPXUU+y5556sX7+eY489lmeffTahuF5//fUO6znuuON45plnkvyk8S1evDi2PWDAgNj2oEGDYttLly5N6TnBG1Xbr18/TjjhBJ555hm++c1vAnD11VdTWRlvmkcREUmEErwUc859CnwQfVkG/DBO8d+02L7PdXfWaUlIQ0MD9957L+ANsx8yZAj9+vXb6vHXv/4V8Hq5Hnroobh1FhUVUVJSQklJCRUVFRxzzDG89NJL7LLLLoTDYc4991zWr1/faWzt1fPyyy+z2267xeppO91Kd7SckPiQQw5ptd18qTTR5LSrzIzLLrsM8ObLS2a+QBERaZ8SvPT4XYvtK83syLYFzOwK4Ojoyyrg+h6IS4AnnngilmxFIhGqqqrafTSPWgW2GhiQiLKyMv7xj38A3nJlf/zjH7sUb2lpaayetWvXdrmetmpra2M9ggceeCBjx46N7Rs+fHhsCbHnnnsuJbPAx9M8RQ2kp8dQRKSvUYKXBs65Z4Bboy8LgGfM7Fkzu8LMrjOzOUDzzUYOuMA5p+n7e0hzsrbtttsyZ86cuI+rr74agFmzZvHFF18kfa7999+fU045BfCWKlu+fHmXYt5333057bTTYvV0d34+gL/97W+x++uaV5JoqXlk8aZNm2I9nomIRCJJx7Jw4cLYdstLxSIi0jVK8NLnB8At0W0DjsRL6i4BmhctrQG+45z7V8+Hl5tq1tXw6t//yz3n3sX0U2/n0Z8+zFdvLoyN2KysrIzN9zZt2jR22223uI/vf//7sdUdmifjTdbll18OeJeGr7vuui5/tssvvxwzo7GxkWuvvbbL9YTDYf70pz/F4ho+fDinnnrqVuVOOeWU2KTKP/7xj3nnnXc6rfuTTz7ZahqXRO48+Nvf/hbbPuGEEzotLyIi8SnBSxPnXNg59wPgQGA6sABvpYpNwEfAdcAuzrmuD9FMgWAw2O4yVtmmblMd9337n/xh+9/z/B+eZe4jc/j0mU946443uPW4W9iwZAMAd999dyzhmDZtWqf1DhgwgKOP9q6k33XXXbE1ZpOx8847x3rxbr311i73vk2cODGWiN12221xByNs2rSJ9evXs379elavXs2XX37JE088wXe+8x2GDx/OJZdcQiQSoby8nKeeeqrdKVz69+/P/fffTzAYpKqqikMOOYQf/ehHfPzxx63KOef47LPP+PGPf8ykSZP4/PPPW+0/4YQTuPPOO1mzZs1Wx33xxRf85Cc/4ZZbvL+Fjj76aHbccccufT+9Ua78fvUVaq/sovaKT/PgpZlz7nXg9UzH0ZFcWBaqZl0NNx5yAxsrNxBuChNu2pKEuYhj1KTRDNhmAM652OXZyZMns8MOOyRU/7Rp0/j3v//N8uXLeeGFFzjqqKOSjvGyyy7joYceoqGhgWuvvZabbrop6Tqa63nwwQdjvXg333xzu+UOPfTQuPXk5eVxzDHH8Pe//52RIztc2pgpU6bwwgsvcPLJJ7N+/XpuvPFGbrzxRnbccUeGDx9OY2MjixcvbnXfXNvv9cMPP+TJJ73B5aNGjWLo0KEEAgGWLl3aKtkdM2ZMh58nW+XC71dfovbKLmqv+NSDJ1nNOcedp9zmJXeN7feu7XD4BABee+01FixYACTWe9fs2GOPjc0R15XBFgC77LJLbCqQznrf4tlpp51ivYG33357QgMSmqckGT9+PN///vd54oknWLduHY8++mjc5K7ZlClTWLx4MXfeeSf77+9N7Thv3jxefvllZs2axdKlS/H5fBx11FE888wzsd64ZmeeeSbbbLMN4A2gmD17Nm+99VYsuZswYQKXXHIJc+bMYcyYMcl8HSIi0gHTzBy5xcwqKyoqKvrKXGJL3l3MLUfeRKgh1GGZo644hkN/8XVCoRChkFcuGAzi9/sTPk9DQwPOOcyM/Px8wLuXrXl0aV5eHj5f/L+XOjp/OuppKZE6k/Hll1/y+eefs3btWgoKChg+fDiTJ09uNTlye5YuXcqnn37K6tWrCYVClJWVseeee8aSPxER2UqXl/ZRgpdj+lqCd885d/Hhv+fiwh3/HBcPKubi//6UgWNTt8SXiIhID+hygqd78CSrff7SZ3GTO4CatTVcv+sfqdh9JPml+QQKguz3nf3Y6ZhdeihKERGRnqUET7JaU11iE/A656j8wLtfLVgYZMepE9MZloiISEZpkIVktWBB8kPkzecjvzQ/DdGIiIj0DkrwJKtte9B2mC+5WxRCDU1se8B2aYpIREQk85Tg9XENDQ00NDRkOowuO/jiQ/DnJT4aFoOx+2/LgG2yczmsbG+vvkbtlV3UXtlF7RWfErw+rq6ujrq6ukyH0WVjD9iWwdsPwRdM7EfZH/Rz+K+npjmq9Mn29upr1F7ZRe2VXdRe8SnBk6xmZpz/7+9RMri00548f36AY685nu2+tn3svXBTmOo11WxetZlQY8dz6YmIiGQTjaKVrFc2tIyfvv4LHvjefXz58uf4/L4tEx+bNxAjvySfE//nJHY7eQ8AVn+2irf+8QZz7/uAUGMIM8Pn97Hzybuy3w8PYPiuIzL4iURERLpHCZ7khJLBJZz/6AVsWLqBt+54g6XvLaGpvol+w/sxedpe7PD1Cfj8PpxzvHz1f3j9/72GwxFpigDeNCrhSJiPHp7Lx49+yKSz9+ToPx+Hz69ObhERyT5K8CSn9B/Vn6OuOKbD/S9d9QJv/u/rhJvaX7c2EvISvg/ufY9IKMJx/+9EzLo8kbiIiEhGqHtC+oyVH63gjZtmEarv/F67UH2IuQ98wKLXFvZAZCIiIqmlBE/6jDf/93VcJPG1l8NNYd64aVYaIxIREUkPJXjSJzTVNfHRw3Njl2AT4SKO+S99Sc3a6jRGJiIiknpK8KRPqF61GRLvvIvx5/nZtHRT6gMSERFJIyV40idEwg66OFYiEkm8109ERKQ3UIInfULJkBIi4eQTtXBjmNJhZWmISEREJH2U4EmfkF+az7jDd8B8yXXjjdijgn4V/dIUlYiISHoowZM+Y/+LDux0ObOWgoVBDvjRQWmMSEREJD2U4EmfMebAsYyfOoFAQefzewfyA4zcezQTjpnYA5GJiIiklhI86TPMjJNvO5XtDxtHoDDYYblgYZBR+4zmzPvP1lJlIiKSlbRUmfQpgfwAp90zjc+fnsfrf3+NZR9UEsgLgEGoIcSQHYdy4I+/xsQTd8YfSPxyroiISG+iBE/6HJ/Px47H7cSOx+3ExiUb2LhkA85Bv4p+DNh2YKbDExER6TYleH1cUVFRpkPIqPLR/Skf3T/TYSSsr7dXtlF7ZRe1V3ZRe8WnBK+Py8vLy3QIkgS1V3ZRe2UXtVd2UXvFpzvIRURERHKMEjwRERGRHKMET0RERCTHKMETERERyTFK8ERERERyjBI8kQQ552isaSTcGM50KCIiInFpmhSRTqz8eCVv3/oWnz72CaGGEADFQ0rY57v7sPu0PSgeVJzhCEVERFoz51ymY5AUMrPKioqKisrKyoTKV1dXA1BSUpLOsLJSuCnMkz95gk8f/wTnHJFQpNX+QEEAF3Ecf+OJ7PyNnXskJrVXdlF7ZRe1V3bpI+1lXT1QPXh9XCgUynQIvZKLOB654GEWvDSfcFP7l2RD9d5398TFj2MGO52Y/iRP7ZVd1F7ZRe2VXdRe8ekePJF2fPzYxyx4aX4siYsn3Bjm3xc9Tu362h6ITEREpHNK8ETa8cZNryeU3DUzvzH3gTmt3nMR3f4gIiKZoUu0Im2s/mw1675Ym9QxoboQb//fW4zebxvevuNdvnjpS0J1IXwBH0MnDmG/C/ZlhyPG48/zpylqERGRLdSDJ9LGuvnr8Ocnn4htXrGZu069h3nPfEaozuv9i4QirPhwJU/+8mn+vv/NrP5sdarDFRER2YoSPJE2IuFI54U6EG4Mt3tptqmuiZp1NUw/6W7WfJlc76CIiEiylOCJtFE6tKTLSZ5ZnBHtzht5++B3H0bTE4mISDopwRNpY+Reo8grzkv+QH/nv04u4ti8ooqlsxObp1BERKQrlOCJtOHz+9jngn0J5Cc5BimY2H17oaYw79zxbhciExERSYwSPJF27PntvSgZWoIvkOCviN+H+RKccDwCqzTYQkRE0kgJnkg7CsoKOPeJb1FWUUagoOOePF/AB36DJKc/iTR1fSCHiIhIZ5TgiXSgbEQZF7x8IQf/cgrFg4sJFATIK8kjryQPX9BHxeQKdjl1NwKlBfEHV7SjeHBxmqIWERHRRMciceWX5rP/xQew3w/3Z9l7lVSvrsaf52fAdgMZuO1ANi2v4uOn5iVVZ7AwyB5n7J6egEVERFCCJ5IQ8xkj9xq11fv9RpQxZv9tWPjaV5DgzCfOOXY+fmKKIxQREdlCl2hFuunrvz0s4RG3gfwAh116KMHCIABVKzbz6TOfM+ehj5j33BfUrKtNZ6giItJHqAdPpJuG7DCYM2acxgPfepBwU5hIqJ0BFAb+oJ8DfrAfe50zmaXvLWPWLW+x+O2lBPL8OK8IocYw4w7djgO+tw/DJg6JHe6co3LOChbNXUq4McyA4f3Z/sAxFPUvbDemmvW1bFpZjYs4SgcXUza0JD0fXkREeiXTjPq5xcwqKyoqKiorE5tId+PGjQCUl5enL6g+Yv3iDbxxy5t89NjH+AN+iKZtoYYQo/YayQE/3J9tDxzL7Hs/4KU/zSTcFG73sq75DF/Ax/F/OooJh49jzmOf8OaM92K9e845fH4f4aYIOxyyLQdduA8Dx/THOceid5fx9v0fsuT9Zfij8/KFm8IM22Ewe5+5K+MPGoOvzYTMkXCEhpomfAEfeYWBpAeMSMf0+5Vd1F7ZpY+0V5f/QVaCl2OSTfAk9RqqG1j81hLqN9UTLAwyfJdhlI8qB2Dec1/wxK+e9ZK7TviDfkZMGsHKL9YSami/vPmMQJ6fU244lo+e/5LPXl5IqIPEMZDvp2KXoZx87VTyCoOs/GIt7z76KZ/NXBTrdcwvyWPS8RPY49gdKB2kkb4iIhmmBE88SvB6r0g4wt8P+j9qN9QlfpDPsIJgp8XMb/iL8giH4v8++/P8jNx5KMXDS/li1mIiEYcLu63KuIjjiIv3Zbejxiceq4iIpFqXEzzdgyfSQxbM/IrGuqbkDoo4XDiCdbLOrQs7QrWNWF78ZDDcGGbJnBX4PltLR1Mthxu93sL/3PQWGOx2pJI8EZFso1G0Ij3kw0c/IVQfSv7AcIKrXoQdifTIu4gjXN95ohluivDC39+iWiN7RUSyjhI8kR6yeXV11w5M5jaKRJPBiMNFOq/XfMacZ75I/PwiItIr6BKtSA/xJ7lebdq55slZOhZuDPP+v+ex72m78NmsxXz66lfUbqonmB9g1E5D2P2o8ZRp2TURkV5HgyxyjAZZ9F4vXvcK7903h3BTgr1szYI+LJjg32JBPxZIMJHMD3R6b1+s2n75GNDUYjSvP+jDRRxjJ43gmJ8eQEFJXmLnFRGRRHV5kIUu0Yr0kEln7JbU1daYRBM2AF8S/xYkMd9dU32oVXIH3j16kbBj0ZwV3P2zZ6ivbkj83CIiklZK8Pq4SCRCJJJkj5J0yYBt+jNyjxH4Akn82vkt8YmHjcSTNp9hySSDcYSbImxeW8vj185MSX25RL9f2UXtlV3UXvEpwevjqqqqqKqqynQYfcYJ/3M0BWX5mD9+cuUP+rCADxJc4xaAgD/xZDCYXK9gZ/WGQxEq561h7ZKNidfbB+j3K7uovbKL2is+JXgiPah0SAnfemga/UeVEygIbH13hQ/8+X6G7DCYaTO+SV5BsNM7MAL5AQZuO4BAUecTIvv8hgV9yV3KTaLH8b0nP0u8XhERSRsleCI9rN+IMr739Hl888bjGbPvaAIFAcwgWBhg/CHbceYd3+S8B89k9OSRnDPjm5QOLvaSwTb8eX58fh+7nbgj37n/NPY6dRf8QR++9noHDQJ5fsYfPJYB2w7An2jS5k/8Um4kFGHx3JWx1y7iaKwPJTQ3n4iIpJZG0eaYZEfR9pHFmrNaJBxh4RtLePueD1j1xVoioTCFZQXsfMwO7HHyLpQNLYmVXfXlOmY/9BGf/me+N8+dGS7iGHfQNux16q6M3HUotRvruf+Xz7FpVXWHa9z6/OatdBHwJX7ZFygqL+DwH+zDu09+xor567w1cQ1G7jCYvY6fwLaThuPz9Z2/K/X7lV3UXtmlj7SX1qIVjxK83JZoe4UaQtRsqMNFHEX9C8krbH35tqk+xEcvfMnbD31Mzfo6/EE/OEeoMcKIHQdRPrKMebMWE+lkbdu2fH7DX5LXbuIYyPNTOqCQUy8/hLJB3tx5K7/awPvPf8GKhesJNUUoLMljpwO3YacDx1BQnP3Truj3K7uovbJLH2kvJXjiUYKX21LdXs45Vi9YT/X6Onx+HwMqSuk3rJTVX23g7p8/SySU5Ai1gMWdh8/nNwpL8znp11/jhenvsW5ZFZE2q2oE8vxEIo59jp3A/t+Y2KoHccPqajauriEScZT2L2DwyH5J9TD2NP1+ZRe1V3bpI+3V5X/gtJKFSB9mZgzdfiBD27w/ZGx/Bo3ux+qFG5KrsJOJkyNhR+3mBu697EWITpTcVqjR6/179+nPqdlYx9fPncQXHyzn3ee/ZE3lJgIBP5g3crekvIC9jhjHzgdsQ6CdkcGRiMOMXp0EioikgxI8EWnX187eg8eueSXxlTcSmE4FwIUdLuzAFz/xCjWF+fT1xayq3MS6VTWEmrzErzEcipXZtLaWVx76iLmvLeLUnx5AYUk+G9fWMGfmIj5+ewkNdV7ZkvICJh08lp33HU1hm0u/VRvqWDBvNfU1jQSCfoZUlDF6+4FKCkUkq+kSbY7RJdrc1tPt9d6Tn/HK9PcJx7lUaz7DgbekWjJJUfOULZ0xIM8fdxJnn98YNKKMkTsOZu6sxZhBONz637ZA0IeLwNRpu7HjniNZVbmJ1579nKUL1nuDSsLR3j6fUVAUZO9DtmP3fUe3GkVctbGOpV9toLExRF5egJFj+tOvf2G7Ma1dU8NXC1YSDjsGDChj9Nj+FBR0PpWNZI7+PcwufaS9dIlWRFJv8nETKBtczMt3vkf1ujrCobA3MhYI5nuXRAv7F1K1oS75Hq9E/7h0xEbjdiQSdqxZVsW61dVEOogjFO2JfP6+uayqrGLOO0tjiWuk5aXisKN6UwMzn/6MpfPXcexZe7B8yUbefGUhlYs24A94l5bNZ4RDESq2KWf/Q7Zj5Jj+AHz52RreemMRa9fU4DMDHGY+IpEIO0wcyr4HbEP/AUWJfXYRkS5SD16OUQ9ebstUeznnWDZvDZ+/sYTajfXkFQYZMWEQEw7chkeufZWln65JvlIfWF6CK2oEfIlNuBzt7Uso2cyP3ysIXq/fsNHlrFi5easewZb8fmPKUTuwbkMdH89dQaiDHk/zgd/v4+TTdmPk6PLOY5QepX8Ps0sfaS/14IlI+pgZIycOYeTEIVvtKx1UDLYm1rOXeKWpia2VBHr7YkKRTpdsCzVFqFywvtNkMBx2vPz0Z1h+gHh3LLoIhCIRHvnXXM761p4MjE4XIyKSan1nxlERSYtdDhmb+MoYLXUy4rbLEr0qEXaJl01guhjnIFIfSqjOUFOEWa8sTOzcIiJdoARPRLpl1MQhFPcrSO6g6GCGhLW3/FqHdaehazDO5dmttDP1S3sWzl9HTXWjV304QtXmBjZtqqepqf3VRUREkqFLtCLSLWbGEd/bi8eufy3uaNtWEhk928xHcklbpmc3CTtI4NZCn9+Y88EyahtCfPrJaiLRnj/nHGPHDmDPyRWM7OUTOYtI76UePBHptrG7Defoi/bFH/R12DPn8xv+oI9RuwyJjcBNSJyVMbbSGyY1TvCybygU4e23lvDxx6sIhSNEIs5b1cPBV1+t55FHP+GJJz7tcMCGiEg8SvBEJCUm7D+as689ggn7j8bn9xEsCJBXECBYEMAf8DHxoDGce/1UTvnNFEaMH9TuyhMt+QI+LN8PyVzKTeZeQJ8l3jOYTAyJ1hm9BTDSziVd57zLtl8t2sCTT36KZjsQkWRpmpQck+w0KTU1NQAUF2s0XzbIlvaqq25g5fz1NNY1eVOqjB9IftGWFSQi4QjvPPUZs5/7klBTmHAojIvgzTHnHMO3G8jBp+3KiiUbmfnIx51e+vX5DWeG8ye2mgbgTZ6caOKWTFm/dTo6Nybo63SwSSDgY+oR45gwYesRzJJa2fL7JZ4+0l5dviShBC/HJJvgiWRSJBxhwZwVrF68kabGEEWl+YybXEH/YaWxMu+/vIBXHvoIM9sq0fNFB19M3G80lYs2sHl9Xbs9Yq2O8RnOwAV8ifW2JbCSRisJzK+XbNkhQ0o4a9rurFhVzUefrWJjVQM+g4EDitht4lAG9m89cXIk4li8bBOfL1xHTV0TAb+P4UNK2Hn8YIoKt15No7EpzMq1NTQ2hQkG/AwdVERBnm7RFukFlOCJRwme5KKq9bXMffUr5rzyFQ11TQD4gz4m7jOKSYdtz+CKMqo31fOv//c6mzfWd9jjFwj6GLn9QALFeSz6fE1sdYv2NC9bRp4/0YGxXi9fopM3J9HT5/MZxQMLqa1vIhJ2sSkHfdFexSGDijn60O3pX17I5wvWMfPdJTQ2RgiFt3w+v99wEdhum3IOO2AM+XkBNlTV8/6nq/hi0YZYXQDhiGP70eVMmjiUQeWtl2KrbwixrqqeUChCfl6AIf0LWx3bknOOxqYIEefID/o7LCciHVKCJx4leJLLnHM0NYZxEUdeQWCry7ENdU188Noi3n9lIU0NYczn/esYDjv6DSxir8O2Y+LeowB45cl5zH1zSWzJsRjzVpvo17+QY6btzr8fmEv15gYicaZK8fmNwqI86l2EcCiBf1N9eMldMgNC8ju+nGtAMOhnt52H8eFnqwnHyUh9PqOsJI/9Jo/kpbeW4JxrN4E1A58ZX99vG8Zt05/VG2qZ8/lavlq+Cb/P+2IjzhHw+dhl+4HstO1Aigq8Xr/quiY+XbSeeYs30Bj9bs1g7LAydtl2AIPLCzM/GEYkOyjBE48SPBHvEuWyBeuormrA5zPKBxUzdFS/rcrVbG7gw7eW8OkHy2moa8If8DFsVDmTDxpDxZj+mBm1NY08ds8HrFlVTSQcaTVI1gx8fh+DhpRw0tl7UF3dyCMPzKWxMbRV72Bz71UEvHvvkk1w4iR4LfkTWKrN5zOcA18CvY0+n7HTuEF8tmQjLuLaXbDE7zMCfh8nHDyW1RvreOPjld4l9TaZo+H1io4ZWsKUPUbiV4+eSGeU4IlHCZ5I6jnnWFG5iffeWMziBesJNUUIBH2M3nYAe+6/DcNHbZmvLhyOsPDLdbz79hLWrK4mEnYE8/zssONg9thzJPfeP7drU58UJDbQwxfw4UtwlRBf0J/whNO+BBJHv9/wB3xxl2sDLyEcObiYw/cchZkRCkeoXFtDTX0IDMoKg1QMLNYlXREleNJMCZ5I7/af/3zJx5+s6nQwSCuGl+Al2OsXyE9wgITP8Cd6H2DQh/k6TxzNR0J1+n3GfjsNo6q+ic+XbQK2TCFoeMnijiPL2XXsAALpWtZOpPfrcoKnYVIiIj1o0qQRfPzJquQOSvKSrnMusXvckkgyI2GX0PLBLpLY+cMRxxufrMTfwSCWUMTx0eINLF5TzdF7jiI/6Pd6UjfU8dXKzdQ2hgj4fAwsy2fc8DIKO0hqG5rCNIUjBPw+8gM+3fsnfYZ68HKMevBEer///ncBH360MrFLtT6Sm3qFxO7Di5VNtLfPwJ/g1Cm+gKXsMrHPYEBpPruOGcB7C9bRGIq0urfPb959jSMHFrHvDkPID/oJhSMsXVvDZ8s2sbmuCTOvd7Ao388OI/qxzZBS8pKZFFskc3SJVjzJJngNDQ0A5OfnpzMsSRG1V3bpqL2cc7z08gI+/ngl4Q5G5/r9RhggL/kBGZlO8Cx6L14ifAEf1kky6IuObI73mXwGhfkBDpo4lDc+X01Tm0Swmd9n+Mw4eKdhDCj12mVjTSMLVlaxtqqeSMQRDPgYObCYMUNKKehgIIpz3oATn3oEM6aP/HuoBE88ySZ4GzduBKC8vDx9QUnKqL2yS2fttWTJRma/V8nixRvx+5sHaTgGDixi0qQKXnrzK0KJTLvSkkEg0UmKk7gHL6n79fyGL4UJHng5brCTdYnNvP8bBvydL0Pn9xkH7TiUeZUb2VTbiHO0GiHsM+/1NoNL2HXMAHxmNIUiLF1Xw6JV1dQ3hWP1jBhQxJghJZS2M4m0pE8f+fdQ9+CJiGSb0aPLGT26nNraRjZtaiDiHCXFefTrVwDA2qo6PvhoZVIDMhK9NJrOskn9LynBHjDnOr+3rzlJi7jOBxyHI46Zn64k6Pe1O/VL81e+ZE0NdY1hxgwuYe6i9YARadExEo44KtfWULmuhooBRey8TX98Zjjn2FTbxLL1NdQ1hACjOD9AxcAiylos29c6fkc44vCZaQSxdJsSPBGRDCsqyqOonf/pT9p5OHMTHHHbfJ9ZotOexFbqSKhwEmVJrmzX+yc6Fo64hBKkiIOwc3Evs0acY82mOtZW1XsTPLeTDrrof5avryMUcWw7pIRPKzfS0BRuNYBkU20jKzbWUpQXYKfR/SktDOKcY3N9E8vW1bJ+c0Os9ryAjxH9ixjav5BgO8l18yViAw0ckXYpwRMR6aXKSvP5xlETePSZeR3eqwdespYX9LPHLsP54NNVcVey8PuMkuI8wgb1jWHi3aXTnDT6Er2MizdNSqIJh/ktLclJMnceRSIutqZxh2WcV6nP4vcgRpxj9cY61lXVx1m+DWoaQsyev4bdxg5kxYZaNtY0bhVzYyjCkrXVLF5bzYSKfgwsLcA5R1VdEys31FJV1xRLBovzAwzvX0R5cV6rZNU5R3V9iNqGEBHn8PuMssK8du8rdM5R2ximocn7mQj4jZKCoCajzmJK8EREerHRFf04/YSd+c/MhaxdXwstlhbzR1ekGFVRxhEHb0dZaT5DBhXz6ttLaGgMEWqRFAb8RsTB9mP6c+j+Y2hoDPP4S19SUxdqtWbtlvI+CvL9HDipgpffWxY3aQRvsEEw4MP5E5x9xUjo3rt0SyYZdK7zK8oRBy7syOskMYo4mLNwHXlBHx11YzZ/j58t28S44Y7Vm+qoa9o6Ka9pCLFwVRV+n48JFf0oCPpZX93Ayo31hCKRVvMLrthQR0GenxH9i2I9iOurG1hf3ei1sUGsa9DVUVYUZFBpQatRx5GIo7phy89NwO+jJD+gy8q9jAZZ5BgNsshtaq/skur2WrOuho8/X8Omqnp8ZgzsX8guOw6lrHTrUbpLllfx+cJ11NQ2EQz4GDakhJ3GDaKwYMtAgEjEsWj5Jj6Yt5oVa2pi7w8ZUMikiUMZO7Icv89Yvqaap19fDLhWSWOzgN8oL83nmAPGMPvzNcxftqnDhNDvM2+EsHU+EKIln0Ggk0EWLXmJU+cMCCY6IMQs4UmXgwFLaIRt8zJvicgLJDZ4pH9xHpvrQ/F7Z4ERAwqpbQxT1xBu9z7EZj6D0YNKCAaM9dUNVEdXHGmZOAKUFAQZWJIXvZS9RSgciSWrAX9i30si+si/hxpFKx4leLlN7ZVdsqm9nHM0hSIEO5gMuKEpzBeLNzDni7VU1zXF3h8xuJjdxw1i1LDS2OCCyjXVzF2wjlXra6OX+IywcxQXBNh1u4EMKivgmfeWtRqs0JmO4mqPN+I2wQTPaPcet/b4DAL+xJLMgN8SvryZl+Bn8ydYp5e0JhZn0O9LqOfNZ17S3FnvrN9nVPQvIuAzahvDVNU30RR2sSzFAcV5fkoLgt2eizCbfr+6QaNopbX6+vqEysUbldbQ0EAyfwAEAgECgdY/UpFIhMbGxoTrAMjLy8PX9i/AUIhQKJRwHWbW7txI2f6Z2pPtnykX2ymbP1NDeOv3mj/TLtsPYuftBhIKO0LhCC4Sio1WbYzOSQYwuCzI1/cYRk19iM21TUScozA/wKB+hQSDXg/i5O0H8v6CdZ1e+gUvsUnmXr1k7htLpmymBzOEwy6heL2RxPEHjzQLhSPk+TpPBiMOGpsiBDpJysIRx7L1NQT9PiJYrGewZSvXNIapaQxTHIDiPF+0fkdDGJoi0Z5Bg6BBQcAoKizY6jwtf58S+f9dtv4bUVCw9WdPlBK8HBSJRJg/f35CZYcNG9bhviVLlsQmkkzE4MGDGTp0aKv3GhsbE46l2fbbb7/VD/W6detYs2ZNwnXk5+czbty4rd7P9s80ePDgrd7P9s+Ui+2Uy5/JzAgGvPvtvvzyq4Q/Uy3gC2/5TLuMGUDQ7+Otz1d79xW26ajw4d0PNqIkwpr6YNxLiFt4f7Amcy9YcjO6pKfe5lveUikSieBLoLcxmWQwksBUNeAleY4I/k7OX90UYd26tQSCQQpK+oFzW9Y7dlAfjlAfNuoiDfQvaj2AZMmSJfTv3x8god+J3vr71FJ7/0bsvPPOScXXkhI8ERHJiAmjyhnRP493PlnE6ho/jREv0cnzO4YVRxhQGMFvkOeHZTVBXJw0yAwC5qJTtCSWMgWS7BlM5t6xZDr70tEvmOqRxC3r7fSzmTegx9dJMmjmI79sIN6thVvfX9ic7DWGIqzd3MCg0vxYG/gCAZzf6w0uKC6lsb6WSLid7uc+TAmeiIhkTF7Ax7DiCMOKO16Xd1BhmIDPsawmSChirRI9v8+beHhISYDBwc00RnwsqSnA4dpNCJsvGpYXQF3Yn/B9gAFf4veL+ZOYKiZdkrm7PplbB1J9174ZuAS+q7BzbKxtpDg/QG1jmNJBw2P7SgcOwcxorKuheuN6Qo2J95TnMg2yyDHNgywS7Xauq6vDzNq9STVb7hlqT7bfB9UeM6Ourg5ofVNxtn+mXGyn5s/U8ibwXPlMLfX0Z3LOsbqqgWXrG6hrDBPwGwNLC9h+eCkBH7HPVN8UZtGaWr5aXUM4Eu1FinY9bTOokLFDiinOD/Bx5WaWr6+Lm+SZwYDiIDUNkcSmf8Eb7ZpoghfwGf5EB3n4SHjErTd4JMGl5RJYBq5ZogM9wIs10e8hmGAPot8X77v1GqjAB20HUfv9fgKBQKtjs+H3qaCgQKNoxaNRtLlN7ZVd1F6ZFYk4NtU20hSOEPD7KCsMtkqQnHMsWr2ZL5ZX0dAUia0O4TPD4ehfnMfEUf0ZUJrP3K/Ws2pjXYdJni86ctdnCc4DGJXoCFrwpl5JPHFMbHQseIlYoklbUjGkIcEzSCgh7lcYJOgzGsMR6kOtk/M8v1EQ8GfLJM4aRSsiItKSz2f0L9m6p6SZmTF2aBljhpSytqqe5Ws3EY5AaUkRw8oLKS3cMmfg7mMHsGJDHQtXbmZzfVPsXjAH+M0YM6SEMUNKqGsM8f7CzkcH+31GXtCf8DrDviQv+yZzhTipleWy5NLz5vqmrXramjWGHY3hEEG/URz0xz6Tc46miIv16PrMCPhSN29fT1OCJyIifZqZMbhfIUHn3btVXl7WbpkRA4oYMaCIzXVN1NR7y3/lB330L9ly838wkMde2w9iXuUmqmobW08IHB3/UV6Sx44V5ZjBh4s3tJoIuC2fwcDSfDbVJX4pMeBLvOcs/iXPNvUm2MsG0bWOM5gYRRxYJ4M8msKOahemKOijKewIbXVF09EY9hL4/IAv6xI9JXgiIiJJKC0Mturda6u4IMie2w+itiHE8vW11DV692cV5wcYPqCIwrwt/+udtK23Hu3y9bVb7hfEu7xcVpTHyIFF9C/JZ1NtI1+s2NTp6hQjBxaxvrqRpjhrF4OXgBUG/QnNQ2g0J2ydFo1pu5pFPOm6UhqJOPydJKWhiKOmMRz3cnbYOWqbwhQF/VmV5CnBExERSYOi/ADbD9+6N7ClgN/HqEEljBxYzOa6JprCEXxmFOUHyA9uGfjQryiPnUf1Z+XGWtZubojemOV1CToH/YrzGNG/iJKCIIPLCvlqdTW1DaGtLmk2jzfpX5zHqIHF1DeFWbqulohzWyWPzRPODCjJozDPz+rNiY1O9SfRgwjJJXjpSLDCCfT2AdQ2hSmOtkmE5vn+onHh9Yb2pgRQgyxyTLKDLJpHEOXl5aUzLEkRtVd2UXtll2xpr3Akwua6JsIRh89nFOe3v+xXXWOINVX1VEcvJ/vNKC/OY1BpQaul3JxzbK4PsX5zAw0hby45n88oL8qjf3FebGBKTUOI1VV1OLf1vXDNvXxDSguoagh12oPYzGfJrzyS8JJ1JDYgw6s3seQx4IvflWlAXhIxJkCjaMWTbIInIiKSKOccNQ0hNkV7G8Fbz7ZfYZDifG8akkjEsXpzA03hSNxBESV5fhrCHc9/2FYyCRt4CVuiI4n9RuJlE0jg8lLXm6dRtCIiIpJeZkZJQZCSgo7vQfT5jKFl+dQ1hamqD9EY2pLEGVCc76e0IEjQ76O6oYnN9Z3PIWcknoBtiTWp4inVGHHkZ3jCayV4IiIiklJmRlFegKK8AOEWU48E2vR+leQHCfh8bK5vItTOgA+feYNTGkMRErzq652fJJOrNORhESCx6aPTQwmeiIiIpI3fZ/jjZFAFQT8FQT+NoQgNoTAR5/CZEfT7yI9OBB0OOjbUNiY8D17SvX1JlU5MKIFRvOmkBE9EREQyLi/ga3ewCHhJYv+iPDbVNRGOM3Yg4PNWG05qUugkyycq0yMclOCJiIhIr+cleUGawo7aphChsDdNieElh0VBPwG/j9rGEA1JXM9NYtxG0smgS2D6lXRRgiciIiJZwczICxh5gY6nsikM+nGEaexssme85C6pOfuSHuihS7QiIiIi3WZmFAX9BP2O+qZwu4Mz8v1GQdBPxDnqQ4lN1eJLcn67JDoG0yLT5885ZnaFmbkkHv+TyXirq6uprq7OZAiSBLVXdlF7ZRe1V3aJ115mRp7fR1lBkLL8AMV5foqCfkry/JQXeKN7fWYEfD4KA/GGgHg9fYEuzGuXzATO6aAevD4uFOp8/iHpPdRe2UXtlV3UXtkl0fbqbBSv3+f1+IWdoym8ZUoXbySv4Tdv4EZjAuv2NjMy34OmBC+9/gEs6KTMuz0RiIiIiLTPzAiY0cEgXgwI+qApwSQvxcuVdYkSvPT6l3PulUwHISIiIt3jN8OiSV5HaZ4PCPaC5A6U4ImIiIgkxGdGvt+IOEe4RaJnBgHrHYldMyV4IiIiIknwmeHL4CoVicj0PYAiIiIikmJK8NLrMjP7wsxqzazGzBaZ2aNm9h0zy890cCIiIpKblOCl16HAOKAQKAK2Ab4B3A4sMLNDMhibiIiI5Cjdg5c+K4DXgC+AGmAAsC9wIN6I6wrgP2Z2knPuiWQqNrPKOLuHOefYuHFjQnW56Hw/bcsHAgFKSkriHtvY2EhtbW1C54mnsLCQ/Pz4HZo1NTU0NTV1+1xlZWX4fPH/rkn0u4snXd9fe+2l76970vn9tW0vfX/dk+7vr7m9qqur9f11Q0/9/DW3V2dy8ftLqJ5u1yBtvQVMAWa6dn76zGw34D5gIuAH7jWzCc65ZT0apYiIiOQsSzQDltQxs8HAe8Co6Fs3O+cuSlHdlRUVFRWVlfE6+bZo/mujvLw8FaeXNFN7ZRe1V3ZRe2WXPtJeXR6q26fuwTOz55JcJzbe45WuxuGcWwNc3uKtk7r94URERESi+lSC18s8DkSi28PNbJsMxiIiIiI5pK/dg/c0sChFdX3RnYOdcxvNbD0wKPrWEGBxt6MSERGRPq9PJXjOuRszHUMbugFSREREUk6XaDPEzPoBA1u8tTpTsYiIiEhu6VM9eL3MCWxJsFcBS1JU75AVK1YwcuTIhAo3j6LuTQskS8fUXtlF7ZVd1F7ZpS+017Jly2Y75/bsyrFK8DLAzAYCV7V46/H25szroqZIJMKyZcsS7REcFn1emaLzS3qpvbKL2iu7qL2yi9orDs2Dl0JmdiGwI3CTc+7LDsrsAtwL7BJ9qw6Y4JxLVQ9eUppXxXDOJdblJxml9souaq/sovbKLmqv+NSDl1oFwI+AH5nZJ8C7wFKgFugP7AN8jS0TF0aAszOV3ImIiEhuUoKXPjtFHx1ZCpzrnPtvD8UjIiIifYQSvNS6E29+vP2BvYGRwGCgHKjBGyn7LvAU8LBzrvurH4uIiIi0oQQvhZxzVcAz0YeIiIhIRmgePBEREZEcowRPREREJMdomhQRERGRHKMePBEREZEcowRPREREJMcowRMRERHJMUrwRERERHKMEjwRERGRHKMET0RERCTHKMETERERyTFK8ERERERyjBK8PsjM9jez283sCzOrMbMNZjbXzP5oZmMyHV9fYGYVZnaimV1tZs+b2Tozc9HHoi7UZ2Z2qpk9YWaLzazezFaa2Wtm9mMzK0nDx+gTzCzPzI40sz+b2SvR77XRzDab2Xwze8DMvmlmCa/trfZKHzObaGbfMrNbzOz1aBttjLbZmuh715nZDknUqfbKADP7Tot/F52ZvZLgcWovtJJFn2JmfuAG4IeAdVCsGrjIOXdXT8XV15jZdcAlcYosds6NSaK+ocADwJQ4xb4CTnXOzU60XgEz+x5wHVCeQPGPgGnOuY86qVPtlUZmthHol0DRMPB34JfOuXCc+tReGWBmw4FPaf2796pzbkonx6m9ohL+i1Nywk3AhdFtB7wAvA0UAEcBuwAlwHQzq3HOPZyRKHNfQTvvrQcGJFuRmRUBTwF7Rt+qAx4BvgCGAt+MPo8FnjOz/ZxzX3Yl6D5qD7b8DyYCzAHeBFYCecBewBF4V0N2AV4xs4Occ5+2V5naq0ctAD7Ea6uVeAndCOAQYEfAD/wUGASc014Faq+Muhnvd68R73etU2qvNpxzevSBB3AkXlLngHrgmHbKXNWizEZgYKbjzsUHcDnwPPBH4ERgJDCmxXe/KIm6rmtx3EJg+zb7i4HnWpR5LdOfP5sewD+AFcBvgBEdlNkDWNLiO35D7ZXRNjsFGNpJmTOBhhbf86Fqr97zwEvEmr/Tlv9feqWT49ReLT9vpgPQo4caGt5r8UP96zjlWv7wX5/puPvKoysJHjAEqG1x3H4dlCuPJinN5Y7K9OfNlgdeD11hAuX2xOvha/6Od1N79e4HcH2L7/gWtVfveAD98XpcHTAD71Jrpwme2mvrhwZZ9AFmNgGYFH25Ge9SbUeuabE9zcw6uldPMu9koDC6/bJz7s32CjnnNgK3tHjrrDTHlTOcc+865+oSKDcb73aHZvu2U0zt1bv8t8X2mHb2q70y4294l1HXAD9P4ji1VxtK8PqG41tsP++cq+6ooHNuJrA2+rICmJzOwKRbjmux/WgnZVvuPzY64EZSa1GL7fbup1R79S4t26i2nf1qrx5mZocD50Zf/sw5ty6Jw9VebSjB6xv2aLH9RgLlW/7ls0eHpSTTkmnXT4BN0e0yYNu0RNS3DW+x3d7/mNRevYSZ5QMXt3hrZjvF1F49yMyKgVujL//jnLsnySrUXm0owesbdmyxPT+B8i3L7NhhKckYMysHhrV4a0G88s67+WRhi7fUrilkZkOA/Vu8NavN/nLUXhllZkEzG2lmp+P9Edt8GX0hcEebsuWovXraNXiXyuuA7ydzoNqrfUrw+oaWP/jLEyi/rINjpfdo2S7VzrmqBI5Ru6bP1UAwuv2m23qaFLVXBpjZ482T5OJNt7EUuJ8tvT3vAAe1c9uK2qsHmdl+wEXRl1c55+ImaO1Qe7VDCV7f0HLW7poEyre8H6U0xbFIaiTbpqB2TQszOxY4P/oyAvyinWJqr96lCjgD2Nc5194fvWqvHmJmecDtePnIR8D/dKEatVc7NNFx31DYYrsxgfINHRwrvUeybQpq15Qzsx2Be9myMsy1zrn27v9Re2XGdLZcLi8CRgGHA9sAdwPHmdn32+nxUXv1nN8BE/H+OPqucy7UhTrUXu1Qgtc31OFN8AiJzQie3+ZY6X1atktCs7yjdk0pMxsJPIt3kzbAk8AVHRRXe2WAc+7fbd8zMx/wPbxlys4EdjezKc65NS2Kqb16gJntAvw6+vJ/nXNvxysfh9qrHbpE2ze0vL+kuMNSWxS12N6c4lgkNZJtU1C7pkx0UMWLeD1BAC8Bp8TpfVB79RLOuYhz7ha2XEqfCPylTTG1V5pFpya5A+/e1WXAb7tRndqrHUrw+oZVLbZHJFC+ooNjpfdo2S4lZpbIPSRq1xQws4F4yd0O0bdmAcc75xo6Pkrt1Qv9L969eACnR6fpaKb2Sr9v4a0UA3BRggMjOqL2aocSvL6h5Yi+7RIo37JMu4umS2Y55zbgLefTLG67RlckaTnXk9q1C6LTMbwA7BJ96x3gaOdcexPlxqi9eh/nXBPwefRlEJjQYp/aK/1Gtth+rHm0c9sHrVccObjN/nJQe3VECV7fMKfF9v4dFWphvw6Old5lTovtztp1ItAvur2Z1nNASQKivQLPsWXZvw+Aqc65RC/vzGmxrfbqHVreh952NYM5LbbVXr3fnBbbai80yKKveBK4Lro91cyKnXPtDiU3swOBwdGXy4HZPRCfdM2TwJHR7ZPwLjl15Bsttp9yzoXTFlUOil6+ewbYJ/rWR8Dh0XUtE6X26kXMrIjWE9y2nS5F7ZVec4C7Eig3DJga3V6F90dWs5YjZtVebTnn9OgDD+B9wEUfv4pT7ukW5f6c6bj7ygNvBvfm731RgscMxRv95fCmGNing3L98P7n1Vz/MZn+vNn0wJtC4eUW39+nwJAu1KP26kUP4MoW3/E8tVfvfABTWny3r8Qpp/Zq89Al2r7jdy22rzSzI9sWMLMrgKOjL6uA63sgLuki59wq4MboSwPuM7NW955EeynuY8s6qW84557uuSizW3QS1keBQ6JvfQEc6pxbnWxdaq/0M7PrzOwX0VHOHZUpNLMrgctavN12FK3aK8uovbZm0YxW+gAz+z/gguhLBzwPvA0U4HVt79Zi3xnOuX/1eJB9gJkNBm5p83YRcFR0uxZvfrVWnHPfbKeuEuBVttwXVgc8jJeIDAFOYcsyPBuA/Z1zn3XzI/QZZnYbW1apALiB1kscdeRj59xzbd9Ue6WXmc0AzgXCeLeXzMG7+b4Ob77CCcBhbLn/Cryly6a5dv5nqPbKPDObwpaBFq8656bEKav2akEJXh8SnXfoRuIv5FwDXOycm94zUfU9ZjYG+CrZ45xz1t77ZjYMeAA4OM7hi4DTnHPvJHvevszMXiH+99qRu5xz53VQp9orTczsTrzpNxJRD/wRb/WRDu/BUntlVjIJXrS82itKgyz6kOg/Yj8ws3uB7wBfw5sXrxFYgnf/3a3OuaSTD8kc59xKMzsEOA2YBuyO99fqJry/XB8BbneJj/aUNFJ7pdWFeD1yhwB7Atvj9djk4f3xuhb4EO+eyvudc2s7q1DtlV3UXluoB09EREQkx2iQhYiIiEiOUYInIiIikmOU4ImIiIjkGCV4IiIiIjlGCZ6IiIhIjlGCJyIiIpJjlOCJiIiI5BgleCIiIiI5RgmeiIiISI5RgiciIiKSY5TgiYiIiOQYJXgiIiIiOUYJnoiIiEiOUYInIiIikmOU4ImIiIjkGCV4IiJdYGaLzMyZ2aJMxyKJMbMp0TZzZnZFpuMRSadApgMQkdxhZq6DXWGgCtgErAXmAu8BzzvnFvZQeCIifYYSPBHpCX6gf/QxBtgT+A7gzOwl4Drn3EuZC09EJLcowRORdLkG2NDidTFQDmwH7AMMAQz4OnCYmd0C/NQ519jDcYqI5BwleCKSLrc55xZ1tNPMvg78FpiCl+j9ACgDzu6J4LrLOTcm0zGIiHREgyxEJCOccy8ChwJ/aPH2WWb2wwyFJCKSM5TgiUjGOM/lwL0t3v6tmRXEO848J5nZfWa20MxqzKzazL40s9vMbK9Ezm9mU83sHjP7Inp8g5ktM7M5ZvaQmX3XzEZ1cGzCo2jN7Egze8zMlptZvZktjtZ/aHT/mBajO2d0UMeMFmXGRN+bHH1/QbTedWb2XzM7z8xS8u97i3O+En2db2bfM7OXo99VyMw2tjnGFx2xep2ZvWpmK82sMfodf2VmD5jZN8zMunB+f/Tz/Tdab72ZLTGzu81s1xR95rxo+zSf+10zG5iKukV6jHNODz300CMlD8C1eIxJ4rgRQGOLY8+NU3Ys8G6bc7X3uBnwd1BHPvBIAnU4YFYHdSyK7l8UJ1YfcEcn9f8Fb+BJ8+sZHdQ1o+V3C/wKCMWp999AXgrb9JXoeee0c66NbY75OMHv9lVgUBLnHwj8N059IWBanLqmtCh7RQdlioHnW5T7L1Ca6d8tPfRI9qF78EQk45xzy83sSeCk6FuHAne1LWdm44DXgKHRt1YDzwILoq8nAicAhXj39OUD57dzyv9pca4Q8B+8qVs2AkXAaGDvaH3d8Xfg2y1ev4mXqNQDuwDHAT/rQr3fBS7FS4qfBD7Em4pmH+AovHsaj4+WuaJLkW+tAHgM2A1YBjwFLMUbGb1Pm7KDos9fAe/gtc9mvO92YjTGIuBrwBNmdpBzLtzJ+QN4SfnB0fM/2eL8JwDj8EZr32Fms51znyf7Ac2sP/A0sF/0rf/f3v3HWl3XcRx/vkO7cYdIlmFqNYNNFkRYd6YJ5crKsDKqVfya1CaJY8ulm/MHlYtysjL6sYGLNjEovS1J2Kh0gDqyWoDLZGn+1gZFwg0oJEHe/fH+ns6Hb+fH93s5l3Pu2euxnfH98Tmf7+fcMy5vPj/en3XAZ9z9YNm6RNqu3RGmXnrp1T0vBtmDl733quS9T9e4fyKwLSlzM9BTo9zpwO+Tcpfk7o8iAiwHdgMTG7TprcC8OveepUEPHjANOEK1Z2lOjTITgOeSckV78JzoJRtXo9ylRLDnRN7B3hZ+pw4sA17T5D03A5Mb3H89EZhX6vxCiecvzX/vRPB3Z1JmRZ26LkzKfC137zQiUK7cXwWc0O6/U3rpNdiX5uCJSKf4Y3Jca97bTOCc7Pj77n6du/8nX8jddxA9Y3uzS9fmikwgevYAVrv79noNcven3f32Am2v5VqiJw1gsbuvqlH/Y0RPopesex/wEXd/Kn/D3e8B7spORxO9oa2y0d0XeJMerey7eaTB/ReBTwE7sku1ellrWePuV+W/d3c/DFwBHMgufaJgfQCY2VnAZqJXFWJ4f25Wr8iwpABPRDrFQHJ8gpmNzt3/YvbnyzQZdnT3XUSPDsAFZjYmuZ3+3htZvpnNZc+7ODt9iZhnV5O7byWGBctY4e4vNLj/8+R4Ssm6G1ncqorc/QAxBArQZ2Y9jcpnvtqgvn8CG7PT15nZmUXaYWYTieBuXHbpG+6+0N3LBt0iHUVz8ESkU/wrd34S0VOFmfUClZWxW919T4H6Kj1IryKCnPuz88eAQ8SQ71wz+yXRM9TKf9D7iPlgAA+4+/4m5dcRvY5F/brJ/b8kx28oUW8jB4AHy74pW/E7mZg3OYrqzwViCBziuziL+G7q2enuf2ryuPzn/muTtp1LDBWfkl26xt3rBuMiw4kCPBHpFCflzvclxxOIIADg/AZ73tZzauXA3feZ2W3AQqqraZ8zs/VET85D3iBBc0Hjk+NHC5QvUib1fJP7abA8qmTd9TzpzRdCAJHGBrgMuBqYVLD+MU3uN/vMUO5zX0C0bxQxB3K+u/+owDNEhgUFeCLSKV6bHB/K9Xodaw6y3tz51cTvv/lED99bgAXZCzN7CvgpMddv1yCeNyY5HqhXKFGkRzLVbFVnGgC3airO3uZFIgce8GNgVsn6G+Y+pPlnhnKf+6Lk+C4Fd9JtFOCJSKd4R3Kcn1+W/q7aQnURQVFb0hOP/W4XmNlSYA7wAWJYtdJLOA64EfiSmc1x97Uln5cm8O2WuVxHCpabTzW420+sur0XeAJ4ETjo7kcAzOwm4CstbmdR9xG9eL3AzCy1yq1taotIyynAE5FOMS05zs/12p0c73H3b7XigR650hYBi7J5fucB04mgbywxbNxvZpPc/ckSVae9dqfULVWuzHCxMPvzMHChu29rUPbk49Ceeh4CbiHmP44Evm1mh939e21sk0jLaBWtiLSdmZ0OfDS5tDFX5AmqPUjnmlnL/3Pq7gfcfaO7X0P04P02u9UDzCtZXRoMFpmDVnSeWkfLguSJ2emDTYI7aPPndvcNREqVStqV75rZle1rkUjrKMATkU6whOrw6A6gP73p7gPAw9npGCKZ75Bx938TCXUrzi5ZxRYi2TDA+8wsv4Akr8wK2k6WzqNsOK/QzMYCU4e2Oc25+71ELsKXs0s/MLPL29gkkZZQgCcibWPhJmB2cnlxrQTGwPLkeEm2rVShZxxLGzMvlSmc5WSrpDIZSSzqqMnM3glcMuiWdZYBqnMOpzT52d9CNeF0W7n7euDTRPocA24zs8+3t1Uix0YBnoi0hZm9H9jA0ZPsV7r7sjpvWUlsdA+RP+1+M6s7xGdmbzazRcRq2PT6bDNbbWZT6wUgZnYacENyaXPDD1PbkuT4RjObnS9gZmcDd3P0ooxhK0teXBmWHQ98PVtV+z9m1mNm3yHSqHTMAhR3Xwd8lpg7aMAKM5vb3laJDJ4WWYjIULnczNLFBr3EpPrxxOb0Y5N7TmwPVbeny90PmdkMItg6g0ie+4iZbSb2nt1N9JadQayInUz8Q70hV9WJxCrPWcBOM/sd8DiRBmQMMRw7HXh1Vv5xIu1HKe7+gJktJ7bQGgGsyuZ3bSLmfE0ihpp7gFuBLyc/i+FsCdVVzjcAHzezTcAuYgu6jxH7Be8C1lJ8m7Ih5+5rzGwW8Z+CEcDtZvaKu/+kzU0TKU0BnogMlesLlHEiAPumu29qWtj9WTPrA+4APkgEcNM4egVuvv4/564dSo7fCMxo8Mg/AJ9091JDtImFRNB5WXb+nuyVWkqkEqkEeIN9Vkdw934zmwJcl116O9U9XiteIIZEpx/HphXi7j8zsxHAKiLIuyNbXdvf5K0iHUUBnogcD68QuwzsBf5BDLVuBX7l7s+Uqcjd/wZ8yMymAp8D3kv02p1MbKe1A9hO9JStdffnc+9fbWbbgA8TedAmAmcSPYwHs/dvJRZ63FPJ2TYY2c4P88ysn9hL993EQoRdRPC4zN3vM7Nzkrft+/+ahhd3v97MNhAB7vlEGpgB4BlgDfBDd99jZh0X4AG4+53ZSu2VRJC3Ogvy7m5z00QKM+2nLCLSXmY2E6gMA17ZYB6iiEghWmQhItJ+FyfHD9ctJSJSkHrwRETayMzeRqw87QF2Am/KhnZFRAZNPXgiIkPEzC41s/lmNrLO/XcB66nmg1uu4E5EWkE9eCIiQ8TMriBWyO4jFn1sB/YTCy3OI1b/VnLgPQr01UnyLCJSilbRiogMvdFEzrt6W6z9Bpih4E5EWkU9eCIiQ8TMeoGLiHxvfURy51OJ3RL+TiRo7nf3X7SrjSLSnRTgiYiIiHQZLbIQERER6TIK8ERERES6jAI8ERERkS6jAE9ERESkyyjAExEREekyCvBEREREuowCPBEREZEuowBPREREpMsowBMRERHpMgrwRERERLqMAjwRERGRLqMAT0RERKTLKMATERER6TIK8ERERES6zH8Bf6IBO9dTn2YAAAAASUVORK5CYII=",
"text/plain": [
""
]
},
"metadata": {
"image/png": {
"height": 302,
"width": 316
},
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"plot_design_scores(scores, cutoff=cutoff)"
]
},
{
"attachments": {},
"cell_type": "markdown",
"id": "91259c25-6cde-43d7-abdc-7041a3e29db9",
"metadata": {},
"source": [
"## Verify design with counterfactual prediction\n",
"\n",
"The designed perturbation can certainly be passed back to the\n",
"[counterfactual](api/cascade.model.CASCADE.counterfactual.rst) method to verify\n",
"whether the target markers are up-regulated, just like what we did in\n",
"[stage 3](counterfactual.ipynb).\n",
"\n",
"One notable difference is that we need to specify the `design` argument of the\n",
"[counterfactual](api/cascade.model.CASCADE.counterfactual.rst) method, which\n",
"tells the model to use the designed interventional scales and biases\n",
"(instead of those from the training set) when making counterfactual\n",
"predictions."
]
},
{
"cell_type": "code",
"execution_count": 25,
"id": "ef02ad10-2cde-41df-9b45-040cde2ee79e",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T10:07:47.756140Z",
"iopub.status.busy": "2025-03-23T10:07:47.755906Z",
"iopub.status.idle": "2025-03-23T10:07:48.126901Z",
"shell.execute_reply": "2025-03-23T10:07:48.126122Z",
"shell.execute_reply.started": "2025-03-23T10:07:47.756121Z"
}
},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"\u001b[32m18:07:48.121\u001b[0m | \u001b[33m\u001b[1mWARNING \u001b[0m | \u001b[33m3810563\u001b[0m:\u001b[36mdata\u001b[0m:\u001b[36mencode_regime\u001b[0m - \u001b[33m\u001b[1mOverwriting existing regime \"ctfact\".\u001b[0m\n"
]
}
],
"source": [
"source.obs[\"my_pert\"] = \"KLF1\"\n",
"encode_regime(source, \"ctfact\", key=\"my_pert\")"
]
},
{
"cell_type": "code",
"execution_count": 26,
"id": "1a509af4-f1ce-48e6-b835-4877e837184b",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T10:07:48.187539Z",
"iopub.status.busy": "2025-03-23T10:07:48.187335Z",
"iopub.status.idle": "2025-03-23T10:07:56.707109Z",
"shell.execute_reply": "2025-03-23T10:07:56.706157Z",
"shell.execute_reply.started": "2025-03-23T10:07:48.187521Z"
}
},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"\u001b[32m18:07:48.188\u001b[0m | \u001b[33m\u001b[1mWARNING \u001b[0m | \u001b[33m3810563\u001b[0m:\u001b[36mdata\u001b[0m:\u001b[36mconfigure_dataset\u001b[0m - \u001b[33m\u001b[1mOverwriting existing `regime` = \"interv\".\u001b[0m\n"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"\u001b[32m18:07:48.370\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[33m3810563\u001b[0m:\u001b[36mcore\u001b[0m:\u001b[36mpredict_mode\u001b[0m - \u001b[1mNumber of topological generations: [68, 88, 70, 101]\u001b[0m\n"
]
},
{
"data": {
"application/vnd.jupyter.widget-view+json": {
"model_id": "a9a049b623ba4c848472dac8057dd989",
"version_major": 2,
"version_minor": 0
},
"text/plain": [
"Predicting: | …"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"configure_dataset(source, use_regime=\"ctfact\")\n",
"ctfact = cascade.counterfactual(source, design=design, sample=True)"
]
},
{
"cell_type": "code",
"execution_count": 27,
"id": "d5ced311-eea9-4979-9c4d-3633beeb94d9",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T10:07:56.708902Z",
"iopub.status.busy": "2025-03-23T10:07:56.708674Z",
"iopub.status.idle": "2025-03-23T10:08:05.456742Z",
"shell.execute_reply": "2025-03-23T10:08:05.456022Z",
"shell.execute_reply.started": "2025-03-23T10:07:56.708879Z"
}
},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"\u001b[32m18:07:56.710\u001b[0m | \u001b[33m\u001b[1mWARNING \u001b[0m | \u001b[33m3810563\u001b[0m:\u001b[36mdata\u001b[0m:\u001b[36mconfigure_dataset\u001b[0m - \u001b[33m\u001b[1mOverwriting existing `regime` = \"ctfact\".\u001b[0m\n"
]
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"\u001b[32m18:07:57.178\u001b[0m | \u001b[1mINFO \u001b[0m | \u001b[33m3810563\u001b[0m:\u001b[36mcore\u001b[0m:\u001b[36mpredict_mode\u001b[0m - \u001b[1mNumber of topological generations: [68, 88, 70, 101]\u001b[0m\n"
]
},
{
"data": {
"application/vnd.jupyter.widget-view+json": {
"model_id": "79b9e0e88921485b900ddf4948848c6a",
"version_major": 2,
"version_minor": 0
},
"text/plain": [
"Predicting: | …"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"configure_dataset(source, use_regime=\"interv\")\n",
"nil = cascade.counterfactual(source, design=design, sample=True)"
]
},
{
"cell_type": "code",
"execution_count": 28,
"id": "01fd956d-7b1b-4cb6-9d3d-76a33cb7e889",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T10:08:05.458023Z",
"iopub.status.busy": "2025-03-23T10:08:05.457785Z",
"iopub.status.idle": "2025-03-23T10:08:05.709447Z",
"shell.execute_reply": "2025-03-23T10:08:05.708687Z",
"shell.execute_reply.started": "2025-03-23T10:08:05.458000Z"
}
},
"outputs": [
{
"data": {
"text/plain": [
"AnnData object with n_obs × n_vars = 10000 × 1064\n",
" obs: 'guide_id', 'gemgroup', 'ncounts', 'knockup', 'my_pert', 'role'\n",
" obsm: 'X_pca', 'covariate'\n",
" layers: 'counts', 'interv', 'ctfact', 'X_ctfact'"
]
},
"execution_count": 28,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"combined = ad.concat({\"nil\": nil, \"ctfact\": ctfact}, label=\"role\", index_unique=\"-\")\n",
"combined.X = np.log1p(combined.X * (1e4 / combined.obs[[\"ncounts\"]].to_numpy()))\n",
"combined"
]
},
{
"cell_type": "code",
"execution_count": 29,
"id": "84c8a39b-b967-4aef-b260-4b2e85c4f534",
"metadata": {
"execution": {
"iopub.execute_input": "2025-03-23T10:08:05.711046Z",
"iopub.status.busy": "2025-03-23T10:08:05.710839Z",
"iopub.status.idle": "2025-03-23T10:08:06.062374Z",
"shell.execute_reply": "2025-03-23T10:08:06.061724Z",
"shell.execute_reply.started": "2025-03-23T10:08:05.711028Z"
}
},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"/rd1/user/caozj/CASCADE/conda/lib/python3.11/site-packages/pandas/core/arraylike.py:399: RuntimeWarning: divide by zero encountered in log10\n",
" result = getattr(ufunc, method)(*inputs, **kwargs)\n"
]
},
{
"data": {
"text/html": [
"\n",
"\n",
"
\n",
" \n",
" \n",
" | \n",
" names | \n",
" scores | \n",
" logfoldchanges | \n",
" pvals | \n",
" pvals_adj | \n",
" pct_nz_group | \n",
" -logfdr | \n",
" is_marker | \n",
"
\n",
" \n",
" \n",
" \n",
" | 0 | \n",
" KLF1 | \n",
" 492.331299 | \n",
" 2.114285 | \n",
" 0.0 | \n",
" 0.0 | \n",
" 1.0 | \n",
" 350.0 | \n",
" False | \n",
"
\n",
" \n",
" | 1 | \n",
" PNMT | \n",
" 438.697296 | \n",
" 2.480975 | \n",
" 0.0 | \n",
" 0.0 | \n",
" 1.0 | \n",
" 350.0 | \n",
" False | \n",
"
\n",
" \n",
" | 2 | \n",
" TMSB10 | \n",
" 366.669769 | \n",
" 1.523769 | \n",
" 0.0 | \n",
" 0.0 | \n",
" 1.0 | \n",
" 350.0 | \n",
" False | \n",
"
\n",
" \n",
" | 3 | \n",
" S100A10 | \n",
" 247.035904 | \n",
" 1.381347 | \n",
" 0.0 | \n",
" 0.0 | \n",
" 1.0 | \n",
" 350.0 | \n",
" False | \n",
"
\n",
" \n",
" | 4 | \n",
" HBG2 | \n",
" 236.489166 | \n",
" 1.407301 | \n",
" 0.0 | \n",
" 0.0 | \n",
" 1.0 | \n",
" 350.0 | \n",
" True | \n",
"
\n",
" \n",
"
\n",
"