{ "cells": [ { "cell_type": "markdown", "id": "4e690b3c-4dec-450e-a7f8-f63987e60cdb", "metadata": {}, "source": [ "[![Open In Colab](https://colab.research.google.com/assets/colab-badge.svg)](https://colab.research.google.com/github/lucascamillomd/pyaging/blob/main/tutorials/tutorial_dnam.ipynb) [![Open In nbviewer](https://img.shields.io/badge/View%20in-nbviewer-orange)](https://nbviewer.jupyter.org/github/lucascamillomd/pyaging/blob/main/tutorials/tutorial_dnam.ipynb)" ] }, { "cell_type": "markdown", "id": "62e1ac68-927d-4ca8-a2ab-bd99a7ee52ab", "metadata": {}, "source": [ "# Illumina Human Methylation Arrays" ] }, { "cell_type": "markdown", "id": "41699e8b-e682-4617-867c-fd0a624ae0ef", "metadata": {}, "source": [ "This tutorial is a brief guide for the implementation of an array of bulk DNA-methylation epigenetic clocks that predict age in humans. In this notebook, we will demonstrate the breadth of epigenetic clock models available in `pyaging` by showing:\n", "\n", "- Horvath's 2013 ElasticNet-based clock ([paper](https://genomebiology.biomedcentral.com/articles/10.1186/gb-2013-14-10-r115));\n", " \n", "- AltumAge, a highly accurate deep-learning based clock ([paper](https://www.nature.com/articles/s41514-022-00085-y));\n", " \n", "- PCGrimAge, a principal-component based version of the GrimAge clock ([paper](https://www.nature.com/articles/s43587-022-00248-2));\n", "\n", "- GrimAge2, the latest version of GrimAge ([paper](https://www.aging-us.com/article/204434/text]));\n", "\n", "- DunedinPACE, a biomarker of the pace of aging ([paper](https://elifesciences.org/articles/73420))." ] }, { "cell_type": "markdown", "id": "d8ddf439-0b19-4a1a-9491-3a0e1ee94447", "metadata": {}, "source": [ "We just need two packages for this tutorial." ] }, { "cell_type": "code", "execution_count": 1, "id": "dc347ae0-41ae-46ac-ba50-08cacd4c9241", "metadata": {}, "outputs": [], "source": [ "import pandas as pd\n", "import pyaging as pya" ] }, { "cell_type": "markdown", "id": "f3eebcc0-d3a2-42a4-b87f-1637d0db2588", "metadata": {}, "source": [ "## Download and load example data" ] }, { "cell_type": "markdown", "id": "096cb9c0-67d1-4a37-958e-13da15b15435", "metadata": {}, "source": [ "Let's download the publicly avaiable dataset GSE139307 with Illumina's 450k array. The CpG coverage of the 450k array should be good enough for most clocks." ] }, { "cell_type": "code", "execution_count": 2, "id": "4b2680c6-0c91-4d6a-bd12-6f9cc06932f7", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "|-----> 🏗️ Starting download_example_data function\n", "|-----------> Data found in pyaging_data/GSE139307.pkl\n", "|-----> 🎉 Done! [0.0007s]\n" ] } ], "source": [ "pya.data.download_example_data('GSE139307')" ] }, { "cell_type": "code", "execution_count": 3, "id": "3192d1bd-f5d6-426c-8c65-29971e46c4b4", "metadata": {}, "outputs": [], "source": [ "df = pd.read_pickle('pyaging_data/GSE139307.pkl')" ] }, { "cell_type": "code", "execution_count": 4, "id": "ddf9e1d5-ddb1-42e2-9dfc-9cf4441b326d", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
datasettissue_typeagegendercg00000029cg00000108cg00000109cg00000165cg00000236cg00000289...ch.X.93511680Fch.X.938089Fch.X.94051109Rch.X.94260649Rch.X.967194Fch.X.97129969Rch.X.97133160Rch.X.97651759Fch.X.97737721Fch.X.98007042R
GSM4137709GSE139307sperm84.0M0.0848110.9206960.8568510.0845670.8386990.247273...0.0617510.0459420.0376310.0564550.2498720.0490220.0856910.0374350.0778200.106234
GSM4137710GSE139307sperm69.0M0.0996260.9190730.8900240.1155410.8525840.198103...0.0750770.0418490.0325730.0897900.2502450.0790950.0797560.0462290.0912560.120241
GSM4137711GSE139307sperm69.0M0.1172280.9202760.8943170.1171270.8392580.213410...0.0686790.0495150.0580970.0799190.2997580.0793050.0898150.0653640.0868640.156005
GSM4137712GSE139307sperm69.0M0.0770960.9102040.9084000.0738850.8616150.163276...0.0700910.0332890.0388360.1082130.2954280.0507310.0999430.0475970.0784800.107480
GSM4137713GSE139307sperm67.0M0.0635240.9116080.8846430.0798770.8646540.176169...0.0823680.0384110.0487870.0886310.3166940.0418730.0793030.0488230.0890100.117903
\n", "

5 rows × 485516 columns

\n", "
" ], "text/plain": [ " dataset tissue_type age gender cg00000029 cg00000108 \\\n", "GSM4137709 GSE139307 sperm 84.0 M 0.084811 0.920696 \n", "GSM4137710 GSE139307 sperm 69.0 M 0.099626 0.919073 \n", "GSM4137711 GSE139307 sperm 69.0 M 0.117228 0.920276 \n", "GSM4137712 GSE139307 sperm 69.0 M 0.077096 0.910204 \n", "GSM4137713 GSE139307 sperm 67.0 M 0.063524 0.911608 \n", "\n", " cg00000109 cg00000165 cg00000236 cg00000289 ... \\\n", "GSM4137709 0.856851 0.084567 0.838699 0.247273 ... \n", "GSM4137710 0.890024 0.115541 0.852584 0.198103 ... \n", "GSM4137711 0.894317 0.117127 0.839258 0.213410 ... \n", "GSM4137712 0.908400 0.073885 0.861615 0.163276 ... \n", "GSM4137713 0.884643 0.079877 0.864654 0.176169 ... \n", "\n", " ch.X.93511680F ch.X.938089F ch.X.94051109R ch.X.94260649R \\\n", "GSM4137709 0.061751 0.045942 0.037631 0.056455 \n", "GSM4137710 0.075077 0.041849 0.032573 0.089790 \n", "GSM4137711 0.068679 0.049515 0.058097 0.079919 \n", "GSM4137712 0.070091 0.033289 0.038836 0.108213 \n", "GSM4137713 0.082368 0.038411 0.048787 0.088631 \n", "\n", " ch.X.967194F ch.X.97129969R ch.X.97133160R ch.X.97651759F \\\n", "GSM4137709 0.249872 0.049022 0.085691 0.037435 \n", "GSM4137710 0.250245 0.079095 0.079756 0.046229 \n", "GSM4137711 0.299758 0.079305 0.089815 0.065364 \n", "GSM4137712 0.295428 0.050731 0.099943 0.047597 \n", "GSM4137713 0.316694 0.041873 0.079303 0.048823 \n", "\n", " ch.X.97737721F ch.X.98007042R \n", "GSM4137709 0.077820 0.106234 \n", "GSM4137710 0.091256 0.120241 \n", "GSM4137711 0.086864 0.156005 \n", "GSM4137712 0.078480 0.107480 \n", "GSM4137713 0.089010 0.117903 \n", "\n", "[5 rows x 485516 columns]" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "df.head()" ] }, { "cell_type": "markdown", "id": "cfbebcff-687e-4199-8b8a-6712577cebc0", "metadata": {}, "source": [ "For PCGrimAge and GrimAge2, both age and sex are features. Therefore, to get the full prediction, let's convert the column `gender` into a column called `female`, with 1 being female and 0 being male." ] }, { "cell_type": "code", "execution_count": 5, "id": "671906c4-ae7e-4808-859a-470e19757477", "metadata": {}, "outputs": [], "source": [ "# needs only numerical data (doesn't work with strings)\n", "df['female'] = (df['gender'] == 'F').astype(int)" ] }, { "cell_type": "markdown", "id": "125e979c-3149-4e00-9a78-820b9c89c0ef", "metadata": {}, "source": [ "Moreover, it is important to note that some probes are duplicated in the EPICv2 array, following the format cg#########_BC11 and cg#########_TC11 for the opposite strands. Given that at this moment most clocks have not been trained with EPICv2 data directly, it is recommended to average these probes. This is particularly the case for DunedinPACE, from which some clock probes were duplicated in the update from EPICv1. To remedy this issue, simply use the following function to aggregate any duplicated probes that may be present." ] }, { "cell_type": "code", "execution_count": 6, "id": "cdb19129-49b3-4b81-b264-73ccbc7943c3", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "|-----> 🏗️ Starting epicv2_probe_aggregation function\n", "|-----> ⚙️ Looking for duplicated probes started\n", "|-----------> in progress: 100.0000%\n", "|-----------> There are no duplicated probes. Returning original data\n", "|-----> 🎉 Done! [8.0469s]\n" ] } ], "source": [ "df = pya.pp.epicv2_probe_aggregation(df)" ] }, { "cell_type": "markdown", "id": "798b8998-dc50-44ac-bc4e-c97c4f1ec183", "metadata": {}, "source": [ "## Convert data to AnnData object" ] }, { "cell_type": "markdown", "id": "f8269bff-fadf-4524-801d-6294655c005b", "metadata": {}, "source": [ "AnnData objects are highly flexible and are thus our preferred method of organizing data for age prediction." ] }, { "cell_type": "code", "execution_count": 7, "id": "037192c8-ef2e-440a-88ae-06bddbab80a8", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "|-----> 🏗️ Starting df_to_adata function\n", "|-----> ⚙️ Create anndata object started\n", "|-----------? Dropping 1 columns with only NAs: ['cg01550828'], etc.\n", "|-----> ⚠️ Create anndata object finished [0.4365s]\n", "|-----> ⚙️ Add metadata to anndata started\n", "|-----------> Adding provided metadata to adata.obs\n", "|-----> ✅ Add metadata to anndata finished [0.0010s]\n", "|-----> ⚙️ Log data statistics started\n", "|-----------> There are 37 observations\n", "|-----------> There are 485513 features\n", "|-----------> Total missing values: 489\n", "|-----------> Percentage of missing values: 0.00%\n", "|-----> ✅ Log data statistics finished [0.0292s]\n", "|-----> ⚙️ Impute missing values started\n", "|-----------> Imputing missing values using knn strategy\n", "|-----> ✅ Impute missing values finished [6.2591s]\n", "|-----> ⚙️ Add imputer strategy to adata.uns started\n", "|-----> ✅ Add imputer strategy to adata.uns finished [0.0003s]\n", "|-----> 🎉 Done! [6.8518s]\n" ] } ], "source": [ "adata = pya.pp.df_to_adata(df, metadata_cols=['gender', 'tissue_type', 'dataset'], imputer_strategy='knn')" ] }, { "cell_type": "markdown", "id": "3f4e479d-0b80-4079-9cd5-79f7d5130ba4", "metadata": {}, "source": [ "Note that the original DataFrame is stored in `X_original` under layers. is This is what the `adata` object looks like:" ] }, { "cell_type": "code", "execution_count": 8, "id": "3259b596-5679-42c1-967b-5297f1612cf4", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "AnnData object with n_obs × n_vars = 37 × 485513\n", " obs: 'gender', 'tissue_type', 'dataset'\n", " var: 'percent_na'\n", " uns: 'imputer_strategy'\n", " layers: 'X_original', 'X_imputed'" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "adata" ] }, { "cell_type": "markdown", "id": "bff7c621-f7b7-4fc0-9c91-6f82a9211e8b", "metadata": {}, "source": [ "## Predict age" ] }, { "cell_type": "markdown", "id": "07f5a0ae-1901-4f49-af64-320974584231", "metadata": {}, "source": [ "We can either predict one clock at once or all at the same time. For convenience, let's simply input all four clocks of interest at once. The function is invariant to the capitalization of the clock name. " ] }, { "cell_type": "code", "execution_count": 9, "id": "5383ef0f-ff8b-4e41-bbb6-7fd7ff6cc6be", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "|-----> 🏗️ Starting predict_age function\n", "|-----> ⚙️ Set PyTorch device started\n", "|-----------> Using device: cpu\n", "|-----> ✅ Set PyTorch device finished [0.0017s]\n", "|-----> 🕒 Processing clock: horvath2013\n", "|-----------> ⚙️ Load clock started\n", "|-----------------> Data found in pyaging_data/horvath2013.pt\n", "|-----------> ✅ Load clock finished [0.0059s]\n", "|-----------> ⚙️ Check features in adata started\n", "|-----------------> All features are present in adata.var_names.\n", "|-----------> ✅ Check features in adata finished [0.0936s]\n", "|-----------> ⚙️ Predict ages with model started\n", "|-----------------> There is no preprocessing necessary\n", "|-----------------> The postprocessing method is anti_log_linear\n", "|-----------------> in progress: 100.0000%\n", "|-----------> ✅ Predict ages with model finished [0.0026s]\n", "|-----------> ⚙️ Add predicted ages and clock metadata to adata started\n", "|-----------> ✅ Add predicted ages and clock metadata to adata finished [0.0017s]\n", "|-----> 🕒 Processing clock: altumage\n", "|-----------> ⚙️ Load clock started\n", "|-----------------> Data found in pyaging_data/altumage.pt\n", "|-----------> ✅ Load clock finished [0.0080s]\n", "|-----------> ⚙️ Check features in adata started\n", "|-----------------> All features are present in adata.var_names.\n", "|-----------> ✅ Check features in adata finished [0.1086s]\n", "|-----------> ⚙️ Predict ages with model started\n", "|-----------------> The preprocessing method is scale\n", "|-----------------> There is no postprocessing necessary\n", "|-----------------> in progress: 100.0000%\n", "|-----------> ✅ Predict ages with model finished [0.0087s]\n", "|-----------> ⚙️ Add predicted ages and clock metadata to adata started\n", "|-----------> ✅ Add predicted ages and clock metadata to adata finished [0.0016s]\n", "|-----> 🕒 Processing clock: pcgrimage\n", "|-----------> ⚙️ Load clock started\n", "|-----------------> Data found in pyaging_data/pcgrimage.pt\n", "|-----------> ✅ Load clock finished [0.2855s]\n", "|-----------> ⚙️ Check features in adata started\n", "|-----------------> All features are present in adata.var_names.\n", "|-----------> ✅ Check features in adata finished [0.1395s]\n", "|-----------> ⚙️ Predict ages with model started\n", "|-----------------> There is no preprocessing necessary\n", "|-----------------> There is no postprocessing necessary\n", "|-----------------> in progress: 100.0000%\n", "|-----------> ✅ Predict ages with model finished [0.1859s]\n", "|-----------> ⚙️ Add predicted ages and clock metadata to adata started\n", "|-----------> ✅ Add predicted ages and clock metadata to adata finished [0.0007s]\n", "|-----> 🕒 Processing clock: grimage2\n", "|-----------> ⚙️ Load clock started\n", "|-----------------> Data found in pyaging_data/grimage2.pt\n", "|-----------> ✅ Load clock finished [0.0032s]\n", "|-----------> ⚙️ Check features in adata started\n", "|-----------------> All features are present in adata.var_names.\n", "|-----------> ✅ Check features in adata finished [0.0921s]\n", "|-----------> ⚙️ Predict ages with model started\n", "|-----------------> There is no preprocessing necessary\n", "|-----------------> The postprocessing method is cox_to_years\n", "|-----------------> in progress: 100.0000%\n", "|-----------> ✅ Predict ages with model finished [0.0030s]\n", "|-----------> ⚙️ Add predicted ages and clock metadata to adata started\n", "|-----------> ✅ Add predicted ages and clock metadata to adata finished [0.0014s]\n", "|-----> 🕒 Processing clock: dunedinpace\n", "|-----------> ⚙️ Load clock started\n", "|-----------------> Data found in pyaging_data/dunedinpace.pt\n", "|-----------> ✅ Load clock finished [0.0084s]\n", "|-----------? ⚠️ Clock 'dunedinpace' is for research purposes only. Please check the clock's documentation or notes for more information.\n", "|-----------> ⚙️ Check features in adata started\n", "|-----------------> All features are present in adata.var_names.\n", "|-----------> ⚠️ Check features in adata finished [0.1133s]\n", "|-----------> ⚙️ Predict ages with model started\n", "|-----------------> The preprocessing method is quantile_normalization_with_gold_standard\n", "|-----------------> There is no postprocessing necessary\n", "|-----------------> in progress: 100.0000%\n", "|-----------> ✅ Predict ages with model finished [0.0741s]\n", "|-----------> ⚙️ Add predicted ages and clock metadata to adata started\n", "|-----------> ✅ Add predicted ages and clock metadata to adata finished [0.0011s]\n", "|-----> 🎉 Done! [4.5616s]\n" ] } ], "source": [ "pya.pred.predict_age(adata, ['Horvath2013', 'AltumAge', 'PCGrimAge', 'GrimAge2', 'DunedinPACE'])" ] }, { "cell_type": "code", "execution_count": 10, "id": "9c0c6926-2944-4274-aefa-eb099b8e6737", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
gendertissue_typedatasethorvath2013altumagepcgrimagegrimage2dunedinpace
GSM4137709MspermGSE13930733.62477637.00721395.50611477.5810571.326327
GSM4137710MspermGSE13930728.82934429.42689983.93424465.9263461.215611
GSM4137711MspermGSE13930728.31654522.79892882.70933463.3583411.271091
GSM4137712MspermGSE13930724.85063018.07917384.26946260.2188801.276866
GSM4137713MspermGSE13930725.94211120.07198584.35698561.2359191.262023
\n", "
" ], "text/plain": [ " gender tissue_type dataset horvath2013 altumage pcgrimage \\\n", "GSM4137709 M sperm GSE139307 33.624776 37.007213 95.506114 \n", "GSM4137710 M sperm GSE139307 28.829344 29.426899 83.934244 \n", "GSM4137711 M sperm GSE139307 28.316545 22.798928 82.709334 \n", "GSM4137712 M sperm GSE139307 24.850630 18.079173 84.269462 \n", "GSM4137713 M sperm GSE139307 25.942111 20.071985 84.356985 \n", "\n", " grimage2 dunedinpace \n", "GSM4137709 77.581057 1.326327 \n", "GSM4137710 65.926346 1.215611 \n", "GSM4137711 63.358341 1.271091 \n", "GSM4137712 60.218880 1.276866 \n", "GSM4137713 61.235919 1.262023 " ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "adata.obs.head()" ] }, { "cell_type": "markdown", "id": "0274db5f-c66d-4b74-a4f2-d5e40ed2824e", "metadata": {}, "source": [ "For curiosity, we can also check if there are any correlations amongst these clocks." ] }, { "cell_type": "code", "execution_count": 11, "id": "bbe696da-df04-4e28-b72a-b70bf6a9bffc", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
horvath2013altumagepcgrimagegrimage2dunedinpace
horvath20131.0000000.6762420.2118810.4591930.354771
altumage0.6762421.0000000.1564560.4400440.164101
pcgrimage0.2118810.1564561.0000000.8594900.061491
grimage20.4591930.4400440.8594901.0000000.183725
dunedinpace0.3547710.1641010.0614910.1837251.000000
\n", "
" ], "text/plain": [ " horvath2013 altumage pcgrimage grimage2 dunedinpace\n", "horvath2013 1.000000 0.676242 0.211881 0.459193 0.354771\n", "altumage 0.676242 1.000000 0.156456 0.440044 0.164101\n", "pcgrimage 0.211881 0.156456 1.000000 0.859490 0.061491\n", "grimage2 0.459193 0.440044 0.859490 1.000000 0.183725\n", "dunedinpace 0.354771 0.164101 0.061491 0.183725 1.000000" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" } ], "source": [ "adata.obs.iloc[:, 3:].corr('pearson')" ] }, { "cell_type": "markdown", "id": "069520e6-0ffe-43d5-a7c1-c3f726b0a1ac", "metadata": {}, "source": [ "Having so much information printed can be overwhelming, particularly when running several clocks at once. In such cases, just set verbose to False." ] }, { "cell_type": "code", "execution_count": 12, "id": "7009e4a4-90c6-4c4f-9bbd-bcb3e3836eee", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
gendertissue_typedatasethorvath2013altumagepcgrimagegrimage2dunedinpace
GSM4137709MspermGSE13930733.62477637.00721395.50578077.5810571.326308
GSM4137710MspermGSE13930728.82934429.42689983.93424465.9263461.215614
GSM4137711MspermGSE13930728.31654522.80555182.70933463.3583411.271033
GSM4137712MspermGSE13930724.85063018.06010784.26946260.2188801.276866
GSM4137713MspermGSE13930725.94211120.07198584.35698561.2359191.262023
\n", "
" ], "text/plain": [ " gender tissue_type dataset horvath2013 altumage pcgrimage \\\n", "GSM4137709 M sperm GSE139307 33.624776 37.007213 95.505780 \n", "GSM4137710 M sperm GSE139307 28.829344 29.426899 83.934244 \n", "GSM4137711 M sperm GSE139307 28.316545 22.805551 82.709334 \n", "GSM4137712 M sperm GSE139307 24.850630 18.060107 84.269462 \n", "GSM4137713 M sperm GSE139307 25.942111 20.071985 84.356985 \n", "\n", " grimage2 dunedinpace \n", "GSM4137709 77.581057 1.326308 \n", "GSM4137710 65.926346 1.215614 \n", "GSM4137711 63.358341 1.271033 \n", "GSM4137712 60.218880 1.276866 \n", "GSM4137713 61.235919 1.262023 " ] }, "execution_count": 12, "metadata": {}, "output_type": "execute_result" } ], "source": [ "pya.data.download_example_data('GSE139307', verbose=False)\n", "df = pd.read_pickle('pyaging_data/GSE139307.pkl')\n", "df['female'] = (df['gender'] == 'F').astype(int)\n", "df = pya.pp.epicv2_probe_aggregation(df, verbose=False)\n", "adata = pya.preprocess.df_to_adata(df, metadata_cols=['gender', 'tissue_type', 'dataset'], imputer_strategy='mean', verbose=False)\n", "pya.pred.predict_age(adata, ['Horvath2013', 'AltumAge', 'PCGrimAge', 'GrimAge2', 'DunedinPACE'], verbose=False)\n", "adata.obs.head()" ] }, { "cell_type": "markdown", "id": "6f246126-0be1-47b2-86f4-c81dcc752da0", "metadata": {}, "source": [ "After age prediction, the clocks are added to `adata.obs`. Moreover, the percent of missing values for each clock and other metadata are included in `adata.uns`." ] }, { "cell_type": "code", "execution_count": 13, "id": "14649c41-167e-4771-a4de-1ed2ae72eb51", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "AnnData object with n_obs × n_vars = 37 × 485513\n", " obs: 'gender', 'tissue_type', 'dataset', 'horvath2013', 'altumage', 'pcgrimage', 'grimage2', 'dunedinpace'\n", " var: 'percent_na'\n", " uns: 'imputer_strategy', 'horvath2013_percent_na', 'horvath2013_missing_features', 'horvath2013_metadata', 'altumage_percent_na', 'altumage_missing_features', 'altumage_metadata', 'pcgrimage_percent_na', 'pcgrimage_missing_features', 'pcgrimage_metadata', 'grimage2_percent_na', 'grimage2_missing_features', 'grimage2_metadata', 'dunedinpace_percent_na', 'dunedinpace_missing_features', 'dunedinpace_metadata'\n", " layers: 'X_original', 'X_imputed'" ] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" } ], "source": [ "adata" ] }, { "cell_type": "markdown", "id": "ba818dde-1561-4184-9c54-3aedc5c8de8e", "metadata": {}, "source": [ "We can also look at which features seem to be missing from each clock (if there are any)." ] }, { "cell_type": "code", "execution_count": 14, "id": "083b2f80-c191-4f84-abcc-d18aa9659e99", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[]" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "adata.uns['dunedinpace_missing_features']" ] }, { "cell_type": "markdown", "id": "123ea2ce-8db1-492d-9d13-c57447030ad8", "metadata": {}, "source": [ "## Get citation" ] }, { "cell_type": "markdown", "id": "99ad630f-a3b1-4cf7-a180-b8b56bd548e1", "metadata": {}, "source": [ "The doi, citation, and some metadata are automatically added to the AnnData object under `adata.uns[CLOCKNAME_metadata]`." ] }, { "cell_type": "code", "execution_count": 15, "id": "7d2d10dc-4ffe-4940-a7f1-2041b933f7b6", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "{'clock_name': 'horvath2013',\n", " 'data_type': 'methylation',\n", " 'species': 'Homo sapiens',\n", " 'year': 2013,\n", " 'approved_by_author': '⌛',\n", " 'citation': 'Horvath, Steve. \"DNA methylation age of human tissues and cell types.\" Genome biology 14.10 (2013): 1-20.',\n", " 'doi': 'https://doi.org/10.1186/gb-2013-14-10-r115',\n", " 'notes': None,\n", " 'version': None}" ] }, "execution_count": 15, "metadata": {}, "output_type": "execute_result" } ], "source": [ "adata.uns['horvath2013_metadata']" ] }, { "cell_type": "code", "execution_count": 16, "id": "e2a5311a-ed7d-4e1b-9fbb-b4ad676ce9da", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "{'clock_name': 'altumage',\n", " 'data_type': 'methylation',\n", " 'species': 'Homo sapiens',\n", " 'year': 2022,\n", " 'approved_by_author': '✅',\n", " 'citation': 'de Lima Camillo, Lucas Paulo, Louis R. Lapierre, and Ritambhara Singh. \"A pan-tissue DNA-methylation epigenetic clock based on deep learning.\" npj Aging 8.1 (2022): 4.',\n", " 'doi': 'https://doi.org/10.1038/s41514-022-00085-y',\n", " 'notes': None,\n", " 'version': None}" ] }, "execution_count": 16, "metadata": {}, "output_type": "execute_result" } ], "source": [ "adata.uns['altumage_metadata']" ] }, { "cell_type": "code", "execution_count": 17, "id": "0bbefc0c-acc8-47db-84dc-5ebe80d08500", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "{'clock_name': 'pcgrimage',\n", " 'data_type': 'methylation',\n", " 'species': 'Homo sapiens',\n", " 'year': 2022,\n", " 'approved_by_author': '⌛',\n", " 'citation': 'Higgins-Chen, Albert T., et al. \"A computational solution for bolstering reliability of epigenetic clocks: Implications for clinical trials and longitudinal tracking.\" Nature aging 2.7 (2022): 644-661.',\n", " 'doi': 'https://doi.org/10.1038/s43587-022-00248-2',\n", " 'notes': None,\n", " 'research_only': None,\n", " 'version': None}" ] }, "execution_count": 17, "metadata": {}, "output_type": "execute_result" } ], "source": [ "adata.uns['pcgrimage_metadata']" ] }, { "cell_type": "code", "execution_count": 18, "id": "9ddfb7c3-83ef-44a6-ace3-ffd5553a5770", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "{'clock_name': 'grimage2',\n", " 'data_type': 'methylation',\n", " 'species': 'Homo sapiens',\n", " 'year': 2022,\n", " 'approved_by_author': '⌛',\n", " 'citation': 'Lu, Ake T., et al. \"DNA methylation GrimAge version 2.\" Aging (Albany NY) 14.23 (2022): 9484.',\n", " 'doi': 'https://doi.org/10.18632/aging.204434',\n", " 'notes': None,\n", " 'version': None}" ] }, "execution_count": 18, "metadata": {}, "output_type": "execute_result" } ], "source": [ "adata.uns['grimage2_metadata']" ] }, { "cell_type": "code", "execution_count": 19, "id": "54497c7f-e1cb-4dd4-815a-182af52155b4", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "{'clock_name': 'dunedinpace',\n", " 'data_type': 'methylation',\n", " 'species': 'Homo sapiens',\n", " 'year': 2022,\n", " 'approved_by_author': '✅',\n", " 'citation': 'Belsky, Daniel W., et al. \"DunedinPACE, a DNA methylation biomarker of the pace of aging.\" Elife 11 (2022): e73420.',\n", " 'doi': 'https://doi.org/10.7554/eLife.73420',\n", " 'notes': \"This model is for research purposes only. Commercial users should contact exclusive DunedinPACE licensee TruDiagnosticTM. The automatic failure if fewer than 80% of the CpG probes are available is not implemented and left to the user's discretion.\",\n", " 'version': None,\n", " 'research_only': True}" ] }, "execution_count": 19, "metadata": {}, "output_type": "execute_result" } ], "source": [ "adata.uns['dunedinpace_metadata']" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.9.17" } }, "nbformat": 4, "nbformat_minor": 5 }