283 lines (282 with data), 10.2 kB
{
"cells": [
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"import numpy as np\n",
"import pandas as pd\n",
"import scanpy as sc # v1.6\",\n",
"import sys\n",
"import matplotlib.pyplot as plt\n",
"import os.path\n",
"import anndata\n",
"from matplotlib import rcParams\n",
"import seaborn as sns\n",
"import numba\n",
"import mnnpy\n",
"import scipy"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"#import full COVID-19 PBMC dataset\n",
"os.chdir('/home/ngr18/covid/')\n",
"covid_total = sc.read_h5ad('covid.h5ad')"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"#Import BAL data GSE145926 (reannotated for DC subsets)\n",
"\n",
"os.chdir('/home/ngr18/covid/external_dataset/BAL')\n",
"bal = sc.read_h5ad('full_bal.h5ad')"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"#Subset PBMC data to myeloid populations and reorder categories (for dotplot visualisations)\n",
"\n",
"blood_myeloid = covid_total[covid_total.obs.full_clustering.isin(['CD83_CD14_mono', 'CD14_mono', \n",
" 'CD16_mono', 'C1_CD16_mono',\n",
" 'DC1', 'DC2', 'DC3', 'ASDC', 'pDC', 'DC_prolif',\n",
" 'Mono_prolif']),:]\n",
"\n",
"blood_myeloid.obs.full_clustering = blood_myeloid.obs.full_clustering.cat.reorder_categories([\n",
"'DC1', 'DC2', 'DC3', 'ASDC', 'pDC', 'DC_prolif', \n",
" 'CD83_CD14_mono', 'CD14_mono', 'CD16_mono', 'C1_CD16_mono', 'Mono_prolif'])"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"#Concatenate PBMC and BAL data for visualisations\n",
"\n",
"myeloid = anndata.concat([bal_myeloid, blood_myeloid], index_unique = None)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"##myeloid figure dotplot - fig 2A (left - RNA)\n",
"\n",
"sc.pl.DotPlot(myeloid, [\n",
"'CLEC9A', 'CADM1', 'CLEC10A','CD1C', 'CD14', 'VCAN',\n",
" 'CCR7', 'LAMP3',\n",
" 'AXL', 'SIGLEC6',\n",
" 'LILRA4', 'ITM2C', 'GZMB',\n",
" 'IL1B', 'IER3', 'LDLR', 'CD83',\n",
" 'S100A12', 'CSF3R', \n",
" 'FCGR3A', 'MS4A7', 'LILRB1', 'CSF1R', 'CDKN1C',\n",
" 'C1QA', 'C1QB', 'C1QC', 'CCR1',\n",
" \n",
" 'MARCO', 'MKI67', 'TOP2A'], \n",
" log = True, \n",
" groupby='full_clustering').style(cmap='Blues',dot_edge_color='black', dot_edge_lw=1).swap_axes(False).show(True)\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"##myeloid figure dotplot - fig 2A (right - protein) - using only PBMC data (no CITEseq data for BAL)\n",
"\n",
"sc.pl.DotPlot(blood_myeloid, ['AB_CD141', 'AB_CLEC9A', 'AB_KIT','AB_BTLA',\n",
" 'AB_CD1C', 'AB_CD101', 'AB_FcERIa',\n",
" 'AB_CD5', \n",
" 'AB_CD123', 'AB_CD45RA', 'AB_CD304',\n",
" 'AB_CD14', 'AB_CD99', 'AB_CD64', \n",
" 'AB_CR1', 'AB_ITGAM',\n",
" 'AB_CD16', 'AB_C5AR1',\n",
" 'AB_CX3CR1'], \n",
" log = True, \n",
" groupby='full_clustering', expression_cutoff=0.15).style(cmap='YlOrRd',dot_edge_color='black', dot_edge_lw=1).swap_axes(False).show(True)\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"#subset data again to blood monocyte and BAL macrophages\n",
"\n",
"monos = myeloid[myeloid.obs.full_clustering.isin(['CD83_CD14_mono', 'CD14_mono', 'CD16_mono', 'C1_CD16_mono', \n",
" 'bal_Mac']),:]"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"#Figure 2C left/upper - analysis of healthy only\n",
"healthy = monos[monos.obs.Status_on_day_collection_summary.isin(['bal_healthy', 'Healthy']),:]\n",
"\n",
"sc.pp.normalize_total(healthy, target_sum=1e4)\n",
"\n",
"sc.pp.log1p(healthy)\n",
"\n",
"sc.pp.highly_variable_genes(healthy, n_top_genes = 3000, flavor = 'seurat')\n",
"\n",
"healthy.raw = healthy\n",
"\n",
"healthy = healthy[:, healthy.var.highly_variable]\n",
"\n",
"sc.pp.scale(healthy, max_value=10)\n",
"\n",
"sc.tl.pca(healthy, svd_solver='arpack')\n",
"\n",
"sc.external.pp.harmony_integrate(healthy, 'sample_id', basis='X_pca', adjusted_basis='X_pca_harmony')\n",
"\n",
"sc.pp.neighbors(healthy, n_neighbors=10, n_pcs=50, use_rep = 'X_pca_harmony')\n",
"\n",
"sc.tl.paga(healthy, groups='full_clustering')\n",
"\n",
"#recolor clusters for consistency with bar graph\n",
"\n",
"healthy.uns['full_clustering_colors'][0] = '#76DDC9'\n",
"healthy.uns['full_clustering_colors'][1] = '#E28686'\n",
"healthy.uns['full_clustering_colors'][2] = '#C2A3E2'\n",
"healthy.uns['full_clustering_colors'][3] = '#991111'\n",
"healthy.uns['full_clustering_colors'][4] = '#ffff00'\n",
"\n",
"\n",
"sc.pl.paga(test, color=['full_clustering'], threshold = 0.11, save = 'healthy_myeloid_PAGA.pdf',\n",
" fontsize = 0, node_size_scale = 10, min_edge_width = 5, frameon = False)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"#Figure 2C left/lower - limit to only covid samples\n",
"covid_mono = covid_mono[covid_mono.obs.Status_on_day_collection_summary.isin(['Asymptomatic', 'Critical', \n",
" 'Mild', 'Moderate', 'Severe',\n",
" 'bal_mild', 'bal_severe']),:]\n",
"\n",
"sc.pp.normalize_total(covid_mono, target_sum=1e4)\n",
"\n",
"sc.pp.log1p(covid_mono)\n",
"\n",
"sc.pp.highly_variable_genes(covid_mono, n_top_genes = 3000, flavor = 'seurat')\n",
"\n",
"covid_mono.raw = covid_mono\n",
"\n",
"covid_mono = covid_mono[:, covid_mono.var.highly_variable]\n",
"\n",
"sc.pp.scale(covid_mono, max_value=10)\n",
"\n",
"sc.tl.pca(covid_mono, svd_solver='arpack')\n",
"\n",
"sc.external.pp.harmony_integrate(covid_mono, 'sample_id', basis='X_pca', adjusted_basis='X_pca_harmony')\n",
"\n",
"sc.pp.neighbors(covid_mono, n_neighbors=10, n_pcs=50, use_rep = 'X_pca_harmony')\n",
"\n",
"sc.tl.paga(covid_mono, groups='full_clustering')\n",
"\n",
"#Again, recolor for consistency with bar graph\n",
"\n",
"covid_mono.uns['full_clustering_colors'][0] = '#76DDC9'\n",
"covid_mono.uns['full_clustering_colors'][1] = '#E28686'\n",
"covid_mono.uns['full_clustering_colors'][2] = '#C2A3E2'\n",
"covid_mono.uns['full_clustering_colors'][3] = '#991111'\n",
"covid_mono.uns['full_clustering_colors'][4] = '#ffff00'\n",
"\n",
"\n",
"sc.pl.paga(covid_mono, color=['full_clustering'], threshold = 0.2, save = 'covid_myeloid_PAGA.pdf',\n",
" fontsize = 0, node_size_scale = 10, min_edge_width = 5, frameon = False)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"#Pseudotime plot\n",
"\n",
"covid_mono.uns['iroot'] = np.flatnonzero(covid_mono.obs['full_clustering'] == 'CD83_CD14_mono')[0]\n",
"\n",
"\n",
"sc.pp.log1p(covid_mono)\n",
"sc.pp.scale(covid_mono)\n",
" \n",
" gene_names = ['NFKBIA', 'KLF6', 'VIM', 'CD14', 'S100A8',\n",
" 'HLA-DPA1', 'HLA-DPB1', 'FCGR3A','CTSC', 'CTSL', \n",
" 'CCR1','RSAD2', 'CCL2', 'CXCL10', 'CCL7', 'TNFSF10'] # monocyte\n",
"\n",
"paths = [('erythrocytes', ['mono1', 'mono2', 'mono3', 'mono4', 'bal_Mac']),\n",
" ('neutrophils', ['mono1', 'mono2', 'mono3', 'mono4', 'bal_Mac']),\n",
" ('monocytes', ['mono1', 'mono2', 'mono3', 'mono4', 'bal_Mac'])]\n",
"\n",
"test.obs['distance'] = test.obs['dpt_pseudotime']\n",
"\n",
"_, axs = plt.subplots(ncols=3, figsize=(6, 5), gridspec_kw={'wspace': 0.05, 'left': 0.12})\n",
"plt.subplots_adjust(left=0.05, right=0.98, top=0.82, bottom=0.2)\n",
"for ipath, (descr, path) in enumerate(paths):\n",
" _, data = sc.pl.paga_path(\n",
" test, path, gene_names,\n",
" show_node_names=False,\n",
" ax=axs[ipath],\n",
" ytick_fontsize=12,\n",
" left_margin=0.15,\n",
" n_avg=50,\n",
" annotations=['distance'],\n",
" show_yticks=True if ipath==0 else False,\n",
" show_colorbar=False,\n",
" color_map='coolwarm',\n",
" color_maps_annotations={'distance': 'viridis'},\n",
" title='{} path'.format(descr),\n",
" return_data=True,\n",
" show=False, normalize_to_zero_one = True)"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "covid_py",
"language": "python",
"name": "covid_py"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.8.5"
}
},
"nbformat": 4,
"nbformat_minor": 4
}
Datasets
Models
