2_SpaGCN_DLPFC.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "dfb6c73e-c6d7-4588-8c12-567ea4eeabf6",
   "metadata": {},
   "outputs": [],
   "source": [
    "import random\n",
    "from pathlib import Path\n",
    "\n",
    "import anndata as ad\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import scanpy as sc\n",
    "import scipy as sp\n",
    "import SpaGCN as spg\n",
    "import torch\n",
    "from PIL import Image\n",
    "\n",
    "# disable DecompressionBombWarning\n",
    "Image.MAX_IMAGE_PIXELS = None"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "859ecd66-6c44-4b32-8ece-ac0da99e822e",
   "metadata": {},
   "outputs": [],
   "source": [
    "data_dir = Path(\"./data/LIBD_DLPFC\")\n",
    "result_dir = Path(\"./results/LIBD_DLPFC\")\n",
    "\n",
    "technology = \"Visium\"\n",
    "seed = 42\n",
    "\n",
    "alpha = 1\n",
    "beta = 49\n",
    "p = 0.5"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "b1338605-5ea2-4c2f-86e8-c386dbe0331d",
   "metadata": {},
   "outputs": [],
   "source": [
    "metadata = pd.read_table(data_dir / \"samples.tsv\").loc[:, [\"directory\", \"n_clusters\"]]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "e14c1ef8-4157-4d16-92f2-d2f1e6d450cd",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Set seed\n",
    "random.seed(seed)\n",
    "torch.manual_seed(seed)\n",
    "np.random.seed(seed)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "31a05e97-f76c-4347-b3f2-cf07c3a28196",
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_anndata(path):\n",
    "    X = sp.io.mmread(path / \"counts.mtx\").tocsr()\n",
    "\n",
    "    observations = pd.read_table(path / \"observations.tsv\", index_col=0)\n",
    "    features = pd.read_table(path / \"features.tsv\", index_col=0)\n",
    "\n",
    "    coordinates = (\n",
    "        pd.read_table(path / \"coordinates.tsv\", index_col=0)\n",
    "        .loc[observations.index, :]\n",
    "        .to_numpy()\n",
    "    )\n",
    "\n",
    "    adata = ad.AnnData(\n",
    "        X=X, obs=observations, var=features, obsm={\"spatial_pixel\": coordinates}\n",
    "    )\n",
    "\n",
    "    adata.uns[\"image\"] = np.array(Image.open(path / \"H_E.tiff\"))\n",
    "\n",
    "    return adata"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "50d5b4f3-09d3-4b0e-9a4f-fbb956f4bcef",
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Processing Br5292_151507\n",
      "Calculateing adj matrix using histology image...\n",
      "Var of c0,c1,c2 =  98.7653707607211 347.04188339945216 67.83401448585366\n",
      "Var of x,y,z =  1123.1580554335508 435.69138627328755 1123.1580554335505\n",
      "Run 1: l [0.01, 1000], p [0.0, 4213.694327626796]\n",
      "Run 2: l [0.01, 500.005], p [0.0, 4180.11669921875]\n",
      "Run 3: l [0.01, 250.0075], p [0.0, 4050.718017578125]\n",
      "Run 4: l [0.01, 125.00874999999999], p [0.0, 3601.731689453125]\n",
      "Run 5: l [0.01, 62.509375], p [0.0, 2493.169921875]\n",
      "Run 6: l [0.01, 31.2596875], p [0.0, 1146.388916015625]\n",
      "Run 7: l [0.01, 15.63484375], p [0.0, 349.5950012207031]\n",
      "Run 8: l [0.01, 7.822421875], p [0.0, 77.59295654296875]\n",
      "Run 9: l [0.01, 3.9162109375], p [0.0, 13.964275360107422]\n",
      "Run 10: l [0.01, 1.9631054687499998], p [0.0, 2.123582363128662]\n",
      "Run 11: l [0.9865527343749999, 1.9631054687499998], p [0.23047912120819092, 2.123582363128662]\n",
      "Run 12: l [0.9865527343749999, 1.4748291015625], p [0.23047912120819092, 0.9038530588150024]\n",
      "recommended l =  1.23069091796875\n",
      "Start at res =  0.4 step =  0.1\n",
      "Initializing cluster centers with louvain, resolution =  0.4\n",
      "Epoch  0\n",
      "Res =  0.4 Num of clusters =  5\n",
      "Initializing cluster centers with louvain, resolution =  0.5\n",
      "Epoch  0\n",
      "Res =  0.5 Num of clusters =  6\n",
      "Res changed to 0.5\n",
      "Initializing cluster centers with louvain, resolution =  0.6\n",
      "Epoch  0\n",
      "Res =  0.6 Num of clusters =  7\n",
      "recommended res =  0.6\n",
      "Initializing cluster centers with louvain, resolution =  0.6\n",
      "Epoch  0\n",
      "Epoch  10\n",
      "Epoch  20\n",
      "Epoch  30\n",
      "Epoch  40\n",
      "Epoch  50\n",
      "Epoch  60\n",
      "Epoch  70\n",
      "Epoch  80\n",
      "Epoch  90\n",
      "Epoch  100\n",
      "delta_label  0.000946521533364884 < tol  0.001\n",
      "Reach tolerance threshold. Stopping training.\n",
      "Total epoch: 100\n",
      "Calculateing adj matrix using xy only...\n",
      "Processing Br5292_151508\n",
      "Calculateing adj matrix using histology image...\n",
      "Var of c0,c1,c2 =  93.7957216228338 334.65744596919296 64.13796797670821\n",
      "Var of x,y,z =  1135.0556850308355 460.842098132226 1135.0556850308353\n",
      "Run 1: l [0.01, 1000], p [0.0, 4371.0581513736515]\n",
      "Run 2: l [0.01, 500.005], p [0.0, 4335.59765625]\n",
      "Run 3: l [0.01, 250.0075], p [0.0, 4199.0146484375]\n",
      "Run 4: l [0.01, 125.00874999999999], p [0.0, 3725.5009765625]\n",
      "Run 5: l [0.01, 62.509375], p [0.0, 2553.836669921875]\n",
      "Run 6: l [0.01, 31.2596875], p [0.0, 1134.8193359375]\n",
      "Run 7: l [0.01, 15.63484375], p [0.0, 334.38372802734375]\n",
      "Run 8: l [0.01, 7.822421875], p [0.0, 72.96975708007812]\n",
      "Run 9: l [0.01, 3.9162109375], p [0.0, 13.101568222045898]\n",
      "Run 10: l [0.01, 1.9631054687499998], p [0.0, 1.9971027374267578]\n",
      "Run 11: l [0.9865527343749999, 1.9631054687499998], p [0.2121880054473877, 1.9971027374267578]\n",
      "Run 12: l [0.9865527343749999, 1.4748291015625], p [0.2121880054473877, 0.8460689783096313]\n",
      "Run 13: l [1.23069091796875, 1.4748291015625], p [0.4699711799621582, 0.8460689783096313]\n",
      "Run 14: l [1.23069091796875, 1.352760009765625], p [0.4699711799621582, 0.6424710750579834]\n",
      "Run 15: l [1.23069091796875, 1.2917254638671873], p [0.4699711799621582, 0.5524231195449829]\n",
      "Run 16: l [1.23069091796875, 1.2612081909179687], p [0.4699711799621582, 0.5102592706680298]\n",
      "Run 17: l [1.2459495544433592, 1.2612081909179687], p [0.489882230758667, 0.5102592706680298]\n",
      "recommended l =  1.2535788726806638\n",
      "Start at res =  0.4 step =  0.1\n",
      "Initializing cluster centers with louvain, resolution =  0.4\n",
      "Epoch  0\n",
      "Res =  0.4 Num of clusters =  4\n",
      "Initializing cluster centers with louvain, resolution =  0.5\n",
      "Epoch  0\n",
      "Res =  0.5 Num of clusters =  5\n",
      "Res changed to 0.5\n",
      "Initializing cluster centers with louvain, resolution =  0.6\n",
      "Epoch  0\n",
      "Res =  0.6 Num of clusters =  6\n",
      "Res changed to 0.6\n",
      "Initializing cluster centers with louvain, resolution =  0.7\n",
      "Epoch  0\n",
      "Res =  0.7 Num of clusters =  6\n",
      "Res changed to 0.7\n",
      "Initializing cluster centers with louvain, resolution =  0.7999999999999999\n",
      "Epoch  0\n",
      "Res =  0.7999999999999999 Num of clusters =  7\n",
      "recommended res =  0.7999999999999999\n",
      "Initializing cluster centers with louvain, resolution =  0.7999999999999999\n",
      "Epoch  0\n",
      "Epoch  10\n",
      "Epoch  20\n",
      "Epoch  30\n",
      "Epoch  40\n",
      "Epoch  50\n",
      "delta_label  0.0009124087591240876 < tol  0.001\n",
      "Reach tolerance threshold. Stopping training.\n",
      "Total epoch: 58\n",
      "Calculateing adj matrix using xy only...\n",
      "Processing Br5292_151509\n",
      "Calculateing adj matrix using histology image...\n",
      "Var of c0,c1,c2 =  93.78375244682724 288.79133117268105 68.38564040125236\n",
      "Var of x,y,z =  1289.1861391829373 476.1081024539383 1289.186139182937\n",
      "Run 1: l [0.01, 1000], p [0.0, 4773.417309686415]\n",
      "Run 2: l [0.01, 500.005], p [0.0, 4730.2060546875]\n",
      "Run 3: l [0.01, 250.0075], p [0.0, 4565.0185546875]\n",
      "Run 4: l [0.01, 125.00874999999999], p [0.0, 4006.449462890625]\n",
      "Run 5: l [0.01, 62.509375], p [0.0, 2689.71923828125]\n",
      "Run 6: l [0.01, 31.2596875], p [0.0, 1129.8826904296875]\n",
      "Run 7: l [0.01, 15.63484375], p [0.0, 307.18756103515625]\n",
      "Run 8: l [0.01, 7.822421875], p [0.0, 63.95890808105469]\n",
      "Run 9: l [0.01, 3.9162109375], p [0.0, 10.807766914367676]\n",
      "Run 10: l [0.01, 1.9631054687499998], p [0.0, 1.5493807792663574]\n",
      "Run 11: l [0.9865527343749999, 1.9631054687499998], p [0.1597144603729248, 1.5493807792663574]\n",
      "Run 12: l [0.9865527343749999, 1.4748291015625], p [0.1597144603729248, 0.6446664333343506]\n",
      "Run 13: l [1.23069091796875, 1.4748291015625], p [0.3555208444595337, 0.6446664333343506]\n",
      "Run 14: l [1.352760009765625, 1.4748291015625], p [0.48764967918395996, 0.6446664333343506]\n",
      "Run 15: l [1.352760009765625, 1.4137945556640625], p [0.48764967918395996, 0.5629539489746094]\n",
      "Run 16: l [1.352760009765625, 1.3832772827148436], p [0.48764967918395996, 0.5245124101638794]\n",
      "recommended l =  1.3680186462402344\n",
      "Start at res =  0.4 step =  0.1\n",
      "Initializing cluster centers with louvain, resolution =  0.4\n",
      "Epoch  0\n",
      "Res =  0.4 Num of clusters =  4\n",
      "Initializing cluster centers with louvain, resolution =  0.5\n",
      "Epoch  0\n",
      "Res =  0.5 Num of clusters =  5\n",
      "Res changed to 0.5\n",
      "Initializing cluster centers with louvain, resolution =  0.6\n",
      "Epoch  0\n",
      "Res =  0.6 Num of clusters =  5\n",
      "Res changed to 0.6\n",
      "Initializing cluster centers with louvain, resolution =  0.7\n",
      "Epoch  0\n",
      "Res =  0.7 Num of clusters =  6\n",
      "Res changed to 0.7\n",
      "Initializing cluster centers with louvain, resolution =  0.7999999999999999\n",
      "Epoch  0\n",
      "Res =  0.7999999999999999 Num of clusters =  7\n",
      "recommended res =  0.7999999999999999\n",
      "Initializing cluster centers with louvain, resolution =  0.7999999999999999\n",
      "Epoch  0\n",
      "Epoch  10\n",
      "Epoch  20\n",
      "Epoch  30\n",
      "Epoch  40\n",
      "Epoch  50\n",
      "Epoch  60\n",
      "Epoch  70\n",
      "Epoch  80\n",
      "Epoch  90\n",
      "Epoch  100\n",
      "Epoch  110\n",
      "Epoch  120\n",
      "delta_label  0.0008352474420547087 < tol  0.001\n",
      "Reach tolerance threshold. Stopping training.\n",
      "Total epoch: 127\n",
      "Calculateing adj matrix using xy only...\n",
      "Processing Br5292_151510\n",
      "Calculateing adj matrix using histology image...\n",
      "Var of c0,c1,c2 =  116.88603475658425 267.75334898637584 77.61160052003804\n",
      "Var of x,y,z =  1281.5097913025436 448.86847099807784 1281.5097913025438\n",
      "Run 1: l [0.01, 1000], p [0.0, 4619.083805075803]\n",
      "Run 2: l [0.01, 500.005], p [0.0, 4577.82177734375]\n",
      "Run 3: l [0.01, 250.0075], p [0.0, 4419.732421875]\n",
      "Run 4: l [0.01, 125.00874999999999], p [0.0, 3881.545654296875]\n",
      "Run 5: l [0.01, 62.509375], p [0.0, 2602.369873046875]\n",
      "Run 6: l [0.01, 31.2596875], p [0.0, 1103.5224609375]\n",
      "Run 7: l [0.01, 15.63484375], p [0.0, 305.1445007324219]\n",
      "Run 8: l [0.01, 7.822421875], p [0.0, 64.29341888427734]\n",
      "Run 9: l [0.01, 3.9162109375], p [0.0, 11.02873420715332]\n",
      "Run 10: l [0.01, 1.9631054687499998], p [0.0, 1.6395759582519531]\n",
      "Run 11: l [0.9865527343749999, 1.9631054687499998], p [0.179815411567688, 1.6395759582519531]\n",
      "Run 12: l [0.9865527343749999, 1.4748291015625], p [0.179815411567688, 0.7000995874404907]\n",
      "Run 13: l [1.23069091796875, 1.4748291015625], p [0.39277541637420654, 0.7000995874404907]\n",
      "Run 14: l [1.23069091796875, 1.352760009765625], p [0.39277541637420654, 0.5339831113815308]\n",
      "Run 15: l [1.2917254638671873, 1.352760009765625], p [0.46034789085388184, 0.5339831113815308]\n",
      "recommended l =  1.3222427368164062\n",
      "Start at res =  0.4 step =  0.1\n",
      "Initializing cluster centers with louvain, resolution =  0.4\n",
      "Epoch  0\n",
      "Res =  0.4 Num of clusters =  5\n",
      "Initializing cluster centers with louvain, resolution =  0.5\n",
      "Epoch  0\n",
      "Res =  0.5 Num of clusters =  6\n",
      "Res changed to 0.5\n",
      "Initializing cluster centers with louvain, resolution =  0.6\n",
      "Epoch  0\n",
      "Res =  0.6 Num of clusters =  6\n",
      "Res changed to 0.6\n",
      "Initializing cluster centers with louvain, resolution =  0.7\n",
      "Epoch  0\n",
      "Res =  0.7 Num of clusters =  6\n",
      "Res changed to 0.7\n",
      "Initializing cluster centers with louvain, resolution =  0.7999999999999999\n",
      "Epoch  0\n",
      "Res =  0.7999999999999999 Num of clusters =  6\n",
      "Res changed to 0.7999999999999999\n",
      "Initializing cluster centers with louvain, resolution =  0.8999999999999999\n",
      "Epoch  0\n",
      "Res =  0.8999999999999999 Num of clusters =  8\n",
      "Step changed to 0.05\n",
      "Initializing cluster centers with louvain, resolution =  0.85\n",
      "Epoch  0\n",
      "Res =  0.85 Num of clusters =  6\n",
      "Res changed to 0.85\n",
      "Initializing cluster centers with louvain, resolution =  0.9\n",
      "Epoch  0\n",
      "Res =  0.9 Num of clusters =  8\n",
      "Step changed to 0.025\n",
      "Initializing cluster centers with louvain, resolution =  0.875\n",
      "Epoch  0\n",
      "Res =  0.875 Num of clusters =  8\n",
      "Step changed to 0.0125\n",
      "Initializing cluster centers with louvain, resolution =  0.8624999999999999\n",
      "Epoch  0\n",
      "Res =  0.8624999999999999 Num of clusters =  7\n",
      "recommended res =  0.8624999999999999\n",
      "Initializing cluster centers with louvain, resolution =  0.8624999999999999\n",
      "Epoch  0\n",
      "Epoch  10\n",
      "Epoch  20\n",
      "Epoch  30\n",
      "Epoch  40\n",
      "Epoch  50\n",
      "Epoch  60\n",
      "Epoch  70\n",
      "delta_label  0.0006473888649115235 < tol  0.001\n",
      "Reach tolerance threshold. Stopping training.\n",
      "Total epoch: 76\n",
      "Calculateing adj matrix using xy only...\n",
      "Processing Br5595_151669\n",
      "Calculateing adj matrix using histology image...\n",
      "Var of c0,c1,c2 =  100.23086214321192 337.8538450184181 80.03365116873323\n",
      "Var of x,y,z =  1078.536154917275 340.0698238839131 1078.5361549172749\n",
      "Run 1: l [0.01, 1000], p [0.0, 3650.8786008844017]\n",
      "Run 2: l [0.01, 500.005], p [0.0, 3623.7607421875]\n",
      "Run 3: l [0.01, 250.0075], p [0.0, 3518.8681640625]\n",
      "Run 4: l [0.01, 125.00874999999999], p [0.0, 3150.131103515625]\n",
      "Run 5: l [0.01, 62.509375], p [0.0, 2209.327880859375]\n",
      "Run 6: l [0.01, 31.2596875], p [0.0, 1042.6387939453125]\n",
      "Run 7: l [0.01, 15.63484375], p [0.0, 338.7726745605469]\n",
      "Run 8: l [0.01, 7.822421875], p [0.0, 77.53675842285156]\n",
      "Run 9: l [0.01, 3.9162109375], p [0.0, 14.540018081665039]\n",
      "Run 10: l [0.01, 1.9631054687499998], p [0.0, 2.3441967964172363]\n",
      "Run 11: l [0.9865527343749999, 1.9631054687499998], p [0.2656189203262329, 2.3441967964172363]\n",
      "Run 12: l [0.9865527343749999, 1.4748291015625], p [0.2656189203262329, 1.0199875831604004]\n",
      "Run 13: l [0.9865527343749999, 1.23069091796875], p [0.2656189203262329, 0.576623797416687]\n",
      "Run 14: l [1.108621826171875, 1.23069091796875], p [0.4052586555480957, 0.576623797416687]\n",
      "Run 15: l [1.1696563720703126, 1.23069091796875], p [0.4868898391723633, 0.576623797416687]\n",
      "Run 16: l [1.1696563720703126, 1.2001736450195313], p [0.4868898391723633, 0.5307332277297974]\n",
      "recommended l =  1.1849150085449218\n",
      "Start at res =  0.4 step =  0.1\n",
      "Initializing cluster centers with louvain, resolution =  0.4\n",
      "Epoch  0\n",
      "Res =  0.4 Num of clusters =  3\n",
      "Initializing cluster centers with louvain, resolution =  0.5\n",
      "Epoch  0\n",
      "Res =  0.5 Num of clusters =  4\n",
      "Res changed to 0.5\n",
      "Initializing cluster centers with louvain, resolution =  0.6\n",
      "Epoch  0\n",
      "Res =  0.6 Num of clusters =  7\n",
      "Step changed to 0.05\n",
      "Initializing cluster centers with louvain, resolution =  0.55\n",
      "Epoch  0\n",
      "Res =  0.55 Num of clusters =  6\n",
      "Step changed to 0.025\n",
      "Initializing cluster centers with louvain, resolution =  0.525\n",
      "Epoch  0\n",
      "Res =  0.525 Num of clusters =  5\n",
      "recommended res =  0.525\n",
      "Initializing cluster centers with louvain, resolution =  0.525\n",
      "Epoch  0\n",
      "Epoch  10\n",
      "Epoch  20\n",
      "Epoch  30\n",
      "Epoch  40\n",
      "Epoch  50\n",
      "Epoch  60\n",
      "Epoch  70\n",
      "delta_label  0.0005462988254575252 < tol  0.001\n",
      "Reach tolerance threshold. Stopping training.\n",
      "Total epoch: 79\n",
      "Calculateing adj matrix using xy only...\n",
      "Processing Br5595_151670\n",
      "Calculateing adj matrix using histology image...\n",
      "Var of c0,c1,c2 =  121.03381450869409 316.2601354125137 59.683037373252404\n",
      "Var of x,y,z =  1101.174662904654 298.2898662014167 1101.1746629046545\n",
      "Run 1: l [0.01, 1000], p [0.0, 3488.2730149559443]\n",
      "Run 2: l [0.01, 500.005], p [0.0, 3462.333984375]\n",
      "Run 3: l [0.01, 250.0075], p [0.0, 3362.085693359375]\n",
      "Run 4: l [0.01, 125.00874999999999], p [0.0, 3010.851806640625]\n",
      "Run 5: l [0.01, 62.509375], p [0.0, 2124.27783203125]\n",
      "Run 6: l [0.01, 31.2596875], p [0.0, 1032.4498291015625]\n",
      "Run 7: l [0.01, 15.63484375], p [0.0, 343.41082763671875]\n",
      "Run 8: l [0.01, 7.822421875], p [0.0, 79.83930206298828]\n",
      "Run 9: l [0.01, 3.9162109375], p [0.0, 15.287084579467773]\n",
      "Run 10: l [0.01, 1.9631054687499998], p [0.0, 2.5060837268829346]\n",
      "Run 11: l [0.9865527343749999, 1.9631054687499998], p [0.2777498960494995, 2.5060837268829346]\n",
      "Run 12: l [0.9865527343749999, 1.4748291015625], p [0.2777498960494995, 1.0853335857391357]\n",
      "Run 13: l [0.9865527343749999, 1.23069091796875], p [0.2777498960494995, 0.6091245412826538]\n",
      "Run 14: l [1.108621826171875, 1.23069091796875], p [0.4260367155075073, 0.6091245412826538]\n",
      "Run 15: l [1.108621826171875, 1.1696563720703126], p [0.4260367155075073, 0.5131330490112305]\n",
      "Run 16: l [1.1391390991210937, 1.1696563720703126], p [0.4684847593307495, 0.5131330490112305]\n",
      "recommended l =  1.1543977355957031\n",
      "Start at res =  0.4 step =  0.1\n",
      "Initializing cluster centers with louvain, resolution =  0.4\n",
      "Epoch  0\n",
      "Res =  0.4 Num of clusters =  4\n",
      "Initializing cluster centers with louvain, resolution =  0.5\n",
      "Epoch  0\n",
      "Res =  0.5 Num of clusters =  6\n",
      "Step changed to 0.05\n",
      "Initializing cluster centers with louvain, resolution =  0.45\n",
      "Epoch  0\n",
      "Res =  0.45 Num of clusters =  5\n",
      "recommended res =  0.45\n",
      "Initializing cluster centers with louvain, resolution =  0.45\n",
      "Epoch  0\n",
      "Epoch  10\n",
      "Epoch  20\n",
      "Epoch  30\n",
      "Epoch  40\n",
      "Epoch  50\n",
      "Epoch  60\n",
      "Epoch  70\n",
      "Epoch  80\n",
      "delta_label  0.0005717552887364208 < tol  0.001\n",
      "Reach tolerance threshold. Stopping training.\n",
      "Total epoch: 88\n",
      "Calculateing adj matrix using xy only...\n",
      "Processing Br5595_151671\n",
      "Calculateing adj matrix using histology image...\n",
      "Var of c0,c1,c2 =  75.55798518644477 209.70027180991067 38.939737645325806\n",
      "Var of x,y,z =  1200.9333700368811 385.02424115414897 1200.9333700368811\n",
      "Run 1: l [0.01, 1000], p [0.0, 4097.581773525745]\n",
      "Run 2: l [0.01, 500.005], p [0.0, 4063.753662109375]\n",
      "Run 3: l [0.01, 250.0075], p [0.0, 3934.537841796875]\n",
      "Run 4: l [0.01, 125.00874999999999], p [0.0, 3499.43359375]\n",
      "Run 5: l [0.01, 62.509375], p [0.0, 2485.08544921875]\n",
      "Run 6: l [0.01, 31.2596875], p [0.0, 1196.14990234375]\n",
      "Run 7: l [0.01, 15.63484375], p [0.0, 344.65362548828125]\n",
      "Run 8: l [0.01, 7.822421875], p [0.0, 66.68467712402344]\n",
      "Run 9: l [0.01, 3.9162109375], p [0.0, 10.531600952148438]\n",
      "Run 10: l [0.01, 1.9631054687499998], p [0.0, 1.4576671123504639]\n",
      "Run 11: l [0.9865527343749999, 1.9631054687499998], p [0.14400064945220947, 1.4576671123504639]\n",
      "Run 12: l [0.9865527343749999, 1.4748291015625], p [0.14400064945220947, 0.5986154079437256]\n",
      "Run 13: l [1.23069091796875, 1.4748291015625], p [0.32640981674194336, 0.5986154079437256]\n",
      "Run 14: l [1.352760009765625, 1.4748291015625], p [0.45054948329925537, 0.5986154079437256]\n",
      "Run 15: l [1.352760009765625, 1.4137945556640625], p [0.45054948329925537, 0.5215022563934326]\n",
      "Run 16: l [1.3832772827148436, 1.4137945556640625], p [0.48526692390441895, 0.5215022563934326]\n",
      "recommended l =  1.398535919189453\n",
      "Start at res =  0.4 step =  0.1\n",
      "Initializing cluster centers with louvain, resolution =  0.4\n",
      "Epoch  0\n",
      "Res =  0.4 Num of clusters =  4\n",
      "Initializing cluster centers with louvain, resolution =  0.5\n",
      "Epoch  0\n",
      "Res =  0.5 Num of clusters =  5\n",
      "recommended res =  0.5\n",
      "Initializing cluster centers with louvain, resolution =  0.5\n",
      "Epoch  0\n",
      "Epoch  10\n",
      "Epoch  20\n",
      "Epoch  30\n",
      "Epoch  40\n",
      "Epoch  50\n",
      "Epoch  60\n",
      "Epoch  70\n",
      "Epoch  80\n",
      "Epoch  90\n",
      "Epoch  100\n",
      "delta_label  0.00024330900243309004 < tol  0.001\n",
      "Reach tolerance threshold. Stopping training.\n",
      "Total epoch: 103\n",
      "Calculateing adj matrix using xy only...\n",
      "Processing Br5595_151672\n",
      "Calculateing adj matrix using histology image...\n",
      "Var of c0,c1,c2 =  114.64310906516253 288.0502933634853 58.59128676561912\n",
      "Var of x,y,z =  1180.8210216668813 374.61722984635765 1180.821021666881\n",
      "Run 1: l [0.01, 1000], p [0.0, 4003.043401517705]\n",
      "Run 2: l [0.01, 500.005], p [0.0, 3970.52490234375]\n",
      "Run 3: l [0.01, 250.0075], p [0.0, 3845.517578125]\n",
      "Run 4: l [0.01, 125.00874999999999], p [0.0, 3415.693603515625]\n",
      "Run 5: l [0.01, 62.509375], p [0.0, 2379.66162109375]\n",
      "Run 6: l [0.01, 31.2596875], p [0.0, 1135.8570556640625]\n",
      "Run 7: l [0.01, 15.63484375], p [0.0, 350.8916015625]\n",
      "Run 8: l [0.01, 7.822421875], p [0.0, 75.0586929321289]\n",
      "Run 9: l [0.01, 3.9162109375], p [0.0, 13.156631469726562]\n",
      "Run 10: l [0.01, 1.9631054687499998], p [0.0, 1.997056007385254]\n",
      "Run 11: l [0.9865527343749999, 1.9631054687499998], p [0.2150135040283203, 1.997056007385254]\n",
      "Run 12: l [0.9865527343749999, 1.4748291015625], p [0.2150135040283203, 0.8504630327224731]\n",
      "Run 13: l [1.23069091796875, 1.4748291015625], p [0.4742389917373657, 0.8504630327224731]\n",
      "Run 14: l [1.23069091796875, 1.352760009765625], p [0.4742389917373657, 0.647018551826477]\n",
      "Run 15: l [1.23069091796875, 1.2917254638671873], p [0.4742389917373657, 0.5568807125091553]\n",
      "Run 16: l [1.23069091796875, 1.2612081909179687], p [0.4742389917373657, 0.5146338939666748]\n",
      "recommended l =  1.2459495544433592\n",
      "Start at res =  0.4 step =  0.1\n",
      "Initializing cluster centers with louvain, resolution =  0.4\n",
      "Epoch  0\n",
      "Res =  0.4 Num of clusters =  3\n",
      "Initializing cluster centers with louvain, resolution =  0.5\n",
      "Epoch  0\n",
      "Res =  0.5 Num of clusters =  5\n",
      "recommended res =  0.5\n",
      "Initializing cluster centers with louvain, resolution =  0.5\n",
      "Epoch  0\n",
      "Epoch  10\n",
      "Epoch  20\n",
      "Epoch  30\n",
      "Epoch  40\n",
      "Epoch  50\n",
      "Epoch  60\n",
      "Epoch  70\n",
      "Epoch  80\n",
      "delta_label  0.00099626400996264 < tol  0.001\n",
      "Reach tolerance threshold. Stopping training.\n",
      "Total epoch: 85\n",
      "Calculateing adj matrix using xy only...\n",
      "Processing Br8100_151673\n",
      "Calculateing adj matrix using histology image...\n",
      "Var of c0,c1,c2 =  95.62291694867722 581.9646063440777 74.15680250956527\n",
      "Var of x,y,z =  942.1320230796399 390.6990007265343 942.1320230796399\n",
      "Run 1: l [0.01, 1000], p [0.0, 3629.7390523246568]\n",
      "Run 2: l [0.01, 500.005], p [0.0, 3605.1669921875]\n",
      "Run 3: l [0.01, 250.0075], p [0.0, 3509.94091796875]\n",
      "Run 4: l [0.01, 125.00874999999999], p [0.0, 3172.689453125]\n",
      "Run 5: l [0.01, 62.509375], p [0.0, 2288.652587890625]\n",
      "Run 6: l [0.01, 31.2596875], p [0.0, 1115.2354736328125]\n",
      "Run 7: l [0.01, 15.63484375], p [0.0, 370.28887939453125]\n",
      "Run 8: l [0.01, 7.822421875], p [0.0, 89.5790023803711]\n",
      "Run 9: l [0.01, 3.9162109375], p [0.0, 17.459369659423828]\n",
      "Run 10: l [0.01, 1.9631054687499998], p [0.0, 2.8228759765625]\n",
      "Run 11: l [0.9865527343749999, 1.9631054687499998], p [0.31799256801605225, 2.8228759765625]\n",
      "Run 12: l [0.9865527343749999, 1.4748291015625], p [0.31799256801605225, 1.2253105640411377]\n",
      "Run 13: l [0.9865527343749999, 1.23069091796875], p [0.31799256801605225, 0.691605806350708]\n",
      "Run 14: l [1.108621826171875, 1.23069091796875], p [0.48566174507141113, 0.691605806350708]\n",
      "Run 15: l [1.108621826171875, 1.1696563720703126], p [0.48566174507141113, 0.5837382078170776]\n",
      "Run 16: l [1.108621826171875, 1.1391390991210937], p [0.48566174507141113, 0.5334899425506592]\n",
      "recommended l =  1.1238804626464844\n",
      "Start at res =  0.4 step =  0.1\n",
      "Initializing cluster centers with louvain, resolution =  0.4\n",
      "Epoch  0\n",
      "Res =  0.4 Num of clusters =  4\n",
      "Initializing cluster centers with louvain, resolution =  0.5\n",
      "Epoch  0\n",
      "Res =  0.5 Num of clusters =  6\n",
      "Res changed to 0.5\n",
      "Initializing cluster centers with louvain, resolution =  0.6\n",
      "Epoch  0\n",
      "Res =  0.6 Num of clusters =  6\n",
      "Res changed to 0.6\n",
      "Initializing cluster centers with louvain, resolution =  0.7\n",
      "Epoch  0\n",
      "Res =  0.7 Num of clusters =  8\n",
      "Step changed to 0.05\n",
      "Initializing cluster centers with louvain, resolution =  0.65\n",
      "Epoch  0\n",
      "Res =  0.65 Num of clusters =  9\n",
      "Step changed to 0.025\n",
      "Initializing cluster centers with louvain, resolution =  0.625\n",
      "Epoch  0\n",
      "Res =  0.625 Num of clusters =  7\n",
      "recommended res =  0.625\n",
      "Initializing cluster centers with louvain, resolution =  0.625\n",
      "Epoch  0\n",
      "Epoch  10\n",
      "Epoch  20\n",
      "Epoch  30\n",
      "Epoch  40\n",
      "Epoch  50\n",
      "Epoch  60\n",
      "Epoch  70\n",
      "Epoch  80\n",
      "Epoch  90\n",
      "Epoch  100\n",
      "delta_label  0.0008244023083264633 < tol  0.001\n",
      "Reach tolerance threshold. Stopping training.\n",
      "Total epoch: 103\n",
      "Calculateing adj matrix using xy only...\n",
      "Processing Br8100_151674\n",
      "Calculateing adj matrix using histology image...\n",
      "Var of c0,c1,c2 =  111.39574416328388 601.3615262083176 67.55499930498358\n",
      "Var of x,y,z =  920.7962138115173 408.49702818834794 920.7962138115174\n",
      "Run 1: l [0.01, 1000], p [0.0, 3663.7520994905253]\n",
      "Run 2: l [0.01, 500.005], p [0.0, 3639.207763671875]\n",
      "Run 3: l [0.01, 250.0075], p [0.0, 3543.926513671875]\n",
      "Run 4: l [0.01, 125.00874999999999], p [0.0, 3204.1669921875]\n",
      "Run 5: l [0.01, 62.509375], p [0.0, 2291.34765625]\n",
      "Run 6: l [0.01, 31.2596875], p [0.0, 1039.6795654296875]\n",
      "Run 7: l [0.01, 15.63484375], p [0.0, 324.6616516113281]\n",
      "Run 8: l [0.01, 7.822421875], p [0.0, 80.1252212524414]\n",
      "Run 9: l [0.01, 3.9162109375], p [0.0, 16.391626358032227]\n",
      "Run 10: l [0.01, 1.9631054687499998], p [0.0, 2.8255667686462402]\n",
      "Run 11: l [0.9865527343749999, 1.9631054687499998], p [0.33460450172424316, 2.8255667686462402]\n",
      "Run 12: l [0.9865527343749999, 1.4748291015625], p [0.33460450172424316, 1.2570574283599854]\n",
      "Run 13: l [0.9865527343749999, 1.23069091796875], p [0.33460450172424316, 0.7187591791152954]\n",
      "recommended l =  1.108621826171875\n",
      "Start at res =  0.4 step =  0.1\n",
      "Initializing cluster centers with louvain, resolution =  0.4\n",
      "Epoch  0\n",
      "Res =  0.4 Num of clusters =  5\n",
      "Initializing cluster centers with louvain, resolution =  0.5\n",
      "Epoch  0\n",
      "Res =  0.5 Num of clusters =  6\n",
      "Res changed to 0.5\n",
      "Initializing cluster centers with louvain, resolution =  0.6\n",
      "Epoch  0\n",
      "Res =  0.6 Num of clusters =  7\n",
      "recommended res =  0.6\n",
      "Initializing cluster centers with louvain, resolution =  0.6\n",
      "Epoch  0\n",
      "Epoch  10\n",
      "Epoch  20\n",
      "Epoch  30\n",
      "Epoch  40\n",
      "Epoch  50\n",
      "Epoch  60\n",
      "Epoch  70\n",
      "Epoch  80\n",
      "delta_label  0.0002722570106180234 < tol  0.001\n",
      "Reach tolerance threshold. Stopping training.\n",
      "Total epoch: 85\n",
      "Calculateing adj matrix using xy only...\n",
      "Processing Br8100_151675\n",
      "Calculateing adj matrix using histology image...\n",
      "Var of c0,c1,c2 =  107.46704917299442 424.21344722384123 61.97015788788262\n",
      "Var of x,y,z =  1024.3871538800652 378.9827598046389 1024.3871538800654\n",
      "Run 1: l [0.01, 1000], p [0.0, 3582.2996941343517]\n",
      "Run 2: l [0.01, 500.005], p [0.0, 3556.43994140625]\n",
      "Run 3: l [0.01, 250.0075], p [0.0, 3456.486572265625]\n",
      "Run 4: l [0.01, 125.00874999999999], p [0.0, 3105.644775390625]\n",
      "Run 5: l [0.01, 62.509375], p [0.0, 2202.875244140625]\n",
      "Run 6: l [0.01, 31.2596875], p [0.0, 1017.9781494140625]\n",
      "Run 7: l [0.01, 15.63484375], p [0.0, 318.6366271972656]\n",
      "Run 8: l [0.01, 7.822421875], p [0.0, 76.90209197998047]\n",
      "Run 9: l [0.01, 3.9162109375], p [0.0, 15.373720169067383]\n",
      "Run 10: l [0.01, 1.9631054687499998], p [0.0, 2.5614852905273438]\n",
      "Run 11: l [0.9865527343749999, 1.9631054687499998], p [0.2892416715621948, 2.5614852905273438]\n",
      "Run 12: l [0.9865527343749999, 1.4748291015625], p [0.2892416715621948, 1.1174530982971191]\n",
      "Run 13: l [0.9865527343749999, 1.23069091796875], p [0.2892416715621948, 0.6306978464126587]\n",
      "Run 14: l [1.108621826171875, 1.23069091796875], p [0.4424600601196289, 0.6306978464126587]\n",
      "Run 15: l [1.108621826171875, 1.1696563720703126], p [0.4424600601196289, 0.5321156978607178]\n",
      "Run 16: l [1.1391390991210937, 1.1696563720703126], p [0.48618268966674805, 0.5321156978607178]\n",
      "recommended l =  1.1543977355957031\n",
      "Start at res =  0.4 step =  0.1\n",
      "Initializing cluster centers with louvain, resolution =  0.4\n",
      "Epoch  0\n",
      "Res =  0.4 Num of clusters =  4\n",
      "Initializing cluster centers with louvain, resolution =  0.5\n",
      "Epoch  0\n",
      "Res =  0.5 Num of clusters =  6\n",
      "Res changed to 0.5\n",
      "Initializing cluster centers with louvain, resolution =  0.6\n",
      "Epoch  0\n",
      "Res =  0.6 Num of clusters =  6\n",
      "Res changed to 0.6\n",
      "Initializing cluster centers with louvain, resolution =  0.7\n",
      "Epoch  0\n",
      "Res =  0.7 Num of clusters =  6\n",
      "Res changed to 0.7\n",
      "Initializing cluster centers with louvain, resolution =  0.7999999999999999\n",
      "Epoch  0\n",
      "Res =  0.7999999999999999 Num of clusters =  8\n",
      "Step changed to 0.05\n",
      "Initializing cluster centers with louvain, resolution =  0.75\n",
      "Epoch  0\n",
      "Res =  0.75 Num of clusters =  6\n",
      "Res changed to 0.75\n",
      "Initializing cluster centers with louvain, resolution =  0.8\n",
      "Epoch  0\n",
      "Res =  0.8 Num of clusters =  8\n",
      "Step changed to 0.025\n",
      "Initializing cluster centers with louvain, resolution =  0.775\n",
      "Epoch  0\n",
      "Res =  0.775 Num of clusters =  6\n",
      "Res changed to 0.775\n",
      "Initializing cluster centers with louvain, resolution =  0.8\n",
      "Epoch  0\n",
      "Res =  0.8 Num of clusters =  8\n",
      "Step changed to 0.0125\n",
      "Initializing cluster centers with louvain, resolution =  0.7875\n",
      "Epoch  0\n",
      "Res =  0.7875 Num of clusters =  8\n",
      "Step changed to 0.00625\n",
      "Initializing cluster centers with louvain, resolution =  0.78125\n",
      "Epoch  0\n",
      "Res =  0.78125 Num of clusters =  6\n",
      "Res changed to 0.78125\n",
      "Initializing cluster centers with louvain, resolution =  0.7875\n",
      "Epoch  0\n",
      "Res =  0.7875 Num of clusters =  8\n",
      "Step changed to 0.003125\n",
      "Initializing cluster centers with louvain, resolution =  0.784375\n",
      "Epoch  0\n",
      "Res =  0.784375 Num of clusters =  6\n",
      "Res changed to 0.784375\n",
      "Exact resolution not found\n",
      "Recommended res =  0.784375\n",
      "Initializing cluster centers with louvain, resolution =  0.784375\n",
      "Epoch  0\n",
      "Epoch  10\n",
      "Epoch  20\n",
      "Epoch  30\n",
      "Epoch  40\n",
      "Epoch  50\n",
      "Epoch  60\n",
      "Epoch  70\n",
      "Epoch  80\n",
      "Epoch  90\n",
      "Epoch  100\n",
      "Epoch  110\n",
      "Epoch  120\n",
      "Epoch  130\n",
      "Epoch  140\n",
      "Epoch  150\n",
      "delta_label  0.00027839643652561246 < tol  0.001\n",
      "Reach tolerance threshold. Stopping training.\n",
      "Total epoch: 151\n",
      "Calculateing adj matrix using xy only...\n",
      "Processing Br8100_151676\n",
      "Calculateing adj matrix using histology image...\n",
      "Var of c0,c1,c2 =  111.729715722502 389.90131980864106 55.72587548303189\n",
      "Var of x,y,z =  933.3009317885663 359.0939566139864 933.3009317885665\n",
      "Run 1: l [0.01, 1000], p [0.0, 3451.314709506626]\n",
      "Run 2: l [0.01, 500.005], p [0.0, 3428.4453125]\n",
      "Run 3: l [0.01, 250.0075], p [0.0, 3339.67919921875]\n",
      "Run 4: l [0.01, 125.00874999999999], p [0.0, 3023.3212890625]\n",
      "Run 5: l [0.01, 62.509375], p [0.0, 2175.862548828125]\n",
      "Run 6: l [0.01, 31.2596875], p [0.0, 1018.6404418945312]\n",
      "Run 7: l [0.01, 15.63484375], p [0.0, 325.3836975097656]\n",
      "Run 8: l [0.01, 7.822421875], p [0.0, 78.885498046875]\n",
      "Run 9: l [0.01, 3.9162109375], p [0.0, 15.651283264160156]\n",
      "Run 10: l [0.01, 1.9631054687499998], p [0.0, 2.612065076828003]\n",
      "Run 11: l [0.9865527343749999, 1.9631054687499998], p [0.30002427101135254, 2.612065076828003]\n",
      "Run 12: l [0.9865527343749999, 1.4748291015625], p [0.30002427101135254, 1.1461548805236816]\n",
      "Run 13: l [0.9865527343749999, 1.23069091796875], p [0.30002427101135254, 0.650269627571106]\n",
      "Run 14: l [1.108621826171875, 1.23069091796875], p [0.4575594663619995, 0.650269627571106]\n",
      "Run 15: l [1.108621826171875, 1.1696563720703126], p [0.4575594663619995, 0.5494418144226074]\n",
      "recommended l =  1.1391390991210937\n",
      "Start at res =  0.4 step =  0.1\n",
      "Initializing cluster centers with louvain, resolution =  0.4\n",
      "Epoch  0\n",
      "Res =  0.4 Num of clusters =  4\n",
      "Initializing cluster centers with louvain, resolution =  0.5\n",
      "Epoch  0\n",
      "Res =  0.5 Num of clusters =  4\n",
      "Res changed to 0.5\n",
      "Initializing cluster centers with louvain, resolution =  0.6\n",
      "Epoch  0\n",
      "Res =  0.6 Num of clusters =  7\n",
      "recommended res =  0.6\n",
      "Initializing cluster centers with louvain, resolution =  0.6\n",
      "Epoch  0\n",
      "Epoch  10\n",
      "Epoch  20\n",
      "Epoch  30\n",
      "Epoch  40\n",
      "Epoch  50\n",
      "Epoch  60\n",
      "Epoch  70\n",
      "delta_label  0.0008670520231213873 < tol  0.001\n",
      "Reach tolerance threshold. Stopping training.\n",
      "Total epoch: 73\n",
      "Calculateing adj matrix using xy only...\n"
     ]
    }
   ],
   "source": [
    "for _, sample in metadata.iterrows():\n",
    "    print(\"Processing \" + sample.directory)\n",
    "\n",
    "    sample_dir = data_dir / sample.directory\n",
    "    out_dir = result_dir / sample.directory\n",
    "    adata = get_anndata(sample_dir)\n",
    "\n",
    "    adj = spg.calculate_adj_matrix(\n",
    "        adata.obs[\"col\"],\n",
    "        adata.obs[\"row\"],\n",
    "        adata.obsm[\"spatial_pixel\"][:, 0],\n",
    "        adata.obsm[\"spatial_pixel\"][:, 1],\n",
    "        image=adata.uns[\"image\"],\n",
    "        alpha=alpha,\n",
    "        beta=beta,\n",
    "        histology=True,\n",
    "    )\n",
    "\n",
    "    spg.prefilter_genes(adata, min_cells=3)\n",
    "    spg.prefilter_specialgenes(adata)\n",
    "    adata.X = adata.X.astype(\"float64\")\n",
    "    sc.pp.normalize_per_cell(adata)\n",
    "    sc.pp.log1p(adata)\n",
    "\n",
    "    clf = spg.SpaGCN()\n",
    "\n",
    "    # Find the l value given p\n",
    "    l = spg.search_l(p, adj)\n",
    "    clf.set_l(l)\n",
    "\n",
    "    # Search for suitable resolution\n",
    "    res = spg.search_res(\n",
    "        adata, adj, l, sample.n_clusters, r_seed=seed, t_seed=seed, n_seed=seed\n",
    "    )\n",
    "\n",
    "    clf.train(adata, adj, res=res)\n",
    "\n",
    "    y_pred, prob = clf.predict()\n",
    "\n",
    "    adata.obs[\"cluster\"] = pd.Series(y_pred, index=adata.obs_names, dtype=\"category\")\n",
    "\n",
    "    # Do cluster refinement(optional)\n",
    "    adj_2d = spg.calculate_adj_matrix(\n",
    "        adata.obs[\"col\"], adata.obs[\"row\"], histology=False\n",
    "    )\n",
    "\n",
    "    refined_pred = spg.refine(\n",
    "        sample_id=adata.obs_names.tolist(),\n",
    "        pred=adata.obs[\"cluster\"].tolist(),\n",
    "        dis=adj_2d,\n",
    "        shape=\"hexagon\",\n",
    "    )\n",
    "\n",
    "    adata.obs[\"refined_cluster\"] = pd.Series(\n",
    "        refined_pred, index=adata.obs_names, dtype=\"category\"\n",
    "    )\n",
    "\n",
    "    # label_df = adata.obs[[\"cluster\"]]\n",
    "    label_df_ref = adata.obs[[\"refined_cluster\"]]\n",
    "\n",
    "    ## Write output\n",
    "    out_dir.mkdir(parents=True, exist_ok=True)\n",
    "\n",
    "    label_df_ref.columns = [\"label\"]\n",
    "    label_df_ref.to_csv(out_dir / \"SpaGCN.tsv\", sep=\"\\t\", index_label=\"\")"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python [conda env:muellni-spagcn]",
   "language": "python",
   "name": "conda-env-muellni-spagcn-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.18"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}