Skip to content

Commit 3860a17

Browse files
committed
feat: add pdac_nac project
1 parent d1cda47 commit 3860a17

File tree

1 file changed

+238
-0
lines changed

1 file changed

+238
-0
lines changed

notebooks/00-pdac_nac.ipynb

+238
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,238 @@
1+
{
2+
"cells": [
3+
{
4+
"cell_type": "code",
5+
"execution_count": 1,
6+
"metadata": {},
7+
"outputs": [],
8+
"source": [
9+
"using Muon\n",
10+
"using RData\n",
11+
"using Revise\n",
12+
"using ISCHIA\n",
13+
"using DataFrames\n",
14+
"using Combinatorics"
15+
]
16+
},
17+
{
18+
"cell_type": "code",
19+
"execution_count": 2,
20+
"metadata": {},
21+
"outputs": [
22+
{
23+
"data": {
24+
"text/plain": [
25+
"AnnData object 216180 ✕ 17764"
26+
]
27+
},
28+
"execution_count": 2,
29+
"metadata": {},
30+
"output_type": "execute_result"
31+
}
32+
],
33+
"source": [
34+
"spatial_object = readh5ad(\"../data/06-pdac_nac-clusters-rmv_unk.h5ad\")\n",
35+
"# adata = readh5ad(\"../data/05-pdac_nac-clusters.h5ad\")\n",
36+
"lr_network = load(\"../data/lr_network.rds\")\n",
37+
"\n",
38+
"# # Remove spots where neoadjuvant_chemo is unknown\n",
39+
"# mask = .!(adata.obs.neoadjuvant_chemo .== \"Unknown\")\n",
40+
"# spatial_object = @view adata[mask, :]\n",
41+
"# # spatial_object = adata[mask, :]\n",
42+
"\n",
43+
"spatial_object"
44+
]
45+
},
46+
{
47+
"cell_type": "code",
48+
"execution_count": 3,
49+
"metadata": {},
50+
"outputs": [
51+
{
52+
"data": {
53+
"text/plain": [
54+
"2-element Vector{String}:\n",
55+
" \"No\"\n",
56+
" \"Yes\""
57+
]
58+
},
59+
"metadata": {},
60+
"output_type": "display_data"
61+
}
62+
],
63+
"source": [
64+
"display(unique(spatial_object.obs.neoadjuvant_chemo))"
65+
]
66+
},
67+
{
68+
"cell_type": "code",
69+
"execution_count": 12,
70+
"metadata": {},
71+
"outputs": [],
72+
"source": [
73+
"gene_names = collect(spatial_object.var_names)\n",
74+
"spatial_object.var.name = gene_names\n",
75+
"\n",
76+
"# Create LR_Pairs column\n",
77+
"lr_network[!, :LR_Pairs] = string.(lr_network.from, \"_\", lr_network.to);\n",
78+
"lr_network = lr_network[:, [:from, :to, :LR_Pairs]];\n",
79+
"\n",
80+
"# Filter lr_network to only include genes in adata\n",
81+
"from_filter = in.(lr_network[!, :from], Ref(gene_names))\n",
82+
"to_filter = in.(lr_network[:, :to], Ref(gene_names))\n",
83+
"all_LR_network = lr_network[from_filter .& to_filter, :];"
84+
]
85+
},
86+
{
87+
"cell_type": "code",
88+
"execution_count": 13,
89+
"metadata": {},
90+
"outputs": [],
91+
"source": [
92+
"# Extract unique genes and common genes\n",
93+
"all_LR_genes = unique(vcat(all_LR_network[:, :from], all_LR_network[:, :to]))\n",
94+
"all_LR_genes_comm = intersect(all_LR_genes, collect(gene_names));\n",
95+
"\n",
96+
"# Create LR.pairs and LR.pairs.AllCombos\n",
97+
"LR_pairs = all_LR_network[:, :LR_Pairs]\n",
98+
"all_combos = [join(combo, \"_\") for combo in combinations(all_LR_genes_comm, 2)];"
99+
]
100+
},
101+
{
102+
"cell_type": "code",
103+
"execution_count": 14,
104+
"metadata": {},
105+
"outputs": [
106+
{
107+
"data": {
108+
"text/plain": [
109+
"\"neoadjuvant_chemo\""
110+
]
111+
},
112+
"execution_count": 14,
113+
"metadata": {},
114+
"output_type": "execute_result"
115+
}
116+
],
117+
"source": [
118+
"# spatial_object = adata\n",
119+
"LR_list = all_LR_genes_comm\n",
120+
"LR_pairs = LR_pairs\n",
121+
"exp_th = 1\n",
122+
"corr_th = 0.2;\n",
123+
"\n",
124+
"cc_column = \"CC_k10\"\n",
125+
"cc_list = [\"CC10\"]\n",
126+
"# Condition = unique(spatial_object.obs[!, \"orig.ident\"])\n",
127+
"condition_list = [\"Yes\", \"No\"]\n",
128+
"condition_column = \"neoadjuvant_chemo\""
129+
]
130+
},
131+
{
132+
"cell_type": "code",
133+
"execution_count": 15,
134+
"metadata": {},
135+
"outputs": [
136+
{
137+
"name": "stdout",
138+
"output_type": "stream",
139+
"text": [
140+
"Running for CC10\n",
141+
"Running for Yes\n",
142+
"Preparing L-R presence/absence matrix\n",
143+
"Calculating L-R pairs correlation\n",
144+
"Preparing for cooccurrence\n"
145+
]
146+
},
147+
{
148+
"name": "stderr",
149+
"output_type": "stream",
150+
"text": [
151+
"\r"
152+
]
153+
},
154+
{
155+
"name": "stdout",
156+
"output_type": "stream",
157+
"text": [
158+
"Cooccurrence calculation starts...\n"
159+
]
160+
},
161+
{
162+
"name": "stderr",
163+
"output_type": "stream",
164+
"text": [
165+
"\u001b[32mCalculate Incidence 100%|████████████████████████████████| Time: 0:02:17\u001b[39mm\n",
166+
"\u001b[32mCalculate Co-occurrences 100%|███████████████████████████| Time: 0:02:15\u001b[39m\n",
167+
"\u001b[32mMain Comp 100%|██████████████████████████████████████████| Time: 0:00:43\u001b[39m\n"
168+
]
169+
},
170+
{
171+
"name": "stdout",
172+
"output_type": "stream",
173+
"text": [
174+
"Cooccurrence calculation ended\n",
175+
"\n",
176+
"Summary of cooccurrence results:\n",
177+
"Of 824970 species pair combinations, 573708 pairs (69.54%) were removed from the analysis because expected co-occurrence was < 1 and\n",
178+
"251262 pairs were analyzed\n",
179+
"\n",
180+
"Cooccurrence Summary:\n",
181+
"\n",
182+
"Species => 1285\n",
183+
"Non-random (%) => 66.3\n",
184+
"Sites => 4494\n",
185+
"Negative => 2665\n",
186+
"Random => 84781\n",
187+
"Positive => 163816\n",
188+
"Unclassifiable => 0\n"
189+
]
190+
},
191+
{
192+
"name": "stderr",
193+
"output_type": "stream",
194+
"text": [
195+
"\u001b[32mCalculate Significantly occurring pairs 2%|█ | ETA: 1 days, 7:24:31\u001b[39mm\r"
196+
]
197+
}
198+
],
199+
"source": [
200+
"for cc in cc_list\n",
201+
" println(\"Running for $cc\")\n",
202+
" for condition in condition_list\n",
203+
" println(\"Running for $condition\")\n",
204+
" lr_result = find_enriched_LR_pairs(\n",
205+
" spatial_object,\n",
206+
" [cc],\n",
207+
" [condition],\n",
208+
" LR_list,\n",
209+
" LR_pairs,\n",
210+
" exp_th,\n",
211+
" corr_th,\n",
212+
" cc_column=cc_column,\n",
213+
" condition_column=condition_column\n",
214+
" )\n",
215+
"\n",
216+
" CSV.write(\"outputs/pdac_nac/$(cc)_lr_enrichment_$(condition).csv\", lr_result[\"enriched_LRs\"])\n",
217+
" CSV.write(\"outputs/pdac_nac/$(cc)_cooccurr_mat_$(condition).csv\", lr_result[\"cooccurrence_table\"].results)\n",
218+
" end\n",
219+
"end"
220+
]
221+
}
222+
],
223+
"metadata": {
224+
"kernelspec": {
225+
"display_name": "Julia 1.9.3",
226+
"language": "julia",
227+
"name": "julia-1.9"
228+
},
229+
"language_info": {
230+
"file_extension": ".jl",
231+
"mimetype": "application/julia",
232+
"name": "julia",
233+
"version": "1.10.4"
234+
}
235+
},
236+
"nbformat": 4,
237+
"nbformat_minor": 2
238+
}

0 commit comments

Comments
 (0)