diff --git a/docs/source/tutorial/notebooks/cbioportal_keeper_demo.ipynb b/docs/source/tutorial/notebooks/cbioportal_keeper_demo.ipynb index 1d133ee..50ec284 100644 --- a/docs/source/tutorial/notebooks/cbioportal_keeper_demo.ipynb +++ b/docs/source/tutorial/notebooks/cbioportal_keeper_demo.ipynb @@ -13,7 +13,7 @@ "3. View attribute and data type summaries for a selected study.\n", "4. Extract molecular data (mutations, CNA, RNA-seq) and load each into netflow's `Keeper` for further analysis.\n", "\n", - "## Requires\n", + "## Requirements\n", "* Python 3.10\n", "* The [Netflow](https://github.com/areElkin/netflow.git) package\n", "\n", @@ -35,8 +35,10 @@ "metadata": {}, "outputs": [], "source": [ + "import numpy as np\n", "import pandas as pd\n", "\n", + "import matplotlib.pyplot as plt\n", "from netflow.prep.extract import CBioPortalClient\n", "from netflow.keepers.keeper import Keeper" ] @@ -51,7 +53,7 @@ }, { "cell_type": "code", - "execution_count": 2, + "execution_count": 3, "id": "ecf7f4fe", "metadata": {}, "outputs": [ @@ -73,12 +75,14 @@ "id": "6626023c", "metadata": {}, "source": [ - "## View all available studies" + "## View all available studies\n", + "\n", + "Get a list of all CBioportal studies, filter them down to studies matching a known keyword, and view a descriptive summary of a selected study. " ] }, { "cell_type": "code", - "execution_count": 3, + "execution_count": 4, "id": "b982cf58", "metadata": {}, "outputs": [ @@ -88,27 +92,17 @@ "text": [ "535 studies available on cBioPortal.\n", "\n", - "First 20 studies:\n", + "First 10 studies:\n", " * Breast Cancer (CPTAC GDC, 2025)\n", " * Glioma (MSK, Clin Cancer Res 2019)\n", " * Osteosarcoma (TARGET GDC, 2025)\n", " * Prostate Cancer MDA PCa PDX (MD Anderson, Clin Cancer Res 2024)\n", " * Uterine Endometrioid Carcinoma (CPTAC GDC, 2025)\n", + " * Kidney Renal Papillary Cell Carcinoma (TCGA, PanCancer Atlas)\n", " * Urothelial Carcinoma (Cornell/Trento, Nat Gen 2016)\n", " * Basal Cell Carcinoma (UNIGE, Nat Genet 2016)\n", " * Endometrial Carcinoma (TCGA GDC, 2025)\n", - " * Kidney Renal Clear Cell Carcinoma (TCGA, PanCancer Atlas)\n", - " * Adenoid Cystic Carcinoma (JHU, Cancer Prev Res 2016)\n", - " * Uterine Corpus Endometrial Carcinoma (TCGA, PanCancer Atlas)\n", - " * Colorectal Adenocarcinoma (TCGA, Firehose Legacy)\n", - " * Lung Adenocarcinoma (CPTAC GDC, 2025)\n", - " * Breast Invasive Carcinoma (Broad, Nature 2012)\n", - " * Retinoblastoma cfDNA (MSK, Cancer Med 2020)\n", - " * Appendiceal Cancer (MSK, J Clin Oncol 2022)\n", - " * Evolution and co-occurrence of PI3K pathway gene mutations in endometrial carcinoma molecular subtypes at the single-cell level\n", - " * Urothelial Cancer Patient Derived Organoids (MSK, 2025)\n", - " * Adrenocortical Carcinoma (TCGA, Firehose Legacy)\n", - " * Ovarian Cancer - MSK SPECTRUM (MSK, Nature 2022)\n" + " * Adenoid Cystic Carcinoma (JHU, Cancer Prev Res 2016)\n" ] } ], @@ -116,14 +110,22 @@ "all_studies = cbio.list_studies()\n", "print(f'{len(all_studies)} studies available on cBioPortal.\\n')\n", "\n", - "print('First 20 studies:')\n", - "for name in list(all_studies)[:20]:\n", + "print('First 10 studies:')\n", + "for name in list(all_studies)[:10]:\n", " print(f' * {name}')" ] }, + { + "cell_type": "markdown", + "id": "c24c8e9e-e472-423f-a5e7-febf2c6e5c52", + "metadata": {}, + "source": [ + "## Search for studies matching keyword" + ] + }, { "cell_type": "code", - "execution_count": 4, + "execution_count": 5, "id": "1803c693", "metadata": {}, "outputs": [ @@ -141,7 +143,6 @@ } ], "source": [ - "# Search for studies by keyword\n", "keyword = 'MSK-IMPACT'\n", "matches = [s for s in all_studies if keyword.lower() in s.lower()]\n", "print(f'Studies matching \"{keyword}\":')\n", @@ -149,6 +150,14 @@ " print(f' {m}')" ] }, + { + "cell_type": "markdown", + "id": "32421483-bbd9-44a5-9ede-99a146e22dc4", + "metadata": {}, + "source": [ + "____" + ] + }, { "cell_type": "markdown", "id": "84d71f27", @@ -159,7 +168,17 @@ }, { "cell_type": "code", - "execution_count": 5, + "execution_count": 6, + "id": "61c0f33a-77c9-46ba-b3c4-985c85e6c451", + "metadata": {}, + "outputs": [], + "source": [ + "pd.set_option('display.max_columns', None)" + ] + }, + { + "cell_type": "code", + "execution_count": 7, "id": "5259d632", "metadata": {}, "outputs": [ @@ -167,23 +186,99 @@ "name": "stdout", "output_type": "stream", "text": [ - "Study ID: msk_impact_2017\n", - "\n", - "Description:\n", - " name description sequencedSampleCount\n", - "MSK-IMPACT Clinical Sequencing Cohort (MSK, Nat Med 2017) Targeted sequencing of 10,000 clinical cases using the MSK-IMPACT assay 10945\n" + "Study ID: msk_impact_2017\n" ] } ], "source": [ "STUDY_NAME = 'MSK-IMPACT Clinical Sequencing Cohort (MSK, Nat Med 2017)'\n", "\n", - "# Study ID and description\n", + "# Get study ID and description\n", "study_id = cbio.get_study_id(STUDY_NAME)\n", - "print(f'Study ID: {study_id}')\n", - "\n", + "print(f'Study ID: {study_id}')" + ] + }, + { + "cell_type": "markdown", + "id": "afa66f65-bfaf-42c1-86b7-eec515e9cc3a", + "metadata": {}, + "source": [ + "### *Study description*" + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "id": "e9837b8b-3daf-4ff6-aeec-a5d27d362a31", + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Study summary:\n" + ] + }, + { + "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", + "
namedescriptioncancerTypeIdreferenceGenomesequencedSampleCount
0MSK-IMPACT Clinical Sequencing Cohort (MSK, Na...Targeted sequencing of 10,000 clinical cases u...mixedhg1910945
\n", + "
" + ], + "text/plain": [ + " name \\\n", + "0 MSK-IMPACT Clinical Sequencing Cohort (MSK, Na... \n", + "\n", + " description cancerTypeId \\\n", + "0 Targeted sequencing of 10,000 clinical cases u... mixed \n", + "\n", + " referenceGenome sequencedSampleCount \n", + "0 hg19 10945 " + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ "desc = cbio.get_study_description(STUDY_NAME)\n", - "print(f'\\nDescription:\\n{desc[[\"name\", \"description\", \"sequencedSampleCount\"]].to_string(index=False)}')" + "print('Study summary:')\n", + "display(desc[['name', 'description', 'cancerTypeId', 'referenceGenome', 'sequencedSampleCount']])" ] }, { @@ -191,12 +286,12 @@ "id": "5b6790f3-8a6d-424f-b3ba-1d13caf2b7c1", "metadata": {}, "source": [ - "### > View all molecular data types available in this study" + "### *List all molecular data profiles available in this study*" ] }, { "cell_type": "code", - "execution_count": 6, + "execution_count": 9, "id": "ef4d589b", "metadata": {}, "outputs": [ @@ -275,12 +370,12 @@ "id": "81208da5-b895-4386-adae-ac0c8a529e9b", "metadata": {}, "source": [ - "### > Cancer-type breakdown" + "### *Cancer-type breakdown*" ] }, { "cell_type": "code", - "execution_count": 7, + "execution_count": 10, "id": "cccc8a64", "metadata": {}, "outputs": [ @@ -400,12 +495,12 @@ "id": "c537c2e5-fd84-49ea-8c54-55c9b737210b", "metadata": {}, "source": [ - "### > Clinical attributes available" + "### *View patient clinical attributes*" ] }, { "cell_type": "code", - "execution_count": 8, + "execution_count": 11, "id": "90ca7b26", "metadata": {}, "outputs": [ @@ -414,14 +509,46 @@ "output_type": "stream", "text": [ "24 clinical attributes available:\n", - "['CANCER_TYPE', 'CANCER_TYPE_DETAILED', 'DNA_INPUT', 'FRACTION_GENOME_ALTERED', 'MATCHED_STATUS', 'METASTATIC_SITE', 'MUTATION_COUNT', 'ONCOTREE_CODE', 'OS_MONTHS', 'OS_STATUS', 'PRIMARY_SITE', 'SAMPLE_CLASS', 'SAMPLE_COLLECTION_SOURCE', 'SAMPLE_COUNT', 'SAMPLE_COVERAGE', 'SAMPLE_TYPE', 'SEX', 'SMOKING_HISTORY', 'SOMATIC_STATUS', 'SPECIMEN_PRESERVATION_TYPE', 'SPECIMEN_TYPE', 'TMB_NONSYNONYMOUS', 'TUMOR_PURITY', 'VITAL_STATUS']\n" + "\n", + "CANCER_TYPE\n", + "CANCER_TYPE_DETAILED\n", + "DNA_INPUT\n", + "FRACTION_GENOME_ALTERED\n", + "MATCHED_STATUS\n", + "METASTATIC_SITE\n", + "MUTATION_COUNT\n", + "ONCOTREE_CODE\n", + "OS_MONTHS\n", + "OS_STATUS\n", + "PRIMARY_SITE\n", + "SAMPLE_CLASS\n", + "SAMPLE_COLLECTION_SOURCE\n", + "SAMPLE_COUNT\n", + "SAMPLE_COVERAGE\n", + "SAMPLE_TYPE\n", + "SEX\n", + "SMOKING_HISTORY\n", + "SOMATIC_STATUS\n", + "SPECIMEN_PRESERVATION_TYPE\n", + "SPECIMEN_TYPE\n", + "TMB_NONSYNONYMOUS\n", + "TUMOR_PURITY\n", + "VITAL_STATUS\n" ] } ], "source": [ "clinical_attrs = cbio.list_study_attributes(STUDY_NAME)\n", - "print(f'{len(clinical_attrs)} clinical attributes available:')\n", - "print(clinical_attrs)" + "print(f'{len(clinical_attrs)} clinical attributes available:\\n')\n", + "print(*clinical_attrs, sep='\\n')" + ] + }, + { + "cell_type": "markdown", + "id": "6e914d7d-4555-4d3f-9db4-a58952830d02", + "metadata": {}, + "source": [ + "----" ] }, { @@ -429,7 +556,20 @@ "id": "5e38ce63", "metadata": {}, "source": [ - "## Extract data" + "## Extract data\n", + "\n", + "#### ***Data types***\n", + " * Patient clinical attributes\n", + " * Mutation \n", + " * CNA \n", + " * RNA-seq expression\n", + " * DNA Methylation\n", + " * Structural variants\n", + "\n", + "Each section below calls the CBioPortalClient method to extract one of the types above and prints a summary of the data. \n", + "*Note*: Data pulls may take a few minutes for large cohorts.\n", + "\n", + "---" ] }, { @@ -437,12 +577,14 @@ "id": "500086fc-f95d-410a-ad57-95c74c913666", "metadata": {}, "source": [ - "### > Clinical data" + "### *Patient clinical attributes*\n", + "\n", + "* `get_clinical_data`() combines sample metadata (primary site, sample type, tumor purity) with patient metadata (OS months, OS status)" ] }, { "cell_type": "code", - "execution_count": 9, + "execution_count": 12, "id": "581e64ce", "metadata": {}, "outputs": [ @@ -500,7 +642,7 @@ " Breast\n", " Primary\n", " 50\n", - " NaN\n", + " <NA>\n", " 0:LIVING\n", " 0.0\n", " \n", @@ -510,7 +652,7 @@ " Breast\n", " Metastasis\n", " 40\n", - " NaN\n", + " <NA>\n", " 1:DECEASED\n", " 1.0\n", " \n", @@ -520,7 +662,7 @@ " Peritoneum\n", " Primary\n", " 30\n", - " 8.71\n", + " <NA>\n", " 1:DECEASED\n", " 1.0\n", " \n", @@ -530,7 +672,7 @@ " Uterus\n", " Metastasis\n", " 40\n", - " 36.75\n", + " <NA>\n", " 0:LIVING\n", " 0.0\n", " \n", @@ -540,7 +682,7 @@ " Uterus\n", " Primary\n", " NaN\n", - " 8.81\n", + " <NA>\n", " 0:LIVING\n", " 0.0\n", " \n", @@ -557,13 +699,13 @@ "P-0000024-T01-IM3 P-0000024 Uterus Metastasis 40 \n", "P-0000025-T01-IM3 P-0000025 Uterus Primary NaN \n", "\n", - " osMonths osStatus osGroup \n", - "sampleId \n", - "P-0000004-T01-IM3 NaN 0:LIVING 0.0 \n", - "P-0000015-T01-IM3 NaN 1:DECEASED 1.0 \n", - "P-0000023-T01-IM3 8.71 1:DECEASED 1.0 \n", - "P-0000024-T01-IM3 36.75 0:LIVING 0.0 \n", - "P-0000025-T01-IM3 8.81 0:LIVING 0.0 " + " osMonths osStatus osGroup \n", + "sampleId \n", + "P-0000004-T01-IM3 0:LIVING 0.0 \n", + "P-0000015-T01-IM3 1:DECEASED 1.0 \n", + "P-0000023-T01-IM3 1:DECEASED 1.0 \n", + "P-0000024-T01-IM3 0:LIVING 0.0 \n", + "P-0000025-T01-IM3 0:LIVING 0.0 " ] }, "metadata": {}, @@ -577,23 +719,32 @@ "display(clin_df.head(5))" ] }, + { + "cell_type": "code", + "execution_count": 13, + "id": "8d1b162b-fd98-48a2-88d8-be17db9a7000", + "metadata": {}, + "outputs": [], + "source": [ + "sample_ids = clin_df.index.to_list()" + ] + }, { "cell_type": "markdown", "id": "377b80d5", "metadata": {}, "source": [ - "### > Extract mutation data and load into Netflow's Keeper\n", + "### *Extract mutation data and load into Netflow's Keeper*\n", "\n", - "* `get_mutation_data()` returns a binary *(genes, samples)* matrix where 1 indicates gene mutations in a sample , 0 otherwise. \n", + "* `get_mutation_data()` returns a binary matrix *(n_genes, n_samples)* where 1 indicates gene mutations in a sample , 0 otherwise. \n", "\n", "*Note*: \n", - "* Data pulls may take a few minutes for a large cohort.\n", - "* Data extraction utilities like `get_mutation_data()` return the raw data (`mut_raw`) and a processed version (`mut_proc`), structured for simplified import to the Keeper. `mut_proc` is a (n_features, n_observations) (or) (n_genes, n_samples) dataframe." + "* Data extraction utilities like `get_mutation_data()` return both the raw data (`mut_raw`) and the processed version (`mut_proc`) that can be readily imported to the Keeper. `mut_proc` is a *(n_features, n_observations)* (or) *(n_genes, n_samples)* dataframe." ] }, { "cell_type": "code", - "execution_count": 10, + "execution_count": 14, "id": "dd0ff708", "metadata": {}, "outputs": [ @@ -721,31 +872,7 @@ }, { "cell_type": "code", - "execution_count": 11, - "id": "3d9e8e9c", - "metadata": {}, - "outputs": [ - { - "name": "stdout", - "output_type": "stream", - "text": [ - "Keeper data keys: ['mutation']\n" - ] - } - ], - "source": [ - "# Initialise Keeper with sample IDs from the mutation matrix\n", - "sample_ids = mut_proc.columns.tolist()\n", - "keeper = Keeper(observation_labels=sample_ids)\n", - "\n", - "# Load mutation data \n", - "keeper.add_data(mut_proc, 'mutation')\n", - "print(f'Keeper data keys: {list(keeper.data.keys())}')" - ] - }, - { - "cell_type": "code", - "execution_count": 12, + "execution_count": 15, "id": "c7ed5675", "metadata": {}, "outputs": [ @@ -860,14 +987,14 @@ "id": "7523d170", "metadata": {}, "source": [ - "### > Extract CNA data and load it into the Keeper\n", + "### *Extract CNA data and load it into the Keeper*\n", "\n", "`get_cna_data()` returns an integer *(genes, samples)* matrix where -2 represents deep deletion, -1 - shallow deletion, 0 - diploid, 1 - shallow gain, 2-amplification." ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 16, "id": "7bbbc797", "metadata": {}, "outputs": [ @@ -875,45 +1002,292 @@ "name": "stdout", "output_type": "stream", "text": [ - "Fetching CNA data...\n" + "Fetching CNA data...\n", + "\n", + "CNA matrix shape: (4487450, 9) (genes \u00d7 samples)\n" ] + }, + { + "data": { + "text/plain": [ + "\tP-0000004-T01-IM3\tP-0000015-T01-IM3\tP-0000023-T01-IM3\tP-0000024-T01-IM3\tP-0000025-T01-IM3\tP-0000025-T02-IM5\tP-0000026-T01-IM3\tP-0000027-T01-IM3\tP-0000030-T01-IM3\tP-0000034-T01-IM3\tP-0000036-T01-IM3\tP-0000037-T01-IM3\tP-0000037-T02-IM3\tP-0000039-T01-IM3\tP-0000041-T01-IM3\tP-0000042-T01-IM3\tP-0000043-T02-IM3\tP-0000047-T01-IM3\tP-0000053-T01-IM3\tP-0000055-T01-IM3\tP-0000056-T01-IM3\tP-0000057-T01-IM3\tP-0000058-T01-IM3\tP-0000059-T01-IM3\tP-0000060-T01-IM3\tP-0000061-T01-IM3\tP-0000062-T01-IM3\tP-0000063-T01-IM3\tP-0000064-T02-IM3\tP-0000065-T01-IM3\tP-0000066-T01-IM3\tP-0000066-T02-IM5\tP-0000067-T01-IM3\tP-0000067-T02-IM5\tP-0000068-T01-IM3\tP-0000069-T01-IM3\tP-0000071-T01-IM3\tP-0000073-T01-IM3\tP-0000075-T01-IM3\tP-0000076-T01-IM3\tP-0000077-T01-IM3\tP-0000079-T01-IM3\tP-0000080-T01-IM3\tP-0000081-T01-IM3\tP-0000082-T01-IM3\tP-0000083-T01-IM3\tP-0000084-T01-IM3\tP-0000085-T02-IM3\tP-0000086-T01-IM3\tP-0000088-T01-IM3\tP-0000092-T01-IM3\tP-0000093-T01-IM3\tP-0000095-T01-IM3\tP-0000096-T01-IM3\tP-0000098-T01-IM3\tP-0000100-T02-IM3\tP-0000102-T02-IM3\tP-0000103-T01-IM3\tP-0000104-T01-IM3\tP-0000105-T01-IM3\tP-0000106-T01-IM3\tP-0000107-T01-IM3\tP-0000108-T01-IM3\tP-0000110-T01-IM3\tP-0000112-T01-IM3\tP-0000113-T01-IM3\tP-0000114-T01-IM3\tP-0000114-T02-IM3\tP-0000115-T01-IM3\tP-0000116-T01-IM3\tP-0000117-T01-IM3\tP-0000118-T01-IM3\tP-0000119-T01-IM3\tP-0000120-T01-IM3\tP-0000121-T01-IM3\tP-0000122-T01-IM3\tP-0000123-T01-IM3\tP-0000125-T01-IM3\tP-0000126-T01-IM3\tP-0000127-T01-IM3\tP-0000127-T02-IM3\tP-0000129-T01-IM3\tP-0000130-T01-IM3\tP-0000131-T01-IM3\tP-0000132-T01-IM3\tP-0000133-T01-IM3\tP-0000134-T01-IM3\tP-0000134-T02-IM3\tP-0000135-T01-IM3\tP-0000136-T01-IM3\tP-0000137-T01-IM3\tP-0000138-T01-IM3\tP-0000138-T02-IM3\tP-0000140-T01-IM3\tP-0000141-T01-IM3\tP-0000142-T01-IM3\tP-0000143-T01-IM3\tP-0000144-T01-IM3\tP-0000145-T01-IM3\tP-0000146-T01-IM3\tP-0000147-T01-IM3\tP-0000148-T01-IM3\tP-0000148-T02-IM3\tP-0000149-T01-IM3\tP-0000149-T02-IM3\tP-0000151-T01-IM3\tP-0000152-T01-IM3\tP-0000152-T02-IM5\tP-0000153-T01-IM3\tP-0000153-T02-IM3\tP-0000154-T01-IM3\tP\n... [truncated \u2014 317,514 chars total]" + ] + }, + "metadata": {}, + "output_type": "display_data" } ], "source": [ "print('Fetching CNA data...')\n", "cna_raw, cna_proc = cbio.get_cna_data(STUDY_NAME)\n", "\n", - "print(f'\\nCNA matrix shape: {cna_raw.shape} (genes × samples)')\n", + "print(f'\\nCNA matrix shape: {cna_raw.shape} (genes \u00d7 samples)')\n", "display(cna_proc.head(5))" ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 17, "id": "5b8aa7c4", "metadata": {}, "outputs": [], + "source": [] + }, + { + "cell_type": "code", + "execution_count": 18, + "id": "310fd59f-4230-423a-9df4-f7611f56f3a1", + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAA90AAAGGCAYAAABmGOKbAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8fJSN1AAAACXBIWXMAAA9hAAAPYQGoP6dpAAA0nUlEQVR4nO3deVRV9f7/8ddBZHA44ATIkhSHSsPyqqk4JUliDjfNuuKQZqQNUBrOWQ5loZim5kB+M9Gb5vBNrbQoQtNukpZmDjmW5njAvgpHSFHh/P7oun+e0NIT2wP4fKy11/J8Pu+9z3uzjixffs7e2+JwOBwCAAAAAABFzsPdDQAAAAAAUFoRugEAAAAAMAmhGwAAAAAAkxC6AQAAAAAwCaEbAAAAAACTELoBAAAAADAJoRsAAAAAAJMQugEAAAAAMAmhGwAAAAAAkxC6AQAAAAAwCaEbAIBS4qefftJTTz2l2rVry8fHR1arVa1atdKMGTN07tw5SVKtWrVksVhksVjk4eEhf39/NWzYUIMGDdLmzZuvelyLxaK4uLhC46+//rosFoueeOIJFRQU6OjRo5owYYKaNWumSpUqqWrVqmrXrp2++OKLQvuOHz9eFotFgYGB+u233wrN16pVS126dLlqP1lZWfLx8ZHFYtGePXtu5EcEAMBNR+gGAKAUWLt2rRo2bKjly5era9eueuutt5SQkKDbbrtNw4cP1+DBg43aRo0a6d///rcWLVqkhIQERURE6OOPP1aLFi0UHx9/Xe83adIkjRkzRv3799c777wjDw8Pffjhh5o8ebLq1q2riRMn6uWXX9bZs2f1wAMPaMGCBVc9TmZmpubOnXtD57pixQpZLBYFBQVp8eLFN7QvAAA3m8XhcDjc3QQAAHDdoUOHdPfdd6tGjRpat26dqlev7jR/8OBBrV27VoMHD1atWrUUFhamNWvWONWcO3dOvXv31urVqzVnzhw988wzxpzFYlFsbKxmzZolSZoyZYpGjBihfv36acGCBfLw+P3/8Hfv3q3AwEBVrVrV2DcvL0+NGjVSTk6Ojh49aoyPHz9eEyZMUKNGjXTy5EkdOnRIvr6+xvy1+pSk++67T1WrVlXNmjW1evVq/fzzz3/jpwcAgLlY6QYAoIRLTExUTk6O5s+fXyhwS1LdunWdVrqvxtfXV//+979VuXJlvfbaa7rW/8lPmzZNI0aMUN++fZ0CtyTdddddToFbkry9vdWpUycdO3ZMZ8+eLXS8sWPHKiMj47pXu48cOaKvvvpK0dHRio6O1qFDh7Rp06br2hcAAHcgdAMAUMJ9/PHHql27tlq2bPm3jlOhQgV1795dx48f148//lhofsaMGRo6dKh69+6t5ORkp8D9Z2w2m8qVK6dy5coVmmvTpo3uv/9+JSYmGted/5n3339f5cuXV5cuXdSsWTPVqVOHr5gDAIo1QjcAACWY3W7X8ePH1bBhwyI5XlhYmKTfb8p2pTVr1mjIkCHq1auXFi1apDJlylzX8Q4ePKiVK1eqR48e19xn3LhxysjIUFJS0l8eb/HixXrooYeMr6L37NlTy5cv16VLl66rHwAAbjZCNwAAJZjdbpckVaxYsUiOV6FCBUkq9FXwjIwMSVJoaOh1B+7ffvtNjz76qHx9fTVp0qRr1rVt21YRERF/udq9Y8cO7dy5U7169TLGevXqpV9//VWfffbZdfUEAMDNRugGAKAEs1qtkgqHZFfl5ORIKhzi+/fvr65du+r111/Xm2+++ZfHyc/PV3R0tH788Uf97//+r4KDg/+0fvz48bLZbH+62v3ee++pfPnyql27tg4ePKiDBw/Kx8dHtWrV4ivmAIBiy9PdDQAAANdZrVYFBwdr165dRXK8y8epW7eu07inp6eWL1+ujh07aujQofL399eAAQOueZyBAwdqzZo1Wrx4se6///6/fN+2bduqXbt2SkxM1NNPP11o3uFw6P3331dubq4aNGhQaD4zM1M5OTnGSj0AAMUFoRsAgBKuS5cumjdvntLT0xUeHu7ycXJycrRq1SqFhISofv36heZ9fHz00UcfKSIiQgMHDpS/v7+6d+9eqG748OFasGCBpk+f7vRV8L8yfvx4tWvXTm+//XahuQ0bNujYsWN65ZVXCvV25swZDRo0SKtXr1bfvn2v+/0AALgZ+Ho5AAAl3IgRI1S+fHk9+eSTxrXXV/rpp580Y8aMPz3GuXPn9Nhjj+n06dMaM2aMLBbLVeusVqtSUlJUt25d9erVS2lpaU7zU6ZM0RtvvKEXX3zxLx9T9kf33Xef2rVrp8mTJ+v8+fNOc5e/Wj58+HA98sgjTtvAgQNVr149vmIOACiWWOkGAKCEq1OnjpYsWaKePXuqfv366tevn8LCwnThwgVt2rRJK1as0OOPP27UHz9+XO+9956k31e3f/zxR61YsUI2m01Dhw7VU0899afvV61aNaWmpqpVq1bq1q2b0tLS1KxZM61atUojRoxQvXr1VL9+feM9LnvggQcUGBj4p8ceN26cIiIinMby8vL0wQcf6IEHHpCPj89V9/vnP/+pGTNmKDMzUwEBAX/6HgAA3EyEbgAASoF//vOf2rFjh6ZMmaIPP/xQc+fOlbe3t+6++25NnTpVAwcONGq3b9+uxx57TBaLRRUrVlRISIi6du2qJ598Us2aNbuu9wsJCdHnn3+uNm3a6MEHH9TGjRv1ww8/SJIOHDigxx57rNA+69ev/8vQ3a5dO913333asGGDMbZ27VplZWWpa9eu19yva9eumjp1qpYuXarnn3/+us4BAICbweJwOBzubgIAAAAAgNKIa7oBAAAAADAJoRsAAAAAAJMQugEAAAAAMAmhGwAAAAAAkxC6AQAAAAAwCaEbAAAAAACT8JzuIlJQUKATJ06oYsWKslgs7m4HAAAAAGAih8Ohs2fPKjg4WB4e117PJnQXkRMnTigkJMTdbQAAAAAAbqKjR4+qRo0a15wndBeRihUrSvr9B261Wt3cDQAAAADATHa7XSEhIUYWvBZCdxG5/JVyq9VK6AYAAACAW8RfXV7MjdQAAAAAADAJoRsAAAAAAJMQugEAAAAAMAmhGwAAAAAAkxC6AQAAAAAwCaEbAAAAAACTELoBAAAAADAJoRsAAAAAAJMQugEAAAAAMAmhGwAAAAAAk7g1dG/cuFFdu3ZVcHCwLBaLVq9e7TTvcDg0duxYVa9eXb6+voqMjNSBAwecak6fPq0+ffrIarXK399fMTExysnJcarZsWOH2rRpIx8fH4WEhCgxMbFQLytWrNCdd94pHx8fNWzYUJ988kmRny8AAAAA4Nbi1tCdm5ure+65R7Nnz77qfGJiombOnKmkpCRt3rxZ5cuXV1RUlM6fP2/U9OnTR7t371ZqaqrWrFmjjRs3atCgQca83W5Xhw4dVLNmTW3dulVTpkzR+PHjNW/ePKNm06ZN6tWrl2JiYvT999+rW7du6tatm3bt2mXeyQMAAAAASj2Lw+FwuLsJSbJYLFq1apW6desm6fdV7uDgYA0dOlTDhg2TJGVnZyswMFDJycmKjo7Wnj171KBBA3377bdq2rSpJCklJUWdOnXSsWPHFBwcrLlz52rMmDGy2Wzy8vKSJI0aNUqrV6/W3r17JUk9e/ZUbm6u1qxZY/TTokULNWrUSElJSdfVv91ul5+fn7Kzs2W1WovqxwIAMFmtUWvd3QIgSTo8qbO7WwAA3IDrzYDF9pruQ4cOyWazKTIy0hjz8/NT8+bNlZ6eLklKT0+Xv7+/EbglKTIyUh4eHtq8ebNR07ZtWyNwS1JUVJT27dunM2fOGDVXvs/lmsvvAwAAAACAKzzd3cC12Gw2SVJgYKDTeGBgoDFns9kUEBDgNO/p6anKlSs71YSGhhY6xuW5SpUqyWaz/en7XE1eXp7y8vKM13a7/UZODwAAAABwCyi2K93FXUJCgvz8/IwtJCTE3S0BAAAAAIqZYhu6g4KCJEkZGRlO4xkZGcZcUFCQMjMzneYvXbqk06dPO9Vc7RhXvse1ai7PX83o0aOVnZ1tbEePHr3RUwQAAAAAlHLFNnSHhoYqKChIaWlpxpjdbtfmzZsVHh4uSQoPD1dWVpa2bt1q1Kxbt04FBQVq3ry5UbNx40ZdvHjRqElNTdUdd9yhSpUqGTVXvs/lmsvvczXe3t6yWq1OGwAAAAAAV3Jr6M7JydH27du1fft2Sb/fPG379u06cuSILBaLhgwZookTJ+qjjz7Szp071a9fPwUHBxt3OK9fv746duyogQMHasuWLfr6668VFxen6OhoBQcHS5J69+4tLy8vxcTEaPfu3Vq2bJlmzJih+Ph4o4/BgwcrJSVFU6dO1d69ezV+/Hh99913iouLu9k/EgAAAABAKeLWG6l99913ioiIMF5fDsL9+/dXcnKyRowYodzcXA0aNEhZWVlq3bq1UlJS5OPjY+yzePFixcXFqX379vLw8FCPHj00c+ZMY97Pz0+ff/65YmNj1aRJE1WtWlVjx451epZ3y5YttWTJEr300kt68cUXVa9ePa1evVphYWE34acAAAAAACitis1zuks6ntMNACUTz+lGccFzugGgZCnxz+kGAAAAAKCkI3QDAAAAAGASQjcAAAAAACYhdAMAAAAAYBJCNwAAAAAAJiF0AwAAAABgEkI3AAAAAAAmIXQDAAAAAGASQjcAAAAAACYhdAMAAAAAYBJCNwAAAAAAJiF0AwAAAABgEkI3AAAAAAAmIXQDAAAAAGASQjcAAAAAACYhdAMAAAAAYBJCNwAAAAAAJiF0AwAAAABgEkI3AAAAAAAmIXQDAAAAAGASQjcAAAAAACYhdAMAAAAAYBJCNwAAAAAAJiF0AwAAAABgEkI3AAAAAAAmIXQDAAAAAGASQjcAAAAAACYhdAMAAAAAYBJCNwAAAAAAJiF0AwAAAABgEkI3AAAAAAAmIXQDAAAAAGASQjcAAAAAACYhdAMAAAAAYBJCNwAAAAAAJiF0AwAAAABgEkI3AAAAAAAmIXQDAAAAAGASQjcAAAAAACYhdAMAAAAAYBJCNwAAAAAAJiF0AwAAAABgEkI3AAAAAAAmIXQDAAAAAGASQjcAAAAAACYhdAMAAAAAYBJCNwAAAAAAJinWoTs/P18vv/yyQkND5evrqzp16ujVV1+Vw+EwahwOh8aOHavq1avL19dXkZGROnDggNNxTp8+rT59+shqtcrf318xMTHKyclxqtmxY4fatGkjHx8fhYSEKDEx8aacIwAAAACg9CrWoXvy5MmaO3euZs2apT179mjy5MlKTEzUW2+9ZdQkJiZq5syZSkpK0ubNm1W+fHlFRUXp/PnzRk2fPn20e/dupaamas2aNdq4caMGDRpkzNvtdnXo0EE1a9bU1q1bNWXKFI0fP17z5s27qecLAAAAAChdLI4rl42LmS5duigwMFDz5883xnr06CFfX1+99957cjgcCg4O1tChQzVs2DBJUnZ2tgIDA5WcnKzo6Gjt2bNHDRo00LfffqumTZtKklJSUtSpUycdO3ZMwcHBmjt3rsaMGSObzSYvLy9J0qhRo7R69Wrt3bv3unq12+3y8/NTdna2rFZrEf8kAABmqTVqrbtbACRJhyd1dncLAIAbcL0ZsFivdLds2VJpaWnav3+/JOmHH37Qf/7zHz344IOSpEOHDslmsykyMtLYx8/PT82bN1d6erokKT09Xf7+/kbglqTIyEh5eHho8+bNRk3btm2NwC1JUVFR2rdvn86cOXPV3vLy8mS32502AAAAAACu5OnuBv7MqFGjZLfbdeedd6pMmTLKz8/Xa6+9pj59+kiSbDabJCkwMNBpv8DAQGPOZrMpICDAad7T01OVK1d2qgkNDS10jMtzlSpVKtRbQkKCJkyYUARnCQAAAAAorYr1Svfy5cu1ePFiLVmyRNu2bdPChQv1xhtvaOHChe5uTaNHj1Z2draxHT161N0tAQAAAACKmWK90j18+HCNGjVK0dHRkqSGDRvql19+UUJCgvr376+goCBJUkZGhqpXr27sl5GRoUaNGkmSgoKClJmZ6XTcS5cu6fTp08b+QUFBysjIcKq5/PpyzR95e3vL29v7758kAAAAAKDUKtYr3b/99ps8PJxbLFOmjAoKCiRJoaGhCgoKUlpamjFvt9u1efNmhYeHS5LCw8OVlZWlrVu3GjXr1q1TQUGBmjdvbtRs3LhRFy9eNGpSU1N1xx13XPWr5QAAAAAAXI9iHbq7du2q1157TWvXrtXhw4e1atUqTZs2Td27d5ckWSwWDRkyRBMnTtRHH32knTt3ql+/fgoODla3bt0kSfXr11fHjh01cOBAbdmyRV9//bXi4uIUHR2t4OBgSVLv3r3l5eWlmJgY7d69W8uWLdOMGTMUHx/vrlMHAAAAAJQCxfrr5W+99ZZefvllPfvss8rMzFRwcLCeeuopjR071qgZMWKEcnNzNWjQIGVlZal169ZKSUmRj4+PUbN48WLFxcWpffv28vDwUI8ePTRz5kxj3s/PT59//rliY2PVpEkTVa1aVWPHjnV6ljcAAAAAADeqWD+nuyThOd0AUDLxnG4UFzynGwBKllLxnG4AAAAAAEoyQjcAAAAAACYhdAMAAAAAYBJCNwAAAAAAJiF0AwAAAABgEkI3AAAAAAAmIXQDAAAAAGASQjcAAAAAACYhdAMAAAAAYBJCNwAAAAAAJiF0AwAAAABgEkI3AAAAAAAmIXQDAAAAAGASQjcAAAAAACYhdAMAAAAAYBJCNwAAAAAAJiF0AwAAAABgEkI3AAAAAAAmIXQDAAAAAGASQjcAAAAAACYhdAMAAAAAYBJCNwAAAAAAJiF0AwAAAABgEkI3AAAAAAAmIXQDAAAAAGASQjcAAAAAACYhdAMAAAAAYBJCNwAAAAAAJiF0AwAAAABgEkI3AAAAAAAmIXQDAAAAAGASQjcAAAAAACYhdAMAAAAAYBJCNwAAAAAAJiF0AwAAAABgEkI3AAAAAAAmIXQDAAAAAGASl0L3zz//XNR9AAAAAABQ6rgUuuvWrauIiAi99957On/+fFH3BAAAAABAqeBS6N62bZvuvvtuxcfHKygoSE899ZS2bNlS1L0BAAAAAFCiuRS6GzVqpBkzZujEiRN69913dfLkSbVu3VphYWGaNm2aTp06VdR9AgAAAABQ4vytG6l5enrq4Ycf1ooVKzR58mQdPHhQw4YNU0hIiPr166eTJ08WVZ8AAAAAAJQ4fyt0f/fdd3r22WdVvXp1TZs2TcOGDdNPP/2k1NRUnThxQg899FBR9QkAAAAAQInj6cpO06ZN04IFC7Rv3z516tRJixYtUqdOneTh8XuGDw0NVXJysmrVqlWUvQIAAAAAUKK4FLrnzp2rJ554Qo8//riqV69+1ZqAgADNnz//bzUHAAAAAEBJ5lLoPnDgwF/WeHl5qX///q4cHgAAAACAUsGla7oXLFigFStWFBpfsWKFFi5c+LebutLx48fVt29fValSRb6+vmrYsKG+++47Y97hcGjs2LGqXr26fH19FRkZWeg/BU6fPq0+ffrIarXK399fMTExysnJcarZsWOH2rRpIx8fH4WEhCgxMbFIzwMAAAAAcOtxKXQnJCSoatWqhcYDAgL0+uuv/+2mLjtz5oxatWqlsmXL6tNPP9WPP/6oqVOnqlKlSkZNYmKiZs6cqaSkJG3evFnly5dXVFSUzp8/b9T06dNHu3fvVmpqqtasWaONGzdq0KBBxrzdbleHDh1Us2ZNbd26VVOmTNH48eM1b968IjsXAAAAAMCtx+JwOBw3upOPj4/27t1b6EZphw8fVv369XXu3LkiaW7UqFH6+uuv9dVXX1113uFwKDg4WEOHDtWwYcMkSdnZ2QoMDFRycrKio6O1Z88eNWjQQN9++62aNm0qSUpJSVGnTp107NgxBQcHa+7cuRozZoxsNpu8vLyM9169erX27t17Xb3a7Xb5+fkpOztbVqu1CM4eAHAz1Bq11t0tAJKkw5M6u7sFAMANuN4M6NJKd0BAgHbs2FFo/IcfflCVKlVcOeRVffTRR2ratKkeffRRBQQE6B//+If+53/+x5g/dOiQbDabIiMjjTE/Pz81b95c6enpkqT09HT5+/sbgVuSIiMj5eHhoc2bNxs1bdu2NQK3JEVFRWnfvn06c+ZMkZ0PAAAAAODW4lLo7tWrl55//nmtX79e+fn5ys/P17p16zR48GBFR0cXWXM///yz5s6dq3r16umzzz7TM888o+eff964btxms0mSAgMDnfYLDAw05mw2mwICApzmPT09VblyZaeaqx3jyvf4o7y8PNntdqcNAAAAAIAruXT38ldffVWHDx9W+/bt5en5+yEKCgrUr1+/Ir2mu6CgQE2bNjWO+Y9//EO7du1SUlKS2++MnpCQoAkTJri1BwAAAABA8ebSSreXl5eWLVumvXv3avHixVq5cqV++uknvfvuu05f0f67qlevrgYNGjiN1a9fX0eOHJEkBQUFSZIyMjKcajIyMoy5oKAgZWZmOs1funRJp0+fdqq52jGufI8/Gj16tLKzs43t6NGjrpwiAAAAAKAUcyl0X3b77bfr0UcfVZcuXVSzZs2i6snQqlUr7du3z2ls//79xnuFhoYqKChIaWlpxrzdbtfmzZsVHh4uSQoPD1dWVpa2bt1q1Kxbt04FBQVq3ry5UbNx40ZdvHjRqElNTdUdd9zhdKf0K3l7e8tqtTptAAAAAABcyaWvl+fn5ys5OVlpaWnKzMxUQUGB0/y6deuKpLkXXnhBLVu21Ouvv65//etf2rJli+bNm2c8ystisWjIkCGaOHGi6tWrp9DQUL388ssKDg5Wt27dJP2+Mt6xY0cNHDhQSUlJunjxouLi4hQdHa3g4GBJUu/evTVhwgTFxMRo5MiR2rVrl2bMmKE333yzSM4DAAAAAHBrcil0Dx48WMnJyercubPCwsJksViKui9J0r333qtVq1Zp9OjReuWVVxQaGqrp06erT58+Rs2IESOUm5urQYMGKSsrS61bt1ZKSop8fHyMmsWLFysuLk7t27eXh4eHevTooZkzZxrzfn5++vzzzxUbG6smTZqoatWqGjt2rNOzvAEAAAAAuFEuPae7atWqWrRokTp16mRGTyUSz+kGgJKJ53SjuOA53QBQspj6nG4vLy/VrVvX5eYAAAAAALgVuBS6hw4dqhkzZsiFRXIAAAAAAG4ZLl3T/Z///Efr16/Xp59+qrvuuktly5Z1ml+5cmWRNAcAAAAAQEnmUuj29/dX9+7di7oXAAAAAABKFZdC94IFC4q6DwAAAAAASh2XrumWpEuXLumLL77Q22+/rbNnz0qSTpw4oZycnCJrDgAAAACAksylle5ffvlFHTt21JEjR5SXl6cHHnhAFStW1OTJk5WXl6ekpKSi7hMAAAAAgBLHpZXuwYMHq2nTpjpz5ox8fX2N8e7duystLa3ImgMAAAAAoCRzaaX7q6++0qZNm+Tl5eU0XqtWLR0/frxIGgMAAAAAoKRzaaW7oKBA+fn5hcaPHTumihUr/u2mAAAAAAAoDVwK3R06dND06dON1xaLRTk5ORo3bpw6depUVL0BAAAAAFCiufT18qlTpyoqKkoNGjTQ+fPn1bt3bx04cEBVq1bV+++/X9Q9AgAAAABQIrkUumvUqKEffvhBS5cu1Y4dO5STk6OYmBj16dPH6cZqAAAAAADcylwK3ZLk6empvn37FmUvAAAAAACUKi6F7kWLFv3pfL9+/VxqBgAAAACA0sSl0D148GCn1xcvXtRvv/0mLy8vlStXjtANAAAAAIBcvHv5mTNnnLacnBzt27dPrVu35kZqAAAAAAD8l0uh+2rq1aunSZMmFVoFBwAAAADgVlVkoVv6/eZqJ06cKMpDAgAAAABQYrl0TfdHH33k9NrhcOjkyZOaNWuWWrVqVSSNAQAAAABQ0rkUurt16+b02mKxqFq1arr//vs1derUougLAAAAAIASz6XQXVBQUNR9AAAAAABQ6hTpNd0AAAAAAOD/c2mlOz4+/rprp02b5spbAAAAAABQ4rkUur///nt9//33unjxou644w5J0v79+1WmTBk1btzYqLNYLEXTJQAAAAAAJZBLobtr166qWLGiFi5cqEqVKkmSzpw5owEDBqhNmzYaOnRokTYJAAAAAEBJ5NI13VOnTlVCQoIRuCWpUqVKmjhxIncvBwAAAADgv1wK3Xa7XadOnSo0furUKZ09e/ZvNwUAAAAAQGngUuju3r27BgwYoJUrV+rYsWM6duyYPvjgA8XExOjhhx8u6h4BAAAAACiRXLqmOykpScOGDVPv3r118eLF3w/k6amYmBhNmTKlSBsEAAAAAKCkcil0lytXTnPmzNGUKVP0008/SZLq1Kmj8uXLF2lzAAAAAACUZC59vfyykydP6uTJk6pXr57Kly8vh8NRVH0BAAAAAFDiuRS6/+///k/t27fX7bffrk6dOunkyZOSpJiYGB4XBgAAAADAf7kUul944QWVLVtWR44cUbly5Yzxnj17KiUlpciaAwAAAACgJHPpmu7PP/9cn332mWrUqOE0Xq9ePf3yyy9F0hgAAAAAACWdSyvdubm5Tivcl50+fVre3t5/uykAAAAAAEoDl0J3mzZttGjRIuO1xWJRQUGBEhMTFRERUWTNAQAAAABQkrn09fLExES1b99e3333nS5cuKARI0Zo9+7dOn36tL7++uui7hEAAAAAgBLJpZXusLAw7d+/X61bt9ZDDz2k3NxcPfzww/r+++9Vp06dou4RAAAAAIAS6YZXui9evKiOHTsqKSlJY8aMMaMnAAAAAABKhRte6S5btqx27NhhRi8AAAAAAJQqLn29vG/fvpo/f35R9wIAAAAAQKni0o3ULl26pHfffVdffPGFmjRpovLlyzvNT5s2rUiaAwAAAACgJLuh0P3zzz+rVq1a2rVrlxo3bixJ2r9/v1ONxWIpuu4AAAAAACjBbih016tXTydPntT69eslST179tTMmTMVGBhoSnMAAAAAAJRkN3RNt8PhcHr96aefKjc3t0gbAgAAAACgtHDpRmqX/TGEm23SpEmyWCwaMmSIMXb+/HnFxsaqSpUqqlChgnr06KGMjAyn/Y4cOaLOnTurXLlyCggI0PDhw3Xp0iWnmi+//FKNGzeWt7e36tatq+Tk5JtwRgAAAACA0uyGQrfFYil0zfbNuob722+/1dtvv627777bafyFF17Qxx9/rBUrVmjDhg06ceKEHn74YWM+Pz9fnTt31oULF7Rp0yYtXLhQycnJGjt2rFFz6NAhde7cWREREdq+fbuGDBmiJ598Up999tlNOTcAAAAAQOlkcdzAcrWHh4cefPBBeXt7S5I+/vhj3X///YXuXr5y5coibTInJ0eNGzfWnDlzNHHiRDVq1EjTp09Xdna2qlWrpiVLluiRRx6RJO3du1f169dXenq6WrRooU8//VRdunTRiRMnjGvPk5KSNHLkSJ06dUpeXl4aOXKk1q5dq127dhnvGR0draysLKWkpFxXj3a7XX5+fsrOzpbVai3S8wcAmKfWqLXubgGQJB2e1NndLQAAbsD1ZsAbWunu37+/AgIC5OfnJz8/P/Xt21fBwcHG68tbUYuNjVXnzp0VGRnpNL5161ZdvHjRafzOO+/UbbfdpvT0dElSenq6GjZs6HSzt6ioKNntdu3evduo+eOxo6KijGMAAAAAAOCKG7p7+YIFC8zq45qWLl2qbdu26dtvvy00Z7PZ5OXlJX9/f6fxwMBA2Ww2o+aPd1e//Pqvaux2u86dOydfX99C752Xl6e8vDzjtd1uv/GTAwAAAACUan/rRmpmO3r0qAYPHqzFixfLx8fH3e04SUhIcFrdDwkJcXdLAAAAAIBipliH7q1btyozM1ONGzeWp6enPD09tWHDBs2cOVOenp4KDAzUhQsXlJWV5bRfRkaGgoKCJElBQUGF7mZ++fVf1Vit1quuckvS6NGjlZ2dbWxHjx4tilMGAAAAAJQixTp0t2/fXjt37tT27duNrWnTpurTp4/x57JlyyotLc3YZ9++fTpy5IjCw8MlSeHh4dq5c6cyMzONmtTUVFmtVjVo0MCoufIYl2suH+NqvL29ZbVanTYAAAAAAK50Q9d032wVK1ZUWFiY01j58uVVpUoVYzwmJkbx8fGqXLmyrFarnnvuOYWHh6tFixaSpA4dOqhBgwZ67LHHlJiYKJvNppdeekmxsbHGXdiffvppzZo1SyNGjNATTzyhdevWafny5Vq7ljvaAgAAAABcV6xD9/V488035eHhoR49eigvL09RUVGaM2eOMV+mTBmtWbNGzzzzjMLDw1W+fHn1799fr7zyilETGhqqtWvX6oUXXtCMGTNUo0YNvfPOO4qKinLHKQEAAAAASokbek43ro3ndANAycRzulFc8JxuAChZTHlONwAAAAAAuH6EbgAAAAAATELoBgAAAADAJIRuAAAAAABMQugGAAAAAMAkhG4AAAAAAExC6AYAAAAAwCSEbgAAAAAATELoBgAAAADAJIRuAAAAAABMQugGAAAAAMAkhG4AAAAAAExC6AYAAAAAwCSEbgAAAAAATELoBgAAAADAJIRuAAAAAABMQugGAAAAAMAkhG4AAAAAAExC6AYAAAAAwCSEbgAAAAAATELoBgAAAADAJIRuAAAAAABMQugGAAAAAMAkhG4AAAAAAExC6AYAAAAAwCSEbgAAAAAATELoBgAAAADAJIRuAAAAAABMQugGAAAAAMAkhG4AAAAAAExC6AYAAAAAwCSEbgAAAAAATELoBgAAAADAJIRuAAAAAABMQugGAAAAAMAkhG4AAAAAAExC6AYAAAAAwCSEbgAAAAAATELoBgAAAADAJIRuAAAAAABMQugGAAAAAMAkhG4AAAAAAExC6AYAAAAAwCSEbgAAAAAATELoBgAAAADAJIRuAAAAAABMQugGAAAAAMAkxTp0JyQk6N5771XFihUVEBCgbt26ad++fU4158+fV2xsrKpUqaIKFSqoR48eysjIcKo5cuSIOnfurHLlyikgIEDDhw/XpUuXnGq+/PJLNW7cWN7e3qpbt66Sk5PNPj0AAAAAQClXrEP3hg0bFBsbq2+++Uapqam6ePGiOnTooNzcXKPmhRde0Mcff6wVK1Zow4YNOnHihB5++GFjPj8/X507d9aFCxe0adMmLVy4UMnJyRo7dqxRc+jQIXXu3FkRERHavn27hgwZoieffFKfffbZTT1fAAAAAEDpYnE4HA53N3G9Tp06pYCAAG3YsEFt27ZVdna2qlWrpiVLluiRRx6RJO3du1f169dXenq6WrRooU8//VRdunTRiRMnFBgYKElKSkrSyJEjderUKXl5eWnkyJFau3atdu3aZbxXdHS0srKylJKScl292e12+fn5KTs7W1artehPHgBgilqj1rq7BUCSdHhSZ3e3AAC4AdebAYv1SvcfZWdnS5IqV64sSdq6dasuXryoyMhIo+bOO+/UbbfdpvT0dElSenq6GjZsaARuSYqKipLdbtfu3buNmiuPcbnm8jGuJi8vT3a73WkDAAAAAOBKJSZ0FxQUaMiQIWrVqpXCwsIkSTabTV5eXvL393eqDQwMlM1mM2quDNyX5y/P/VmN3W7XuXPnrtpPQkKC/Pz8jC0kJORvnyMAAAAAoHQpMaE7NjZWu3bt0tKlS93diiRp9OjRys7ONrajR4+6uyUAAAAAQDHj6e4GrkdcXJzWrFmjjRs3qkaNGsZ4UFCQLly4oKysLKfV7oyMDAUFBRk1W7ZscTre5bubX1nzxzueZ2RkyGq1ytfX96o9eXt7y9vb+2+fGwAAAACg9CrWK90Oh0NxcXFatWqV1q1bp9DQUKf5Jk2aqGzZskpLSzPG9u3bpyNHjig8PFySFB4erp07dyozM9OoSU1NldVqVYMGDYyaK49xuebyMQAAAAAAcEWxXumOjY3VkiVL9OGHH6pixYrGNdh+fn7y9fWVn5+fYmJiFB8fr8qVK8tqteq5555TeHi4WrRoIUnq0KGDGjRooMcee0yJiYmy2Wx66aWXFBsba6xUP/3005o1a5ZGjBihJ554QuvWrdPy5cu1di13tAUAAAAAuK5Yr3TPnTtX2dnZateunapXr25sy5YtM2refPNNdenSRT169FDbtm0VFBSklStXGvNlypTRmjVrVKZMGYWHh6tv377q16+fXnnlFaMmNDRUa9euVWpqqu655x5NnTpV77zzjqKiom7q+QIAAAAASpcS9Zzu4ozndANAycRzulFc8JxuAChZSuVzugEAAAAAKEkI3QAAAAAAmITQDQAAAACASQjdAAAAAACYhNANAAAAAIBJCN0AAAAAAJiE0A0AAAAAgEkI3QAAAAAAmITQDQAAAACASQjdAAAAAACYhNANAAAAAIBJCN0AAAAAAJiE0A0AAAAAgEkI3QAAAAAAmITQDQAAAACASQjdAAAAAACYhNANAAAAAIBJCN0AAAAAAJiE0A0AAAAAgEkI3QAAAAAAmITQDQAAAACASQjdAAAAAACYhNANAAAAAIBJCN0AAAAAAJiE0A0AAAAAgEkI3QAAAAAAmITQDQAAAACASTzd3QBurlqj1rq7BUCSdHhSZ3e3AAAAAJiOlW4AAAAAAExC6AYAAAAAwCSEbgAAAAAATELoBgAAAADAJIRuAAAAAABMQugGAAAAAMAkhG4AAAAAAExC6AYAAAAAwCSEbgAAAAAATELoBgAAAADAJIRuAAAAAABMQugGAAAAAMAkhG4AAAAAAExC6AYAAAAAwCSEbgAAAAAATELoBgAAAADAJIRuAAAAAABMQugGAAAAAMAkhO4/mD17tmrVqiUfHx81b95cW7ZscXdLAAAAAIASitB9hWXLlik+Pl7jxo3Ttm3bdM899ygqKkqZmZnubg0AAAAAUAIRuq8wbdo0DRw4UAMGDFCDBg2UlJSkcuXK6d1333V3awAAAACAEojQ/V8XLlzQ1q1bFRkZaYx5eHgoMjJS6enpbuwMAAAAAFBSebq7geLi119/VX5+vgIDA53GAwMDtXfv3kL1eXl5ysvLM15nZ2dLkux2u7mN/k0Feb+5uwVAUvH/u4JbB78XUVzwexEASpbLv7cdDsef1hG6XZSQkKAJEyYUGg8JCXFDN0DJ4zfd3R0AQPHC70UAKJnOnj0rPz+/a84Tuv+ratWqKlOmjDIyMpzGMzIyFBQUVKh+9OjRio+PN14XFBTo9OnTqlKliiwWi+n9wj3sdrtCQkJ09OhRWa1Wd7eDWxyfRxQXfBZRXPBZRHHBZ/HW4HA4dPbsWQUHB/9pHaH7v7y8vNSkSROlpaWpW7dukn4P0mlpaYqLiytU7+3tLW9vb6cxf3//m9ApigOr1covUBQbfB5RXPBZRHHBZxHFBZ/F0u/PVrgvI3RfIT4+Xv3791fTpk3VrFkzTZ8+Xbm5uRowYIC7WwMAAAAAlECE7iv07NlTp06d0tixY2Wz2dSoUSOlpKQUurkaAAAAAADXg9D9B3FxcVf9Ojkg/X5Zwbhx4wpdWgC4A59HFBd8FlFc8FlEccFnEVeyOP7q/uYAAAAAAMAlHu5uAAAAAACA0orQDQAAAACASQjdAAAAAACYhNANuOjw4cOKiYlRaGiofH19VadOHY0bN04XLlxwd2u4Bb322mtq2bKlypUrJ39/f3e3g1vI7NmzVatWLfn4+Kh58+basmWLu1vCLWjjxo3q2rWrgoODZbFYtHr1ane3hFtQQkKC7r33XlWsWFEBAQHq1q2b9u3b5+62UAwQugEX7d27VwUFBXr77be1e/duvfnmm0pKStKLL77o7tZwC7pw4YIeffRRPfPMM+5uBbeQZcuWKT4+XuPGjdO2bdt0zz33KCoqSpmZme5uDbeY3Nxc3XPPPZo9e7a7W8EtbMOGDYqNjdU333yj1NRUXbx4UR06dFBubq67W4ObcfdyoAhNmTJFc+fO1c8//+zuVnCLSk5O1pAhQ5SVleXuVnALaN68ue69917NmjVLklRQUKCQkBA999xzGjVqlJu7w63KYrFo1apV6tatm7tbwS3u1KlTCggI0IYNG9S2bVt3twM3YqUbKELZ2dmqXLmyu9sAANNduHBBW7duVWRkpDHm4eGhyMhIpaenu7EzACgesrOzJYl/G4LQDRSVgwcP6q233tJTTz3l7lYAwHS//vqr8vPzFRgY6DQeGBgom83mpq4AoHgoKCjQkCFD1KpVK4WFhbm7HbgZoRv4g1GjRslisfzptnfvXqd9jh8/ro4dO+rRRx/VwIED3dQ5ShtXPosAAMD9YmNjtWvXLi1dutTdraAY8HR3A0BxM3ToUD3++ON/WlO7dm3jzydOnFBERIRatmypefPmmdwdbiU3+lkEbqaqVauqTJkyysjIcBrPyMhQUFCQm7oCAPeLi4vTmjVrtHHjRtWoUcPd7aAYIHQDf1CtWjVVq1btumqPHz+uiIgINWnSRAsWLJCHB18eQdG5kc8icLN5eXmpSZMmSktLM25YVVBQoLS0NMXFxbm3OQBwA4fDoeeee06rVq3Sl19+qdDQUHe3hGKC0A246Pjx42rXrp1q1qypN954Q6dOnTLmWOXBzXbkyBGdPn1aR44cUX5+vrZv3y5Jqlu3ripUqODe5lBqxcfHq3///mratKmaNWum6dOnKzc3VwMGDHB3a7jF5OTk6ODBg8brQ4cOafv27apcubJuu+02N3aGW0lsbKyWLFmiDz/8UBUrVjTub+Hn5ydfX183dwd34pFhgIuSk5Ov+Q9L/lrhZnv88ce1cOHCQuPr169Xu3btbn5DuGXMmjVLU6ZMkc1mU6NGjTRz5kw1b97c3W3hFvPll18qIiKi0Hj//v2VnJx88xvCLclisVx1fMGCBX95uRhKN0I3AAAAAAAm4QJUAAAAAABMQugGAAAAAMAkhG4AAAAAAExC6AYAAAAAwCSEbgAAAAAATELoBgAAAADAJIRuAAAAAABMQugGAAAAAMAkhG4AAAAAAExC6AYAAIXYbDY999xzql27try9vRUSEqKuXbsqLS1NklSrVi1ZLBZ98803TvsNGTJE7dq1K3S8Y8eOycvLS2FhYTejfQAAig1CNwAAcHL48GE1adJE69at05QpU7Rz506lpKQoIiJCsbGxRp2Pj49Gjhx5XcdMTk7Wv/71L9ntdm3evNms1gEAKHY83d0AAAAoXp599llZLBZt2bJF5cuXN8bvuusuPfHEE8brQYMGKSkpSZ988ok6dep0zeM5HA4tWLBAc+bMUY0aNTR//nw1b97c1HMAAKC4YKUbAAAYTp8+rZSUFMXGxjoF7sv8/f2NP4eGhurpp5/W6NGjVVBQcM1jrl+/Xr/99psiIyPVt29fLV26VLm5uWa0DwBAsUPoBgAAhoMHD8rhcOjOO++8rvqXXnpJhw4d0uLFi69ZM3/+fEVHR6tMmTIKCwtT7dq1tWLFiqJqGQCAYo3QDQAADA6H44bqq1WrpmHDhmns2LG6cOFCofmsrCytXLlSffv2Ncb69u2r+fPn/+1eAQAoCbimGwAAGOrVqyeLxaK9e/de9z7x8fGaM2eO5syZU2huyZIlOn/+vNM13A6HQwUFBdq/f79uv/32IukbAIDiipVuAABgqFy5sqKiojR79uyrXnedlZVVaKxChQp6+eWX9dprr+ns2bNOc/Pnz9fQoUO1fft2Y/vhhx/Upk0bvfvuu2adBgAAxQahGwAAOJk9e7by8/PVrFkzffDBBzpw4ID27NmjmTNnKjw8/Kr7DBo0SH5+flqyZIkxtn37dm3btk1PPvmkwsLCnLZevXpp4cKFunTp0s06LQAA3ILQDQAAnNSuXVvbtm1TRESEhg4dqrCwMD3wwANKS0vT3Llzr7pP2bJl9eqrr+r8+fPG2Pz589WgQYOr3pSte/fuyszM1CeffGLaeQAAUBxYHDd6xxQAAAAAAHBdWOkGAAAAAMAkhG4AAAAAAExC6AYAAAAAwCSEbgAAAAAATELoBgAAAADAJIRuAAAAAABMQugGAAAAAMAkhG4AAAAAAExC6AYAAAAAwCSEbgAAAAAATELoBgAAAADAJIRuAAAAAABM8v8ASzzff9Qze/QAAAAASUVORK5CYII=", + "text/plain": [ + "
" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "fig, ax = plt.subplots(figsize=(10, 4))\n", + "\n", + "alterFreq = np.argmax((cna_proc!=0).sum(axis=1))\n", + "topAltered = cna_proc.index[alterFreq]\n", + "geneAlt = cna_proc.loc[topAltered]\n", + "bins = np.arange(geneAlt.min(), geneAlt.max() + 2) - 0.5\n", + "geneAlt.hist(ax=ax, bins=bins, grid=False)\n", + "\n", + "ax.set_ylabel('Frequency')\n", + "ax.set_xlabel('CNA')\n", + "ax.set_title('CDK2NA')\n", + "plt.tight_layout()\n", + "plt.show()" + ] + }, + { + "cell_type": "markdown", + "id": "9b48de3e-4050-489e-b4bc-97a1275f21e6", + "metadata": {}, + "source": [ + "## *Extract RNA seq data*" + ] + }, + { + "cell_type": "code", + "execution_count": 19, + "id": "a85014a9-11fa-40ec-a731-03730fb32068", + "metadata": {}, + "outputs": [], + "source": [ + "STUDY_NAME_2 = 'Breast Invasive Carcinoma (TCGA, PanCancer Atlas)'\n", + "rna_types, rna_ids = cbio.list_study_data(STUDY_NAME_2)\n", + "rna_raw, rna_proc = cbio.get_rna_seq_data(STUDY_NAME_2)" + ] + }, + { + "cell_type": "code", + "execution_count": 20, + "id": "eaacc165-2566-422a-b83a-d227fb2fc629", + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "RNA-seq df shape: (20338, 1082)\n" + ] + }, + { + "data": { + "text/plain": [ + "\tTCGA-3C-AAAU-01\tTCGA-3C-AALI-01\tTCGA-3C-AALJ-01\tTCGA-3C-AALK-01\tTCGA-4H-AAAK-01\tTCGA-5L-AAT0-01\tTCGA-5T-A9QA-01\tTCGA-A1-A0SB-01\tTCGA-A1-A0SD-01\tTCGA-A1-A0SE-01\tTCGA-A1-A0SF-01\tTCGA-A1-A0SG-01\tTCGA-A1-A0SH-01\tTCGA-A1-A0SI-01\tTCGA-A1-A0SJ-01\tTCGA-A1-A0SK-01\tTCGA-A1-A0SM-01\tTCGA-A1-A0SN-01\tTCGA-A1-A0SO-01\tTCGA-A1-A0SP-01\tTCGA-A1-A0SQ-01\tTCGA-A2-A04N-01\tTCGA-A2-A04P-01\tTCGA-A2-A04Q-01\tTCGA-A2-A04R-01\tTCGA-A2-A04T-01\tTCGA-A2-A04U-01\tTCGA-A2-A04V-01\tTCGA-A2-A04W-01\tTCGA-A2-A04X-01\tTCGA-A2-A04Y-01\tTCGA-A2-A0CK-01\tTCGA-A2-A0CL-01\tTCGA-A2-A0CM-01\tTCGA-A2-A0CO-01\tTCGA-A2-A0CP-01\tTCGA-A2-A0CQ-01\tTCGA-A2-A0CR-01\tTCGA-A2-A0CS-01\tTCGA-A2-A0CT-01\tTCGA-A2-A0CU-01\tTCGA-A2-A0CV-01\tTCGA-A2-A0CW-01\tTCGA-A2-A0CX-01\tTCGA-A2-A0CZ-01\tTCGA-A2-A0D0-01\tTCGA-A2-A0D1-01\tTCGA-A2-A0D2-01\tTCGA-A2-A0D3-01\tTCGA-A2-A0D4-01\tTCGA-A2-A0EM-01\tTCGA-A2-A0EN-01\tTCGA-A2-A0EO-01\tTCGA-A2-A0EP-01\tTCGA-A2-A0EQ-01\tTCGA-A2-A0ER-01\tTCGA-A2-A0ES-01\tTCGA-A2-A0ET-01\tTCGA-A2-A0EU-01\tTCGA-A2-A0EV-01\tTCGA-A2-A0EW-01\tTCGA-A2-A0EX-01\tTCGA-A2-A0EY-01\tTCGA-A2-A0ST-01\tTCGA-A2-A0SU-01\tTCGA-A2-A0SV-01\tTCGA-A2-A0SW-01\tTCGA-A2-A0SX-01\tTCGA-A2-A0SY-01\tTCGA-A2-A0T0-01\tTCGA-A2-A0T1-01\tTCGA-A2-A0T2-01\tTCGA-A2-A0T3-01\tTCGA-A2-A0T4-01\tTCGA-A2-A0T5-01\tTCGA-A2-A0T6-01\tTCGA-A2-A0T7-01\tTCGA-A2-A0YC-01\tTCGA-A2-A0YD-01\tTCGA-A2-A0YE-01\tTCGA-A2-A0YF-01\tTCGA-A2-A0YG-01\tTCGA-A2-A0YH-01\tTCGA-A2-A0YI-01\tTCGA-A2-A0YJ-01\tTCGA-A2-A0YK-01\tTCGA-A2-A0YL-01\tTCGA-A2-A0YM-01\tTCGA-A2-A0YT-01\tTCGA-A2-A1FV-01\tTCGA-A2-A1FW-01\tTCGA-A2-A1FX-01\tTCGA-A2-A1FZ-01\tTCGA-A2-A1G0-01\tTCGA-A2-A1G1-01\tTCGA-A2-A1G4-01\tTCGA-A2-A1G6-01\tTCGA-A2-A259-01\tTCGA-A2-A25A-01\tTCGA-A2-A25B-01\tTCGA-A2-A25C-01\tTCGA-A2-A25D-01\tTCGA-A2-A25E-01\tTCGA-A2-A25F-01\tTCGA-A2-A3KC-01\tTCGA-A2-A3KD-01\tTCGA-A2-A3XS-01\tTCGA-A2-A3XT-01\tTCGA-A2-A3XU-01\tTCGA-A2-A3XV-01\tTCGA-A2-A3XW-01\tTCGA-A2-A3XX-01\tTCGA-A2-A3XY-01\tTCGA-A2-A3XZ-01\tTCGA-A2-A3Y0-01\tTCGA-A2-A4RW-01\tTCGA-A2-A4RX-01\tTCGA-A2-A4RY-01\tTCGA-A2-A4S0-01\tTCGA-A2-A4S1-01\tTCGA-A2-A4S2-01\tTCGA-A2-A4S3-01\tTCGA-A7-A0CD-01\tTCGA-A7-A0CE-01\tTCGA-A7-A0CG-01\n... [truncated \u2014 43,142 chars total]" + ] + }, + "execution_count": 20, + "metadata": {}, + "output_type": "execute_result" + } + ], "source": [ - "# Align to keeper's sample order before loading\n", + "print(f'RNA-seq df shape: {rna_proc.shape}')\n", + "rna_proc.head(3)" + ] + }, + { + "cell_type": "markdown", + "id": "e67b5945-a98e-49ce-9bb4-f831b82465e6", + "metadata": {}, + "source": [ + "## *Extract DNA methylation data*" + ] + }, + { + "cell_type": "code", + "execution_count": 21, + "id": "c1342471-e475-4967-a7f8-3498996ad680", + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Methylation df shape: (1643, 264)\n" + ] + }, + { + "data": { + "text/plain": [ + "\tTCGA-A1-A0SD-01\tTCGA-A2-A04N-01\tTCGA-A2-A04P-01\tTCGA-A2-A04Q-01\tTCGA-A2-A04T-01\tTCGA-A2-A04U-01\tTCGA-A2-A04V-01\tTCGA-A2-A04W-01\tTCGA-A2-A04X-01\tTCGA-A2-A04Y-01\tTCGA-A2-A0CL-01\tTCGA-A2-A0CM-01\tTCGA-A2-A0CP-01\tTCGA-A2-A0CQ-01\tTCGA-A2-A0CS-01\tTCGA-A2-A0CU-01\tTCGA-A2-A0CV-01\tTCGA-A2-A0CW-01\tTCGA-A2-A0D0-01\tTCGA-A2-A0D2-01\tTCGA-A2-A0D3-01\tTCGA-A2-A0D4-01\tTCGA-A2-A0EM-01\tTCGA-A2-A0EO-01\tTCGA-A2-A0EQ-01\tTCGA-A2-A0ER-01\tTCGA-A2-A0ES-01\tTCGA-A2-A0ET-01\tTCGA-A2-A0EV-01\tTCGA-A2-A0EW-01\tTCGA-A2-A0EX-01\tTCGA-A2-A0EY-01\tTCGA-A2-A0T3-01\tTCGA-A7-A0CD-01\tTCGA-A7-A0CE-01\tTCGA-A7-A0CG-01\tTCGA-A7-A0CH-01\tTCGA-A7-A0CJ-01\tTCGA-A7-A0DA-01\tTCGA-A7-A0DB-01\tTCGA-A8-A06N-01\tTCGA-A8-A06O-01\tTCGA-A8-A06P-01\tTCGA-A8-A06Q-01\tTCGA-A8-A06R-01\tTCGA-A8-A06T-01\tTCGA-A8-A06U-01\tTCGA-A8-A06Y-01\tTCGA-A8-A06Z-01\tTCGA-A8-A076-01\tTCGA-A8-A079-01\tTCGA-A8-A07B-01\tTCGA-A8-A07E-01\tTCGA-A8-A07F-01\tTCGA-A8-A07G-01\tTCGA-A8-A07I-01\tTCGA-A8-A07J-01\tTCGA-A8-A07L-01\tTCGA-A8-A07O-01\tTCGA-A8-A07P-01\tTCGA-A8-A07R-01\tTCGA-A8-A07U-01\tTCGA-A8-A07W-01\tTCGA-A8-A081-01\tTCGA-A8-A082-01\tTCGA-A8-A083-01\tTCGA-A8-A084-01\tTCGA-A8-A085-01\tTCGA-A8-A086-01\tTCGA-A8-A08B-01\tTCGA-A8-A08G-01\tTCGA-A8-A08H-01\tTCGA-A8-A08I-01\tTCGA-A8-A08J-01\tTCGA-A8-A08L-01\tTCGA-A8-A08R-01\tTCGA-A8-A08S-01\tTCGA-A8-A08T-01\tTCGA-A8-A08X-01\tTCGA-A8-A08Z-01\tTCGA-A8-A090-01\tTCGA-A8-A091-01\tTCGA-A8-A092-01\tTCGA-A8-A093-01\tTCGA-A8-A094-01\tTCGA-A8-A095-01\tTCGA-A8-A096-01\tTCGA-A8-A097-01\tTCGA-A8-A099-01\tTCGA-A8-A09A-01\tTCGA-A8-A09B-01\tTCGA-A8-A09C-01\tTCGA-A8-A09D-01\tTCGA-A8-A09G-01\tTCGA-A8-A09I-01\tTCGA-A8-A09K-01\tTCGA-A8-A09M-01\tTCGA-A8-A09N-01\tTCGA-A8-A09Q-01\tTCGA-A8-A09R-01\tTCGA-A8-A09T-01\tTCGA-A8-A09X-01\tTCGA-A8-A09Z-01\tTCGA-A8-A0A1-01\tTCGA-A8-A0A2-01\tTCGA-A8-A0A4-01\tTCGA-A8-A0A7-01\tTCGA-A8-A0A9-01\tTCGA-AN-A04D-01\tTCGA-AN-A0AJ-01\tTCGA-AN-A0AK-01\tTCGA-AN-A0AL-01\tTCGA-AN-A0AM-01\tTCGA-AN-A0AS-01\tTCGA-AN-A0AT-01\tTCGA-AN-A0FF-01\tTCGA-AN-A0FJ-01\tTCGA-AN-A0FK-01\tTCGA-AN-A0FL-01\tTCGA-AN-A0FN-01\tTCGA-AN-A0FS-01\tTCGA-AN-A0FT-01\tTCGA-AN-A0FV-01\tTCGA-AN-A0FW-01\tTCGA-AN-A0FX-01\n... [truncated \u2014 10,103 chars total]" + ] + }, + "execution_count": 22, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "STUDY_NAME_3 = 'Breast Invasive Carcinoma (TCGA, Cell 2015)'\n", + "__, methyl_proc = cbio.get_methylation_data(STUDY_NAME_3)\n", + "\n", + "print(f'Methylation df shape: {methyl_proc.shape}')\n", + "methyl_proc.head(3)" + ] + }, + { + "cell_type": "markdown", + "id": "7c1c7182-0a9b-487d-aa3d-e99f26bb8519", + "metadata": {}, + "source": [ + "## *Extract structural variant data*" + ] + }, + { + "cell_type": "code", + "execution_count": 23, + "id": "3af06a4a-c45e-4866-b6ce-8fbc5073d4a4", + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Structural variant df shape: (3283, 852)\n" + ] + } + ], + "source": [ + "target_col = 'tumorPairedEndReadCount'\n", + "__ , sv_proc = cbio.get_structural_variant_data(STUDY_NAME_2, target_col)\n", + "print(f'Structural variant df shape: {sv_proc.shape}')" + ] + }, + { + "cell_type": "code", + "execution_count": 24, + "id": "bb9b1cf1-0bc9-493c-af2a-ca525d82a88e", + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Top 5 genes by samples with structural variants:\n", + "site1HugoSymbol\n", + "BCAS3 21\n", + "FBXL20 21\n", + "SHANK2 18\n", + "RARA 18\n", + "Unknown 17\n", + "\n", + "84 samples with at least one event\n", + "tumorPairedEndReadCount for top 5 genes:\n", + "sampleId TCGA-3C-AALI-01 TCGA-A1-A0SJ-01 TCGA-A1-A0SM-01 TCGA-A2-A04V-01 TCGA-A2-A04W-01 TCGA-A2-A0CZ-01 TCGA-A2-A0D1-01 TCGA-A2-A0D4-01\n", + "site1HugoSymbol \n", + "BCAS3 16.0 8.0 \n", + "FBXL20 1000.0 31.0 \n", + "SHANK2 11.0 \n", + "RARA 42.0 24.0 10.0 98.0\n", + "Unknown \n" + ] + } + ], + "source": [ + "# Display genes with the most structural variant events\n", + "var_gene = (sv_proc != 'missing').sum(axis=1).sort_values(ascending=False)\n", + "top_genes = var_gene.head(5)\n", + "print(\"Top 5 genes by samples with structural variants:\")\n", + "print(top_genes.to_string())\n", + "\n", + "# Display target_col for the top genes \n", + "top_gene_names = top_genes.index\n", + "subset = sv_proc.loc[top_gene_names].replace('missing', pd.NA).dropna(axis=1, how='all')\n", + "print(f\"{subset.shape[1]} samples with at least one event\")\n", + "print(f\"{target_col} for top 5 genes:\")\n", + "print(subset.iloc[:, :8].to_string())" + ] + }, + { + "cell_type": "markdown", + "id": "74770acb-eb82-4bdf-bc01-0a14ef66b418", + "metadata": {}, + "source": [ + "## Loading CBioPortal data into the Keeper" + ] + }, + { + "cell_type": "markdown", + "id": "b3f7d23e-a6f8-43f9-b14b-f5328649720f", + "metadata": {}, + "source": [ + "### Add data from a study into the Keeper" + ] + }, + { + "cell_type": "code", + "execution_count": 25, + "id": "3cb5337c-3b5b-4c5d-bd09-3f86951dab12", + "metadata": {}, + "outputs": [], + "source": [ + "# Initialise Keeper with sample IDs \n", + "keeper = Keeper(observation_labels=sample_ids)\n", + "\n", + "# Load mutation data \n", + "mut_aligned = mut_proc.reindex(columns=sample_ids)\n", + "keeper.add_data(mut_aligned, 'mutation')\n", + "# Load CNA data\n", "cna_aligned = cna_proc.reindex(columns=sample_ids)\n", - "keeper.add_data(cna_aligned, 'cna')\n", - "print(f'Keeper data keys: {list(keeper.data.keys())}')" + "keeper.add_data(cna_aligned, 'cna')" ] }, { "cell_type": "markdown", - "id": "e0c0020c", + "id": "9a04ab91-7fe8-4210-9821-f11676c6ac33", "metadata": {}, "source": [ - "## Keeper summary" + "### Data summary" ] }, { "cell_type": "code", - "execution_count": null, + "execution_count": 26, "id": "3600e824", "metadata": {}, - "outputs": [], + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Keeper observation count: 10945\n", + "Data modalities loaded: ['mutation', 'cna']\n", + "\n", + " [mutation] shape: (413, 10945) dtype: float64\n", + " [cna] shape: (410, 10945) dtype: int64\n" + ] + } + ], "source": [ "print(f'Keeper observation count: {len(keeper.observation_labels)}')\n", "print(f'Data modalities loaded: {list(keeper.data.keys())}\\n')\n", @@ -939,9 +1313,9 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.10.20" + "version": "3.10.19" } }, "nbformat": 4, "nbformat_minor": 5 -} +} \ No newline at end of file diff --git a/netflow/keepers/keeper.py b/netflow/keepers/keeper.py index a86d8ce..648a36a 100644 --- a/netflow/keepers/keeper.py +++ b/netflow/keepers/keeper.py @@ -25,7 +25,7 @@ # POSER, compute_rw_transitions, # import netflow.InfoNet as InfoNet from .._logging import _gen_logger, set_verbose -from ..prep.preprocessing import PCA_tx, log1p_tx, rand_offset_tx +from ..prep.preprocessing import PCA_tx, log1p_tx, rand_offset_tx, feature_correlations from ..probe import summary # from importlib import reload @@ -464,7 +464,11 @@ def standardize(self, **kwargs): data_z = StandardScaler(**kwargs).fit_transform(data.T).T data_z = pd.DataFrame(data=data_z, columns=data.columns, index=data.index) return data_z - + + def feature_correlations(self): + data = self.to_frame() + corr = feature_correlations(data) + return corr class DistanceKeeper: """ A class to store and handle multiple distances. diff --git a/netflow/prep/extract.py b/netflow/prep/extract.py index 61105be..51eace1 100644 --- a/netflow/prep/extract.py +++ b/netflow/prep/extract.py @@ -1,4 +1,5 @@ import json +import math import pandas as pd import requests from bravado.client import SwaggerClient @@ -22,9 +23,8 @@ class CBioPortalClient: 'sampleListId': {'type': 'string'}, 'sampleIds': {'type': 'array', 'items': {'type': 'string'}}, 'entrezGeneIds': {'type': 'array', 'items': {'type': 'integer'}}, - } + }} } - } }, '/api/molecular-profiles/{molecularProfileId}/molecular-data/fetch': { 'post': { @@ -37,8 +37,7 @@ class CBioPortalClient: 'sampleListId': {'type': 'string'}, 'sampleIds': {'type': 'array', 'items': {'type': 'string'}}, 'entrezGeneIds': {'type': 'array', 'items': {'type': 'integer'}}, - } - } + }} } }, '/api/molecular-profiles/{molecularProfileId}/discrete-copy-number/fetch': { @@ -52,10 +51,22 @@ class CBioPortalClient: 'sampleListId': {'type': 'string'}, 'sampleIds': {'type': 'array', 'items': {'type': 'string'}}, 'entrezGeneIds': {'type': 'array', 'items': {'type': 'integer'}}, - } - } + }} } }, + '/api/structural-variant/fetch': { + 'post': { + 'name': 'structuralVariantFilter', + 'in': 'body', + 'required': True, + 'schema': { + 'type': 'object', + 'properties': { + 'molecularProfileIds': {'type': 'array', 'items': {'type': 'string'}}, + 'entrezGeneIds': {'type': 'array', 'items': {'type': 'integer'}}, + }} + }, + } } SAMPLE_ATTR_MAP = { @@ -131,10 +142,85 @@ def _read_response(self, response_wrapper): """ raw_response = response_wrapper.future.result() result_json = json.loads(raw_response.text) - result_df = pd.json_normalize(result_json) + if not result_json: + result_df = pd.DataFrame() + else: + result_df = pd.json_normalize(result_json) return result_df - - + + + def _read_response_paged(self, bravado_callable, gene_ids, filterName, gene_batch_size=1000, **kwargs): + """ Fetch large datasets in batches of `gene_batch_size` and combining responses into a single DataFrame. + + Parameters + ---------- + bravado_callable : callable + A valid Bravado function (e.g.,`self.cbioportal.Structural_Variants.fetchStructuralVariantsUsingPOST`). + gene_ids : `list` of `int` + List of Entrez gene IDs to query to be split into batches of size `gene_batch_size`. + filterName : `str` + Name of the request parameter expected by `bravado_callable` + gene_batch_size : `int` + Number of Entrez gene IDs to include per request (Default = 1000). + **kwargs + Additional parameters required by `bravado_callable`. + + Returns + ------- + result_df : `pandas.DataFrame` + Concatenated DataFrame of batch responses. + """ + chunks = [] + for i in range(0, len(gene_ids), gene_batch_size): + batch = gene_ids[i: i + gene_batch_size] + batch_kwargs = dict(kwargs) + mf = dict(batch_kwargs.pop(filterName, {})) + mf['entrezGeneIds'] = batch + batch_kwargs[filterName] = mf + resp = bravado_callable(**batch_kwargs) + chunk = self._read_response(resp) + if not chunk.empty: + chunks.append(chunk) + + return pd.concat(chunks, ignore_index=True) if chunks else pd.DataFrame() + + + def _get_col(self, study_info, col_name): + """Return column as list. + Parameters + ---------- + study_info : `pandas.DataFrame` + Information returned by metadata extraction methods. + col_name: `str` + Name of a column of `study_info`. + Returns + ------- + result : `list` + Column values as list or an empty list if the column is absent. + """ + return study_info[col_name].tolist() if col_name in study_info.columns else [] + + + @staticmethod + def _flatten_dotted_columns(df): + """Flatten nested column names returned by pd.json_normalize. + Parameters + ------- + df : `pandas.DataFrame` + Input dataframe. + Returns + ------- + df : `pandas.DataFrame` + Output dataframe with flattened column names. + """ + + for col in df.columns: + if '.' in col: + parts = col.split('.') + new_col = parts[0] + ''.join(p[0].upper() + p[1:] for p in parts[1:] if p) + df.rename(columns={col:new_col}, inplace=True) + return df + def list_studies(self): """ Return names of all studies on cBioPortal. Returns @@ -167,7 +253,7 @@ def get_study_id(self, target_name, verbose=False): resp = self.cbioportal.Studies.getAllStudiesUsingGET() all_studies = self._read_response(resp) - study_list = all_studies['name'].to_list() + study_list = self._get_col(all_studies, 'name') if target_name not in study_list: if verbose: print(f'Available studies: {study_list}') @@ -193,6 +279,7 @@ def get_study_description(self, target_name): target_study_id = self.get_study_id(target_name) resp = self.cbioportal.Studies.getStudyUsingGET(studyId=target_study_id) study_info = self._read_response(resp) + study_info = self._flatten_dotted_columns(study_info) return study_info @@ -211,7 +298,7 @@ def list_study_attributes(self, target_name): target_study_id = self.get_study_id(target_name) resp = self.cbioportal.Clinical_Attributes.getAllClinicalAttributesInStudyUsingGET(studyId=target_study_id) clinical_attributes = self._read_response(resp) - attr_list = clinical_attributes['clinicalAttributeId'].to_list() + attr_list = self._get_col(clinical_attributes, 'clinicalAttributeId') return attr_list @@ -226,14 +313,14 @@ def list_study_data(self, target_name): data_list : `list` of `str` List of molecular data types (names) available in the study (e.g., ["COPY_NUMBER_ALTERATION", "MRNA_EXPRESSION", ...]). id_list : `list` of `str` - List of molecular data IDs is returned. + List of associated molecular profile IDs. """ target_study_id = self.get_study_id(target_name) resp = self.cbioportal.Molecular_Profiles.getAllMolecularProfilesInStudyUsingGET(studyId=target_study_id) mp_info = self._read_response(resp) - data_list = mp_info['molecularAlterationType'].to_list() - id_list = mp_info['molecularProfileId'].to_list() + data_list = self._get_col(mp_info, 'molecularAlterationType') + id_list = self._get_col(mp_info, 'molecularProfileId') return data_list, id_list @@ -261,20 +348,26 @@ def get_study_attribute(self, target_name, attribute): if attribute == 'data_types': resp = self.cbioportal.Molecular_Profiles.getAllMolecularProfilesInStudyUsingGET(studyId=target_study_id) molec_profile_info = self._read_response(resp) - molec_profile_list = molec_profile_info['name'].to_list() - molec_profile_ids = molec_profile_info.loc[molec_profile_info['name'].isin(molec_profile_list),'molecularProfileId'] - dtypes_dict = dict(zip(molec_profile_list, molec_profile_ids)) - molec_profile_df = pd.DataFrame(list(dtypes_dict.items()), columns=['name', 'id']) - ans = molec_profile_df + molec_profile_list = self._get_col(molec_profile_info, 'name') + if len(molec_profile_list)==0 or molec_profile_info.empty: + ans = pd.DataFrame(columns=['name', 'id']) + else: + molec_profile_ids = molec_profile_info.loc[molec_profile_info['name'].isin(molec_profile_list),'molecularProfileId'] + dtypes_dict = dict(zip(molec_profile_list, molec_profile_ids)) + molec_profile_df = pd.DataFrame(list(dtypes_dict.items()), columns=['name', 'id']) + ans = molec_profile_df elif attribute == 'cancer_types': resp = self.cbioportal.Clinical_Data.getAllClinicalDataInStudyUsingGET( studyId=target_study_id, attributeId="CANCER_TYPE") cancer_type_info = self._read_response(resp) - cancer_type_df = cancer_type_info['value'].value_counts().reset_index() - cancer_type_df.columns = ['type', 'sampleCount'] - cancer_type_df = cancer_type_df.sort_values(by='sampleCount', ascending=False) - ans = cancer_type_df.reset_index(drop=True) + if cancer_type_info.empty or 'value' not in cancer_type_info.columns: + ans = pd.DataFrame(columns=['type', 'sampleCount']) + else: + cancer_type_df = cancer_type_info['value'].value_counts().reset_index() + cancer_type_df.columns = ['type', 'sampleCount'] + cancer_type_df = cancer_type_df.sort_values(by='sampleCount', ascending=False) + ans = cancer_type_df.reset_index(drop=True) else: raise ValueError(f'Invalid input Attribute {attribute} requested.') return ans @@ -303,7 +396,7 @@ def get_pt_attribute(self, target_study_id, attr): if attr == 'patient_id': resp = self.cbioportal.Patients.getAllPatientsInStudyUsingGET(studyId=target_study_id) pt_info = self._read_response(resp) - vals = pt_info['patientId'].to_list() + vals = self._get_col(pt_info, 'patientId') pts = vals else: if attr not in valid_attrs: @@ -313,12 +406,15 @@ def get_pt_attribute(self, target_study_id, attr): attributeId=attr.upper(), clinicalDataType='PATIENT') pt_info = self._read_response(resp) - vals = pt_info['value'] - if attr == 'os_months': - vals = vals.astype(float) - vals = vals.to_list() - pts = pt_info['patientId'].to_list() - + if pt_info.empty or 'value' not in pt_info.columns: + vals = [] + pts = [] + else: + vals = pt_info['value'] + if attr == 'os_months': + vals = vals.to_numeric(vals, errors='coerce') + vals = vals.to_list() + pts = self._get_col(pt_info, 'patientId') return vals, pts @@ -372,21 +468,25 @@ def get_sample_attribute(self, target_study_id, attr): if attr == 'sample_id': resp = self.cbioportal.Samples.getAllSamplesInStudyUsingGET(studyId=target_study_id) sample_info = self._read_response(resp) - vals = sample_info['sampleId'].to_list() + vals = self._get_col(sample_info, 'sampleId') samples = vals elif attr == 'patient_id': resp = self.cbioportal.Samples.getAllSamplesInStudyUsingGET(studyId=target_study_id) sample_info = self._read_response(resp) - vals = sample_info['patientId'].to_list() - samples = sample_info['sampleId'].to_list() + vals = self._get_col(sample_info, 'patientId') + samples = self._get_col(sample_info, 'sampleId') else: if attr not in valid_attrs: raise ValueError(f"No attribute '{attr}'. Valid options: {sorted(valid_attrs)}") resp = self.cbioportal.Clinical_Data.getAllClinicalDataInStudyUsingGET( studyId=target_study_id, attributeId=attr.upper()) sample_info = self._read_response(resp) - vals = sample_info['value'].to_list() - samples = sample_info['sampleId'].to_list() + if sample_info.empty or 'value' not in sample_info.columns: + vals = [] + samples = [] + else: + vals = self._get_col(sample_info, 'value') + samples = self._get_col(sample_info, 'sampleId') return vals, samples @@ -409,8 +509,10 @@ def get_sample_info(self, target_name, attr_list): sample_id_list, __ = self.get_sample_attribute(target_study_id, 'sample_id') sample_df = pd.DataFrame(sample_id_list, columns=['sampleId']) + attr_list = list(attr_list) #Create copy if 'patient_id' not in attr_list: attr_list = ['patient_id'] + attr_list + for attr in attr_list: try: vals, samples = self.get_sample_attribute(target_study_id, attr) @@ -442,13 +544,17 @@ def map_samples_to_pts(sample_df, pt_df): return summary_df - def map_entrezid_to_hugosymbol(self, input_df): + def map_entrezid_to_hugosymbol(self, input_df, col_in='entrezGeneId', col_out='hugoGeneSymbol'): """ Map gene IDs to HUGO nomenclature. Parameters ---------- input_df : `pandas.DataFrame` Input dataframe containing `entrezGeneId` column. - + col_in : `str` + Name of column to be mapped to hugo gene symbol. Default: `entrezGeneId`. + Set to 'site1EntrezGeneId' or 'site2EntrezGeneId' to match site1 or site2 IDs for structural variants. + col_out : `str` + Name of output column containing hugo gene symbol.Default: `hugoGeneSymbol`. Returns ------- input_df : `pandas.DataFrame` @@ -459,10 +565,23 @@ def map_entrezid_to_hugosymbol(self, input_df): resp = self.cbioportal.Genes.getAllGenesUsingGET() gene_info = self._read_response(resp) gene_map = gene_info.set_index('entrezGeneId')['hugoGeneSymbol'].to_dict() - input_df['hugoGeneSymbol'] = input_df['entrezGeneId'].map(gene_map).fillna("Unknown") + input_df[col_out] = input_df[col_in].map(gene_map).fillna("Unknown") return input_df + def get_gene_ids(self): + """Return all Entrez gene IDs present in a molecular profile. + Returns + ------- + gene_list : `list` of `int` + List of Entrez gene IDs. + """ + resp = self.cbioportal.Genes.getAllGenesUsingGET() + df = self._read_response(resp) + gene_list = self._get_col(df, 'entrezGeneId') + return gene_list + + def get_clinical_data(self, target_name): """ Extract patient clincal data for a given study. Parameters @@ -479,7 +598,11 @@ def get_clinical_data(self, target_name): pt_df = self.get_pt_info(target_name, ['os_months', 'os_status']) clin_df = CBioPortalClient.map_samples_to_pts(sample_df, pt_df) clin_df = clin_df.set_index('sampleId') - clin_df['osGroup'] = clin_df['osStatus'].str.contains('DECEASED', na=False).astype(float) + + if 'osStatus' in clin_df.columns: + clin_df['osGroup'] = clin_df['osStatus'].str.contains('DECEASED', na=False).astype(float) + else: + clin_df['osGroup'] = float('nan') return clin_df @@ -612,7 +735,7 @@ def get_methylation_data(self, target_name): methyl_raw = self.map_entrezid_to_hugosymbol(methyl_raw) methyl_df = methyl_raw.copy() - methyl_df['value'] = pd.to_numeric(methyl_raw['value'], errors='coerce') + methyl_df['value'] = pd.to_numeric(methyl_df['value'], errors='coerce') # Aggregate methylation by mean & reshape to (n_feats, n_obs) methyl_df = ( @@ -625,81 +748,101 @@ def get_methylation_data(self, target_name): return methyl_raw, methyl_df -def get_rna_seq_data(self, target_name): - """Get RNA-seq expression data for all samples in a study. + def get_rna_seq_data(self, target_name): + """Get RNA-seq expression data for all samples in a study. - Parameters - ---------- - target_name : `str` - Name of a study exactly as displayed on cBioPortal. - - Returns - ------- - rna_raw : `pandas.DataFrame` - DataFrame containing RNA-seq expression values for all samples in the study. - rna_df : `pandas.DataFrame` - rna_raw processed for loading into the keeper. - - """ - target_study_id = self.get_study_id(target_name) + Parameters + ---------- + target_name : `str` + Name of a study exactly as displayed on cBioPortal. - sel_type = 'MRNA_EXPRESSION' - rna_types, rna_ids = self.list_study_data(target_name) - if sel_type not in rna_types: - raise ValueError(f"No profile matching {sel_type} found in study {target_name}.") - else: - match_idx = rna_types.index(sel_type) - mpid = rna_ids[match_idx] #TBD: Check if multiple IDs can have the same study data type + Returns + ------- + rna_raw : `pandas.DataFrame` + DataFrame containing RNA-seq expression values for all samples in the study. + rna_df : `pandas.DataFrame` + rna_raw processed for loading into the keeper. - resp = self.cbioportal.Molecular_Data.fetchAllMolecularDataInMolecularProfileUsingPOST( - molecularProfileId=mpid, - molecularDataFilter={'sampleListId': target_study_id + '_all'}, - projection='SUMMARY') - rna_raw = self._read_response(resp) - rna_raw = self.map_entrezid_to_hugosymbol(rna_raw) + """ + target_study_id = self.get_study_id(target_name) - rna_df = rna_raw.copy() - rna_df['value'] = pd.to_numeric(rna_df['value'], errors='coerce') + sel_type = 'MRNA_EXPRESSION' + rna_types, rna_ids = self.list_study_data(target_name) + if sel_type not in rna_types: + raise ValueError(f"No profile matching {sel_type} found in study {target_name}.") + else: + match_idx = rna_types.index(sel_type) + mpid = rna_ids[match_idx] #TBD: Check if multiple IDs can have the same study data type + + gene_ids = self.get_gene_ids() + rna_raw = self._read_response_paged( + self.cbioportal.Molecular_Data.fetchAllMolecularDataInMolecularProfileUsingPOST, + gene_ids, + gene_batch_size=500, filterName='molecularDataFilter', + molecularProfileId=mpid, + molecularDataFilter={'sampleListId': target_study_id + '_all'}, + projection='SUMMARY') + rna_raw = self.map_entrezid_to_hugosymbol(rna_raw) - # Aggregate by mean, reshape to (n_feats, n_obs) - rna_df = ( - rna_df.groupby(['hugoGeneSymbol', 'sampleId'])['value'] - .mean() - .unstack(level='sampleId') - ) - rna_df.index.name = 'hugoGeneSymbol' - rna_df.columns.name = None + rna_df = rna_raw.copy() + rna_df['value'] = pd.to_numeric(rna_df['value'], errors='coerce') - return rna_raw, rna_df + # Aggregate by mean, reshape to (n_feats, n_obs) + rna_df = ( + rna_df.groupby(['hugoGeneSymbol', 'sampleId'])['value'] + .mean() + .unstack(level='sampleId') + ) + rna_df.index.name = 'hugoGeneSymbol' + rna_df.columns.name = None + return rna_raw, rna_df -def get_structural_variant_data(self, target_name): - """Get structural variant (fusion) data for all samples in a study. - Parameters - ---------- - target_name : `str` - Name of a study exactly as displayed on cBioPortal. + def get_structural_variant_data(self, target_name, target_col, agg_func='first', fill_val='missing'): + """Get structural variant (fusion) data for all samples in a study. - Returns - ------- - sv_df : `pandas.DataFrame` - DataFrame containing structural variant data for all samples in the study. - """ - target_study_id = self.get_study_id(target_name) - - sel_type = 'STRUCTURAL_VARIANT' - sv_types, sv_ids = self.list_study_data(target_name) - if sel_type not in sv_types: - raise ValueError(f"No profile matching {sel_type} found in study {target_name}.") - else: - match_idx = sv_types.index(sel_type) - mpid = sv_ids[match_idx] - - resp = self.cbioportal.Structural_Variants.fetchStructuralVariantsUsingPOST( - molecularProfileId=mpid, - structuralVariantFilter={'sampleListId': target_study_id + '_all'}) - sv_df = self._read_response(resp) - sv_df = self.map_entrezid_to_hugosymbol(sv_df) - return sv_df + Parameters + ---------- + target_name : `str` + Name of a study exactly as displayed on cBioPortal. + target_col : `str` + Column from the raw structural variant table to return + (e.g., 'site2HugoSymbol', 'variantClass', 'eventInfo'). + agg_func : `str` or callable + Aggregation function passed to `pivot_table` to handle duplicates (Default = 'first'). + fill_val : `str`, `int`, or `float` + Value used to fill missing entries in the pivoted matrix (Default = 'missing'). + Returns + ------- + sv_raw : `pandas.DataFrame` + DataFrame containing structural variant data for all smaples. + sv_df : `pandas.DataFrame` of shape (n_site1_genes, n_samples) + sv_raw processed for loading into the keeper. + """ + sel_type = 'STRUCTURAL_VARIANT' + sv_types, sv_ids = self.list_study_data(target_name) + if sel_type not in sv_types: + raise ValueError(f"No profile matching {sel_type} found in study {target_name}.") + else: + match_idx = sv_types.index(sel_type) + mpid = sv_ids[match_idx] + + gene_ids = self.get_gene_ids() + sv_raw = self._read_response_paged(self.cbioportal.Structural_Variants.fetchStructuralVariantsUsingPOST, + gene_ids, gene_batch_size=500, filterName='structuralVariantFilter', + structuralVariantFilter={'molecularProfileIds': [mpid]}) + if sv_raw.empty: + raise ValueError(f"No structural variant data returned for study '{target_name}'.") + + if target_col not in sv_raw.columns: + raise ValueError(f"Column '{target_col}' not found. Available columns: {sv_raw.columns.tolist()}") + + sv_df = sv_raw.copy() + sv_df = self.map_entrezid_to_hugosymbol(sv_df, col_in='site1EntrezGeneId', col_out='site1HugoSymbol') + sv_df = self.map_entrezid_to_hugosymbol(sv_df, col_in='site2EntrezGeneId', col_out='site2HugoSymbol') + sv_df = sv_df.pivot_table(index='site1HugoSymbol', columns='sampleId', + values=target_col, aggfunc=agg_func).fillna(fill_val) + + return sv_raw, sv_df \ No newline at end of file diff --git a/netflow/prep/preprocessing.py b/netflow/prep/preprocessing.py index 7574863..32f4761 100644 --- a/netflow/prep/preprocessing.py +++ b/netflow/prep/preprocessing.py @@ -1,6 +1,7 @@ import numpy as np import pandas as pd import sklearn.decomposition +from scipy.stats import spearmanr def MAF_to_binary_matrix(mut, vc_col='Variant_Classification', @@ -140,3 +141,22 @@ def rand_offset_tx(data, scale, center=0., rand_seed=None): data_unique = data + noise return data_unique + +def feature_correlations(data_arr): + """ Compute correlations between all feature pairs. + + Parameters: + ---------- + data_arr: `np.array` (n_features, n_observations) + Input feature array. + + Returns: + ---------- + corr_arr: `np.ndarray` (n_features, n_features) + Array of correlations between every pair of features. + """ + + result = spearmanr(data_arr, axis=1) + corr_arr = result.correlation + + return corr_arr \ No newline at end of file diff --git a/netflow/probe/summary.py b/netflow/probe/summary.py index f121c45..a92ec13 100644 --- a/netflow/probe/summary.py +++ b/netflow/probe/summary.py @@ -597,75 +597,3 @@ def feature_graph_order_correlation_local(poser, graph_nw, data_df, obs_labels, return corr_dict -def ordered_features_correlation_global(poser, graph_nw, data_df, obs_labels, weights=None): - """ Compute correlation between feature pairs, sorted by global node order. - - Parameters: - ---------- - poser: `netflow.pose.POSER` - The object used to construct the POSE. - graph_nw: `networkx.Graph` - The POSE graph. - data_df: `pandas.DataFrame` (n_features, n_observations) - Feature matrix. - obs_labels: `list` of `str` - List of observation labels corresponding to node IDs. - weights: {`None`, `pandas.DataFrame`, (n, n)} - Dataframe of edge weights between nodes (observations). If `None` unweighted hop count is used. - - Returns: - ---------- - corr_arr: `np.ndarray` (n_features, n_features) - Array of correlations between every pair of features sorted by the global ordering of nodes. - Node order is based on weighted distance if provided. Otherwise, based on hop distance if `weights` is `None` - """ - - node_ord_dict = get_global_node_order(poser, graph_nw, weights=weights) - node_ids = list(node_ord_dict.keys()) - obs_ord_labels = [obs_labels[node_id] for node_id in node_ids] - feat_arr = np.array(data_df.loc[:, obs_ord_labels]) - result = spearmanr(feat_arr, axis=1) - corr_arr = result.correlation - - return corr_arr - - -def ordered_features_correlation_branch(poser, graph_nw, data_df, obs_labels, weights=None, min_branch_size=3): - """ Compute correlations between pairs of features on each branch. - - Parameters: - ---------- - poser: `netflow.pose.POSER` - The object used to construct the POSE. - graph_nw: `networkx.Graph` - The POSE graph. - data_df: `pandas.DataFrame` (n_features, n_observations) - Feature matrix. - obs_labels: `list` of `str` - List of observation labels corresponding to node IDs. - weights: {`None`, `pandas.DataFrame`, (n, n)} - Dataframe of edge weights between nodes (observations). If `None` unweighted hop count is used. - min_branch_size: {`None`, `int`} - Skip branches with <= ``min_branch_size`` observations. - - Returns: - ---------- - corr_dict: `dict` - Dictionary of correlations between feature pairs sorted by node order on each branch. - Node order is based on weighted distance if provided. Otherwise, based on hop distance if `weights` is `None` - """ - - branch_ord_dict = get_branch_node_order(poser, graph_nw, weights=weights, min_branch_size=min_branch_size) - corr_dict = {} - for branch_id, ord_dict in branch_ord_dict.items(): - - branch_ord_node_ids = list(ord_dict.keys()) - observation_subset = [obs_labels[node_id] for node_id in branch_ord_node_ids] - feat_arr = np.array(data_df.loc[:, observation_subset]) - result = spearmanr(feat_arr, axis=1) - corr_dict[branch_id] = result.correlation - - return corr_dict - - -