Repository: pykale/pykale (language: python)

Issue:
Construct an example on brain fMRI classification
Source: https://github.com/kundaMwiza/fMRI-site-adaptation
Attempt: https://github.com/shan3333/Mwiza_Pykale


Discussion / hints:

Fix patch:
diff --git a/.gitignore b/.gitignore
index 1e41ec778..2cc314717 100644
--- a/.gitignore
+++ b/.gitignore
@@ -5,6 +5,7 @@
 *.pth
 *.zip
 *.ckpt
+*.json
 
 # Distribution / packaging
 build/
@@ -23,8 +24,10 @@ examples/data/
 examples/*/data/
 examples/*/outputs/
 examples/*/lightning_logs/
+examples/*/tb_logs/
+examples/*/logs/
 examples/*/demo*/
-examples/action_dann_lightn/configs_xianyuan/
+examples/*/configs_experiment/
 
 # Logs
 log-*.txt
diff --git a/docs/requirements.txt b/docs/requirements.txt
index 591b9bb7f..a1dbacb05 100644
--- a/docs/requirements.txt
+++ b/docs/requirements.txt
@@ -4,9 +4,9 @@
 # Learn from https://github.com/rusty1s/pytorch_geometric/blob/master/docs/requirements.txt
 # Remove if switching to DGL
 cython
-https://download.pytorch.org/whl/cpu/torch-1.7.0%2Bcpu-cp37-cp37m-linux_x86_64.whl
-https://pytorch-geometric.com/whl/torch-1.7.0+cpu/torch_scatter-2.0.6-cp37-cp37m-linux_x86_64.whl
-https://pytorch-geometric.com/whl/torch-1.7.0+cpu/torch_sparse-0.6.8-cp37-cp37m-linux_x86_64.whl
+https://data.pyg.org/whl/torch-1.9.0%2Bcpu/torch_scatter-2.0.7-cp37-cp37m-linux_x86_64.whl
+https://data.pyg.org/whl/torch-1.9.0%2Bcpu/torch_sparse-0.6.10-cp37-cp37m-linux_x86_64.whl
+https://download.pytorch.org/whl/cpu/torch-1.9.1%2Bcpu-cp37-cp37m-linux_x86_64.whl
 
 ipython<9.0
 
diff --git a/docs/source/kale.embed.rst b/docs/source/kale.embed.rst
index 04cb5b2e3..f8628370c 100644
--- a/docs/source/kale.embed.rst
+++ b/docs/source/kale.embed.rst
@@ -14,14 +14,13 @@ kale.embed.attention\_cnn module
    :undoc-members:
    :show-inheritance:
 
-kale.embed.seq\_nn module
---------------------------------
+kale.embed.factorization module
+----------------------
 
-.. automodule:: kale.embed.seq_nn
+.. automodule:: kale.embed.factorization
    :members:
    :undoc-members:
    :show-inheritance:
-   :exclude-members:
 
 kale.embed.gcn module
 ----------------------------
@@ -50,21 +49,22 @@ kale.embed.image\_cnn module
    :undoc-members:
    :show-inheritance:
 
-kale.embed.mpca module
-----------------------
+kale.embed.positional\_encoding module
+--------------------------------------
 
-.. automodule:: kale.embed.mpca
+.. automodule:: kale.embed.positional_encoding
    :members:
    :undoc-members:
    :show-inheritance:
 
-kale.embed.positional\_encoding module
---------------------------------------
+kale.embed.seq\_nn module
+--------------------------------
 
-.. automodule:: kale.embed.positional_encoding
+.. automodule:: kale.embed.seq_nn
    :members:
    :undoc-members:
    :show-inheritance:
+   :exclude-members:
 
 kale.embed.video\_feature\_extractor module
 -------------------------------------------
diff --git a/examples/autism_detection/__init__.py b/examples/autism_detection/__init__.py
new file mode 100644
index 000000000..e69de29bb
diff --git a/examples/autism_detection/config.py b/examples/autism_detection/config.py
new file mode 100644
index 000000000..683c291a9
--- /dev/null
+++ b/examples/autism_detection/config.py
@@ -0,0 +1,38 @@
+"""
+Default configurations for classification on resting-state fMRI of ABIDE
+"""
+
+from yacs.config import CfgNode as CN
+
+# -----------------------------------------------------------------------------
+# Config definition
+# -----------------------------------------------------------------------------
+
+_C = CN()
+
+# -----------------------------------------------------------------------------
+# Dataset
+# -----------------------------------------------------------------------------
+_C.DATASET = CN()
+_C.DATASET.ROOT = "../data"
+_C.DATASET.PIPELINE = "cpac"  # options: {‘cpac’, ‘css’, ‘dparsf’, ‘niak’}
+_C.DATASET.ATLAS = "rois_cc200"
+# options: {rois_aal, rois_cc200, rois_cc400, rois_dosenbach160, rois_ez, rois_ho, rois_tt}
+_C.DATASET.SITE_IDS = None  # list of site ids to use, if None, use all sites
+_C.DATASET.TARGET = "NYU"  # target site ids, e.g. "UM_1", "UCLA_1", "USM"
+# ---------------------------------------------------------------------------- #
+# Machine learning pipeline
+# ---------------------------------------------------------------------------- #
+_C.MODEL = CN()
+_C.MODEL.KERNEL = "rbf"
+_C.MODEL.ALPHA = 0.01
+_C.MODEL.LAMBDA_ = 1.0
+# ---------------------------------------------------------------------------- #
+# Misc options
+# ---------------------------------------------------------------------------- #
+_C.OUTPUT = CN()
+_C.OUTPUT.ROOT = "./outputs"  # output_dir
+
+
+def get_cfg_defaults():
+    return _C.clone()
diff --git a/examples/autism_detection/configs/tutorial.yaml b/examples/autism_detection/configs/tutorial.yaml
new file mode 100644
index 000000000..0a220b52f
--- /dev/null
+++ b/examples/autism_detection/configs/tutorial.yaml
@@ -0,0 +1,7 @@
+MODEL:
+  KERNEL: "linear"
+  ALPHA: 1.0
+  LAMBDA_: 1.0
+DATASET:
+  ATLAS: "rois_cc200"
+  SITE_IDS: ['NYU', "UM_1", "UCLA_1", "USM"]
diff --git a/examples/autism_detection/tutorial.ipynb b/examples/autism_detection/tutorial.ipynb
new file mode 100644
index 000000000..50bbc06ba
--- /dev/null
+++ b/examples/autism_detection/tutorial.ipynb
@@ -0,0 +1,489 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "id": "c2ce4959",
+   "metadata": {
+    "pycharm": {
+     "name": "#%% md\n"
+    }
+   },
+   "source": [
+    "# PyKale Tutorial: Domain Adaptation (Generalization) for Autism Detection with Multi-site Brain Imaging Data\n",
+    "| [Open in Colab](https://colab.research.google.com/github/sz144/pykale/blob/brain-example/examples/autism_detection/tutorial.ipynb) (click `Runtime` → `Run all (Ctrl+F9)` |  [Launch Binder](https://mybinder.org/v2/gh/pykale/pykale/HEAD?filepath=examples%2Fautism_detection%2Ftutorial.ipynb) (click `Run` → `Run All Cells`) |"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "ed8a9c14",
+   "metadata": {
+    "pycharm": {
+     "name": "#%% md\n"
+    }
+   },
+   "source": [
+    "## Overview\n",
+    "\n",
+    "- Pre-processing:\n",
+    "    - [Data loading](#Data-Preparation)\n",
+    "    - [Construct brain networks](#Extracting-Brain-Networks-Features)\n",
+    "- Machine learning pipeline:\n",
+    "    - [Baseline: Ridge classifier](#Baseline)\n",
+    "    - [Domain adaptation](#Domain-Adaptation)\n",
+    "    - [Domain generalization](#Domain-Generalization)\n",
+    "\n",
+    "**Reference:**\n",
+    "\n",
+    "[1] Cameron Craddock, Yassine Benhajali, Carlton Chu, Francois Chouinard, Alan Evans, András Jakab, Budhachandra Singh Khundrakpam, John David Lewis, Qingyang Li, Michael Milham, Chaogan Yan, Pierre Bellec (2013). The Neuro Bureau Preprocessing Initiative: open sharing of preprocessed neuroimaging data and derivatives. In *Neuroinformatics 2013*, Stockholm, Sweden.\n",
+    "\n",
+    "[2] Abraham, A., Pedregosa, F., Eickenberg, M., Gervais, P., Mueller, A., Kossaifi, J., ... & Varoquaux, G. (2014). Machine learning for neuroimaging with scikit-learn. *Frontiers in neuroinformatics*, 14.\n",
+    "\n",
+    "[3] Zhou, S., Li, W., Cox, C.R., & Lu, H. (2020). Side Information Dependence as a Regularizer for Analyzing Human Brain Conditions across Cognitive Experiments. *AAAI 2020*, New York, USA. [[Link](https://ojs.aaai.org//index.php/AAAI/article/view/6179)]\n",
+    "\n",
+    "[4] Zhou, S. (2022). Interpretable Domain-Aware Learning for Neuroimage Classification (Doctoral dissertation, University of Sheffield). [[Link](https://etheses.whiterose.ac.uk/31044/1/PhD_thesis_ShuoZhou_170272834.pdf)]"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "3386eedc",
+   "metadata": {
+    "pycharm": {
+     "name": "#%% md\n"
+    }
+   },
+   "source": [
+    "## Setup"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "2ddd53bb",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "if 'google.colab' in str(get_ipython()):\n",
+    "    print('Running on CoLab')\n",
+    "    !pip uninstall --yes imgaug && pip uninstall --yes albumentations && pip install git+https://github.com/aleju/imgaug.git\n",
+    "    !pip install git+https://github.com/pykale/pykale.git\n",
+    "    !git clone https://github.com/pykale/pykale.git\n",
+    "    %cd pykale/examples/autism_detection\n",
+    "else:\n",
+    "    print('Not running on CoLab')"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "45a65434",
+   "metadata": {
+    "pycharm": {
+     "name": "#%% md\n"
+    }
+   },
+   "source": [
+    "This imports required modules."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "78684d23",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    },
+    "scrolled": true
+   },
+   "outputs": [],
+   "source": [
+    "import os\n",
+    "\n",
+    "%matplotlib inline\n",
+    "import matplotlib.pyplot as plt\n",
+    "import numpy as np\n",
+    "from nilearn.datasets import fetch_abide_pcp\n",
+    "import pandas as pd\n",
+    "from config import get_cfg_defaults\n",
+    "\n",
+    "import sys\n",
+    "\n",
+    "from kale.utils.download import download_file_by_url\n",
+    "from kale.interpret import visualize"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "d5547762",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "cfg_path = \"configs/tutorial.yaml\" # Path to `.yaml` config file\n",
+    "\n",
+    "cfg = get_cfg_defaults()\n",
+    "cfg.merge_from_file(cfg_path)\n",
+    "cfg.freeze()\n",
+    "print(cfg)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "5bef71fc",
+   "metadata": {
+    "pycharm": {
+     "name": "#%% md\n"
+    }
+   },
+   "source": [
+    "## Data Preparation\n",
+    "\n",
+    "### Fetch ABIDE fMRI timeseries"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "f1d1b6cb",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "root_dir = cfg.DATASET.ROOT\n",
+    "pipeline = cfg.DATASET.PIPELINE  # fmri pre-processing pipeline\n",
+    "atlas = cfg.DATASET.ATLAS\n",
+    "site_ids = cfg.DATASET.SITE_IDS\n",
+    "abide = fetch_abide_pcp(data_dir=root_dir, pipeline=pipeline, \n",
+    "                        band_pass_filtering=True, global_signal_regression=False, \n",
+    "                        derivatives=atlas, quality_checked=False,\n",
+    "                        SITE_ID=site_ids, \n",
+    "                        verbose=0)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "84d813a3",
+   "metadata": {
+    "pycharm": {
+     "name": "#%% md\n"
+    }
+   },
+   "source": [
+    "### Read Phenotypic data"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "eeb46db6",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "pheno_file = os.path.join(cfg.DATASET.ROOT, \"ABIDE_pcp/Phenotypic_V1_0b_preprocessed1.csv\")\n",
+    "pheno_info = pd.read_csv(pheno_file, index_col=0)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "a9df053f",
+   "metadata": {
+    "pycharm": {
+     "name": "#%% md\n"
+    }
+   },
+   "source": [
+    "View Phenotypic data"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "c871fdae",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "pheno_info.head()"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "12a19241",
+   "metadata": {
+    "pycharm": {
+     "name": "#%% md\n"
+    }
+   },
+   "source": [
+    "### Read timeseries from files"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "5f9a9afa",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "data_dir = os.path.join(root_dir, \"ABIDE_pcp/%s/filt_noglobal\" % pipeline)\n",
+    "use_idx = []\n",
+    "time_series = []\n",
+    "for i in pheno_info.index:\n",
+    "    data_file_name = \"%s_%s.1D\" % (pheno_info.loc[i, \"FILE_ID\"], atlas)\n",
+    "    data_path = os.path.join(data_dir, data_file_name)\n",
+    "    if os.path.exists(data_path):\n",
+    "        time_series.append(np.loadtxt(data_path, skiprows=0))\n",
+    "        use_idx.append(i)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "d665b774",
+   "metadata": {
+    "pycharm": {
+     "name": "#%% md\n"
+    }
+   },
+   "source": [
+    "Use \"DX_GROUP\" (autism vs control) as labels, and \"SITE_ID\" as covariates"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "03d09b3c",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "pheno = pheno_info.loc[use_idx, [\"SITE_ID\", \"DX_GROUP\"]].reset_index(drop=True)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "f2f325ec",
+   "metadata": {
+    "pycharm": {
+     "name": "#%% md\n"
+    }
+   },
+   "source": [
+    "## Extracting Brain Networks Features"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "ad4037a9",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "from nilearn.connectome import ConnectivityMeasure\n",
+    "\n",
+    "correlation_measure = ConnectivityMeasure(kind='correlation', vectorize=True)\n",
+    "brain_networks = correlation_measure.fit_transform(time_series)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "9e31d359",
+   "metadata": {
+    "pycharm": {
+     "name": "#%% md\n"
+    }
+   },
+   "source": [
+    "## Machine Learning for Multi-site Data\n",
+    "\n",
+    "### Cross validation Pipeline"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "25df0e1e",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "from sklearn.metrics import accuracy_score\n",
+    "from sklearn.preprocessing import LabelEncoder\n",
+    "import torch\n",
+    "from torch.nn.functional import one_hot\n",
+    "\n",
+    "def cross_validation(x, y, covariates, estimator, domain_adaptation=False, domain_generalization=False):\n",
+    "    results = {\"Target\": [], \"Num_samples\": [], \"Accuracy\": []}\n",
+    "    unique_covariates = np.unique(covariates)\n",
+    "    n_covariates = len(unique_covariates)\n",
+    "    le = LabelEncoder()\n",
+    "    covariate_mat = one_hot(torch.as_tensor(le.fit_transform(covariates)))\n",
+    "    \n",
+    "    for tgt in unique_covariates:\n",
+    "        idx_tgt = np.where(covariates == tgt)\n",
+    "        idx_src = np.where(covariates != tgt)\n",
+    "        x_tgt = brain_networks[idx_tgt]\n",
+    "        x_src = brain_networks[idx_src]\n",
+    "        y_tgt = y[idx_tgt]\n",
+    "        y_src = y[idx_src]        \n",
+    "        \n",
+    "        if domain_generalization:\n",
+    "            estimator.fit(x_src, y_src, covariate_mat[idx_src])\n",
+    "        elif domain_adaptation:\n",
+    "            estimator.fit(np.concatenate((x_src, x_tgt)), y_src, \n",
+    "                          np.concatenate((covariate_mat[idx_src], covariate_mat[idx_tgt])))\n",
+    "        else:            \n",
+    "            estimator.fit(x_src, y_src)\n",
+    "        y_pred = estimator.predict(x_tgt)\n",
+    "        results[\"Accuracy\"].append(accuracy_score(y_tgt, y_pred))\n",
+    "        results[\"Target\"].append(tgt)\n",
+    "        results[\"Num_samples\"].append(x_tgt.shape[0])\n",
+    "    \n",
+    "    mean_acc = sum([results[\"Num_samples\"][i] * results[\"Accuracy\"][i] for i in range(n_covariates)])\n",
+    "    mean_acc /= x.shape[0]\n",
+    "    results[\"Target\"].append(\"Average\")\n",
+    "    results[\"Num_samples\"].append(x.shape[0])\n",
+    "    results[\"Accuracy\"].append(mean_acc)\n",
+    "    \n",
+    "    return pd.DataFrame(results)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "1ec04d89",
+   "metadata": {
+    "pycharm": {
+     "name": "#%% md\n"
+    }
+   },
+   "source": [
+    "### Baseline"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "1fa57354",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "from sklearn.linear_model import RidgeClassifier\n",
+    "\n",
+    "estimator = RidgeClassifier()\n",
+    "res_df = cross_validation(brain_networks, pheno[\"DX_GROUP\"].values, pheno[\"SITE_ID\"], estimator)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "c9c3a33e",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "res_df"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "93fe6a0c",
+   "metadata": {
+    "pycharm": {
+     "name": "#%% md\n"
+    }
+   },
+   "source": [
+    "### Domain Adaptation"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "d1230e54",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "from kale.pipeline.multi_domain_adapter import CoIRLS\n",
+    "estimator = CoIRLS(kernel=cfg.MODEL.KERNEL, lambda_=cfg.MODEL.LAMBDA_, alpha=cfg.MODEL.ALPHA)\n",
+    "res_df = cross_validation(brain_networks, pheno[\"DX_GROUP\"].values, pheno[\"SITE_ID\"], \n",
+    "                          estimator, domain_adaptation=True)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "8d49c39b",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "res_df"
+   ]
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "Python 3 (ipykernel)",
+   "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.7"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}
diff --git a/examples/toy_domain_adaptation/__init__.py b/examples/toy_domain_adaptation/__init__.py
new file mode 100644
index 000000000..e69de29bb
diff --git a/examples/toy_domain_adaptation/main.py b/examples/toy_domain_adaptation/main.py
new file mode 100644
index 000000000..416f72dd9
--- /dev/null
+++ b/examples/toy_domain_adaptation/main.py
@@ -0,0 +1,99 @@
+import matplotlib.pyplot as plt
+import numpy as np
+import seaborn as sns
+from sklearn.datasets import make_blobs
+from sklearn.linear_model import RidgeClassifier
+from sklearn.metrics import accuracy_score
+from sklearn.preprocessing import OneHotEncoder
+
+from kale.interpret.visualize import distplot_1d
+from kale.pipeline.multi_domain_adapter import CoIRLS
+
+
+def main():
+    np.random.seed(29118)
+    # Generate toy data
+    n_samples = 200
+
+    xs, ys = make_blobs(n_samples, centers=[[0, 0], [0, 2]], cluster_std=[0.3, 0.35])
+    xt, yt = make_blobs(n_samples, centers=[[2, -2], [2, 0.2]], cluster_std=[0.35, 0.4])
+
+    # visualize toy data
+    colors = ["c", "m"]
+    x_all = [xs, xt]
+    y_all = [ys, yt]
+    labels = ["source", "Target"]
+    plt.figure(figsize=(8, 5))
+    for i in range(2):
+        idx_pos = np.where(y_all[i] == 1)
+        idx_neg = np.where(y_all[i] == 0)
+        plt.scatter(
+            x_all[i][idx_pos, 0],
+            x_all[i][idx_pos, 1],
+            c=colors[i],
+            marker="o",
+            alpha=0.4,
+            label=labels[i] + " positive",
+        )
+        plt.scatter(
+            x_all[i][idx_neg, 0],
+            x_all[i][idx_neg, 1],
+            c=colors[i],
+            marker="x",
+            alpha=0.4,
+            label=labels[i] + " negative",
+        )
+    plt.legend()
+    plt.title("Source domain and target domain blobs data", fontsize=14, fontweight="bold")
+    plt.show()
+
+    clf = RidgeClassifier(alpha=1.0)
+    clf.fit(xs, ys)
+
+    yt_pred = clf.predict(xt)
+    print("Accuracy on target domain: {:.2f}".format(accuracy_score(yt, yt_pred)))
+
+    # visualize decision scores of non-adaptation classifier
+    ys_score = clf.decision_function(xs)
+    yt_score = clf.decision_function(xt)
+    title = "Ridge classifier decision score distribution"
+    title_kwargs = {"fontsize": 14, "fontweight": "bold"}
+    hist_kwargs = {"kde": True, "alpha": 0.7}
+    plt_labels = ["Source", "Target"]
+    distplot_1d(
+        [ys_score, yt_score],
+        labels=plt_labels,
+        xlabel="Decision Scores",
+        title=title,
+        title_kwargs=title_kwargs,
+        hist_kwargs=hist_kwargs,
+    ).show()
+
+    # domain adaptation
+    clf_ = CoIRLS(lambda_=1)
+    # encoding one-hot domain covariate matrix
+    covariates = np.zeros(n_samples * 2)
+    covariates[:n_samples] = 1
+    enc = OneHotEncoder(handle_unknown="ignore")
+    covariates_mat = enc.fit_transform(covariates.reshape(-1, 1)).toarray()
+
+    x = np.concatenate((xs, xt))
+    clf_.fit(x, ys, covariates_mat)
+    yt_pred_ = clf_.predict(xt)
+    print("Accuracy on target domain: {:.2f}".format(accuracy_score(yt, yt_pred_)))
+
+    ys_score_ = clf_.decision_function(xs).detach().numpy().reshape(-1)
+    yt_score_ = clf_.decision_function(xt).detach().numpy().reshape(-1)
+    title = "Domain adaptation classifier decision score distribution"
+    distplot_1d(
+        [ys_score_, yt_score_],
+        labels=plt_labels,
+        xlabel="Decision Scores",
+        title=title,
+        title_kwargs=title_kwargs,
+        hist_kwargs=hist_kwargs,
+    ).show()
+
+
+if __name__ == "__main__":
+    main()
diff --git a/examples/toy_domain_adaptation/tutorial.ipynb b/examples/toy_domain_adaptation/tutorial.ipynb
new file mode 100644
index 000000000..32a1773ff
--- /dev/null
+++ b/examples/toy_domain_adaptation/tutorial.ipynb
@@ -0,0 +1,320 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "id": "73cb63d5",
+   "metadata": {
+    "pycharm": {
+     "name": "#%% md\n"
+    }
+   },
+   "source": [
+    "# PyKale Tutorial: Domain Adaptation on Toy Data\n",
+    "| [Open in Colab](https://colab.research.google.com/github/sz144/pykale/blob/brain-example/examples/toy_domain_adaptation/tutorial.ipynb) (click `Runtime` → `Run all (Ctrl+F9)` |  [Launch Binder](https://mybinder.org/v2/gh/sz144/pykale/HEAD?filepath=examples%2Ftoy_domain_adaptation%2Ftutorial.ipynb) (click `Run` → `Run All Cells`) |"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "aede3243",
+   "metadata": {
+    "pycharm": {
+     "name": "#%% md\n"
+    }
+   },
+   "source": [
+    "### Setup"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "ebd39507",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "# import seaborn first to avoid seaborn import error caused by newer scipy version, to be solved later\n",
+    "import seaborn as sns"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "771b7124",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "if 'google.colab' in str(get_ipython()):\n",
+    "    print('Running on CoLab')\n",
+    "    !pip uninstall --yes imgaug && pip uninstall --yes albumentations && pip install git+https://github.com/aleju/imgaug.git\n",
+    "    !pip install git+https://github.com/pykale/pykale.git\n",
+    "    !git clone https://github.com/pykale/pykale.git\n",
+    "    %cd pykale/examples/toy_domain_adaptation\n",
+    "else:\n",
+    "    print('Not running on CoLab')"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "57b52723",
+   "metadata": {
+    "pycharm": {
+     "name": "#%% md\n"
+    }
+   },
+   "source": [
+    "### Generate toy data"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "61582582",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "%matplotlib inline\n",
+    "import matplotlib.pyplot as plt\n",
+    "import numpy as np\n",
+    "\n",
+    "from sklearn.datasets import make_moons, make_blobs"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "3d7ed283",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "n_samples = 1000\n",
+    "\n",
+    "xs, ys = make_blobs(n_samples, centers=[[0, 0], [0, 2]], cluster_std=[0.3, 0.35])\n",
+    "xt, yt = make_blobs(n_samples, centers=[[2, -2], [2, 0.2]], cluster_std=[0.35, 0.4])"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "dc371f25",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "colors = [\"c\", \"m\"]\n",
+    "x_all = [xs, xt]\n",
+    "y_all = [ys, yt]\n",
+    "labels = [\"source\", \"Target\"]\n",
+    "for i in range(2):\n",
+    "    idx_pos = np.where(y_all[i] == 1)\n",
+    "    idx_neg = np.where(y_all[i] == 0)\n",
+    "    plt.scatter(x_all[i][idx_pos, 0], x_all[i][idx_pos, 1], c=colors[i], marker=\"o\", alpha=0.4, \n",
+    "                label=labels[i] + \" positive\")\n",
+    "    plt.scatter(x_all[i][idx_neg, 0], x_all[i][idx_neg, 1], c=colors[i], marker=\"x\", alpha=0.4, \n",
+    "                label=labels[i] + \" negative\")\n",
+    "plt.legend()\n",
+    "plt.title('Source domain and target domain blobs data',fontsize=14,fontweight='bold')"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "dea7787b",
+   "metadata": {
+    "pycharm": {
+     "name": "#%% md\n"
+    }
+   },
+   "source": [
+    "### Classification"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "41607742",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "from sklearn.preprocessing import OneHotEncoder\n",
+    "from sklearn.linear_model import RidgeClassifier\n",
+    "from sklearn.metrics import accuracy_score\n",
+    "from kale.interpret.visualize import distplot_1d\n",
+    "from kale.pipeline.multi_domain_adapter import CoIRLS"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "85769302",
+   "metadata": {
+    "pycharm": {
+     "name": "#%% md\n"
+    }
+   },
+   "source": [
+    "#### Training a standard Ridge classifier"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "681c00ac",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "clf = RidgeClassifier(alpha=1.0)\n",
+    "clf.fit(xs, ys)\n",
+    "\n",
+    "yt_pred = clf.predict(xt)\n",
+    "print('Accuracy on target domain: {:.2f}'.format(accuracy_score(yt, yt_pred)))"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "cfc2a82d",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "ys_score = clf.decision_function(xs)\n",
+    "yt_score = clf.decision_function(xt)\n",
+    "\n",
+    "title = \"Ridge classifier decision score distribution\"\n",
+    "title_kwargs = {\"fontsize\": 14, \"fontweight\": \"bold\"}\n",
+    "hist_kwargs = {\"kde\": True, \"alpha\": 0.7}\n",
+    "plt_labels = [\"Source\", \"Target\"]\n",
+    "distplot_1d(\n",
+    "    [ys_score, yt_score],\n",
+    "    title=title,\n",
+    "    xlabel=\"Decision Scores\",\n",
+    "    labels=plt_labels,\n",
+    "    hist_kwargs=hist_kwargs,\n",
+    "    title_kwargs=title_kwargs,\n",
+    ").show()"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "2bb477cc",
+   "metadata": {
+    "pycharm": {
+     "name": "#%% md\n"
+    }
+   },
+   "source": [
+    "#### Training a domain adaptation classifier"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "1bd49414",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "clf_ = CoIRLS()\n",
+    "# encoding one-hot domain covariate matrix\n",
+    "covariates = np.zeros(n_samples * 2)\n",
+    "covariates[:n_samples] = 1\n",
+    "enc = OneHotEncoder(handle_unknown=\"ignore\")\n",
+    "covariates_mat = enc.fit_transform(covariates.reshape(-1, 1)).toarray()"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "4ba2cc5d",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "x = np.concatenate((xs, xt))\n",
+    "clf_.fit(x, ys, covariates_mat)\n",
+    "yt_pred_ = clf_.predict(xt)\n",
+    "print(\"Accuracy on target domain: {:.2f}\".format(accuracy_score(yt, yt_pred_)))"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "4afdbb4c",
+   "metadata": {
+    "pycharm": {
+     "name": "#%%\n"
+    }
+   },
+   "outputs": [],
+   "source": [
+    "ys_score_ = clf_.decision_function(xs).detach().numpy().reshape(-1)\n",
+    "yt_score_ = clf_.decision_function(xt).detach().numpy().reshape(-1)\n",
+    "plt.figure(figsize=(10, 5))\n",
+    "title = \"Domain adaptation classifier decision score distribution\"\n",
+    "distplot_1d(\n",
+    "    [ys_score_, yt_score_],\n",
+    "    title=title,\n",
+    "    xlabel=\"Decision Scores\",\n",
+    "    labels=plt_labels,\n",
+    "    hist_kwargs=hist_kwargs,\n",
+    "    title_kwargs=title_kwargs,\n",
+    ").show()"
+   ]
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "Python 3 (ipykernel)",
+   "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.7"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}
diff --git a/examples/video_loading/README.md b/examples/video_loading/README.md
index 14c9db86d..ed5cd89d3 100644
--- a/examples/video_loading/README.md
+++ b/examples/video_loading/README.md
@@ -49,20 +49,24 @@ for image in frames:
 ```
 ![alt text](https://github.com/RaivoKoot/images/blob/main/Action_Video.jpg "Action Video")
 # Table of Contents
-- [1. Requirements](#1-requirements)
-- [2. Custom Dataset](#2-custom-dataset)
-- [3. Video Frame Sampling Method](#3-video-frame-sampling-method)
-- [4. Alternate Video Frame Sampling Methods](#4-alternate-video-frame-sampling-methods)
-- [5. Using VideoFrameDataset for Training](#5-using-videoframedataset-for-training)
-- [6. Allowing Multiple Labels per Sample](#6-allowing-multiple-labels-per-sample)
-- [7. Conclusion](#7-conclusion)
-- [8. Acknowledgements](#8-acknowledgements)
+- [Efficient Video Dataset Loading, Preprocessing, and Augmentation](#efficient-video-dataset-loading-preprocessing-and-augmentation)
+    - [Overview: This example demonstrates the use of `kale.loaddata.videos.VideoFrameDataset`](#overview-this-example-demonstrates-the-use-of-kaleloaddatavideosvideoframedataset)
+    - [QuickDemo (main.py)](#quickdemo-mainpy)
+- [Table of Contents](#table-of-contents)
+    - [1. Requirements](#1-requirements)
+    - [2. Custom Dataset](#2-custom-dataset)
+    - [3. Video Frame Sampling Method](#3-video-frame-sampling-method)
+    - [4. Alternate Video Frame Sampling Methods](#4-alternate-video-frame-sampling-methods)
+    - [5. Using VideoFrameDataset for training](#5-using-videoframedataset-for-training)
+    - [6. Allowing Multiple Labels per Sample](#6-allowing-multiple-labels-per-sample)
+    - [7. Conclusion](#7-conclusion)
+    - [8. Acknowledgements](#8-acknowledgements)
 
 ### 1. Requirements
 ```
 # Without these three, VideoFrameDataset will not work.
-torchvision >= 0.8.0
-torch >= 1.7.0
+torchvision >= 0.10.0
+torch >= 1.9.0
 python >= 3.6
 ```
 ### 2. Custom Dataset
@@ -172,10 +176,10 @@ list gives label_n for the whole batch.
 A demo of this can be found at the end in `main.py`. It uses the dummy dataset in directory `demo_dataset_multilabel`,
 that is also automatically downloaded.
 
-### 5. Conclusion
+### 7. Conclusion
 A proper code-based explanation on how to use VideoFrameDataset for training is provided in `main.py`
 
-### 6. Acknowledgements
+### 8. Acknowledgements
 We thank the authors of TSN for their [codebase](https://github.com/yjxiong/tsn-pytorch), from which we took VideoFrameDataset and adapted it
 for general use and compatibility.
 ```
diff --git a/kale/embed/mpca.py b/kale/embed/factorization.py
similarity index 60%
rename from kale/embed/mpca.py
rename to kale/embed/factorization.py
index 513388fa6..3e3509a87 100644
--- a/kale/embed/mpca.py
+++ b/kale/embed/factorization.py
@@ -1,21 +1,21 @@
 # =============================================================================
-# Author: Shuo Zhou, szhou20@sheffield.ac.uk
+# Author: Shuo Zhou, shuo.zhou@sheffield.ac.uk
 #         Haiping Lu, h.lu@sheffield.ac.uk or hplu@ieee.org
 # =============================================================================
 
-"""Python implementation of Multilinear Principal Component Analysis (MPCA)
-
-Reference:
-    Haiping Lu, K.N. Plataniotis, and A.N. Venetsanopoulos, "MPCA: Multilinear Principal Component Analysis of Tensor
-    Objects", IEEE Transactions on Neural Networks, Vol. 19, No. 1, Page: 18-39, January 2008. For initial Matlab
-    implementation, please go to https://uk.mathworks.com/matlabcentral/fileexchange/26168.
+"""Python implementation of a tensor factorization algorithm Multilinear Principal Component Analysis (MPCA)
+    and a matrix factorization algorithm Maximum Independence Domain Adaptation (MIDA）
 """
 import logging
 import warnings
 
 import numpy as np
+from numpy.linalg import multi_dot
 from scipy import linalg
 from sklearn.base import BaseEstimator, TransformerMixin
+from sklearn.metrics.pairwise import pairwise_kernels
+from sklearn.preprocessing import KernelCenterer, LabelBinarizer
+from sklearn.utils.validation import check_is_fitted
 
 # import tensorly as tl
 from tensorly.base import fold, unfold
@@ -77,10 +77,16 @@ class MPCA(BaseEstimator, TransformerMixin):
     Attributes:
         proj_mats (list of arrays): A list of transposed projection matrices, shapes (P_1, I_1), ...,
             (P_N, I_N), where P_1, ..., P_N are output tensor shape for each sample.
-        idx_order (array-like): The ordering index of projected (and vectorised) features in decreasing variance.
+        idx_order (array-like): The ordering index of projected (and vectorized) features in decreasing variance.
         mean_ (array-like): Per-feature empirical mean, estimated from the training set, shape (I_1, I_2, ..., I_N).
         shape_in (tuple): Input tensor shapes, i.e. (I_1, I_2, ..., I_N).
         shape_out (tuple): Output tensor shapes, i.e. (P_1, P_2, ..., P_N).
+
+    Reference:
+        Haiping Lu, K.N. Plataniotis, and A.N. Venetsanopoulos, "MPCA: Multilinear Principal Component Analysis of
+        Tensor Objects", IEEE Transactions on Neural Networks, Vol. 19, No. 1, Page: 18-39, January 2008. For initial
+        Matlab implementation, please go to https://uk.mathworks.com/matlabcentral/fileexchange/26168.
+
     Examples:
         >>> import numpy as np
         >>> from kale.embed.mpca import MPCA
@@ -214,8 +220,9 @@ def transform(self, x):
             x_projected = unfold(x_projected, mode=0)
             x_projected = x_projected[:, self.idx_order]
             if isinstance(self.n_components, int):
-                if self.n_components > np.prod(self.shape_out):
-                    self.n_components = np.prod(self.shape_out)
+                n_features = int(np.prod(self.shape_out))
+                if self.n_components > n_features:
+                    self.n_components = n_features
                     warn_msg = "n_components exceeds the maximum number, all features will be returned."
                     logging.warning(warn_msg)
                     warnings.warn(warn_msg)
@@ -228,7 +235,7 @@ def inverse_transform(self, x):
 
         Args:
             x (array-like tensor): Data to be reconstructed, shape (n_samples, P_1, P_2, ..., P_N), if
-                self.vectorize == False, where P_1, P_2, ..., P_N are the reduced dimensions of of corresponding
+                self.vectorize == False, where P_1, P_2, ..., P_N are the reduced dimensions of corresponding
                 mode (1, 2, ..., N), respectively. If self.vectorize == True, shape (n_samples, self.n_components)
                 or shape (n_samples, P_1 * P_2 * ... * P_N).
 
@@ -258,3 +265,149 @@ def inverse_transform(self, x):
         x_rec = x_rec + self.mean_
 
         return x_rec
+
+
+class MIDA(BaseEstimator, TransformerMixin):
+    """Maximum independence domain adaptation
+    Args:
+        n_components (int): Number of components to keep.
+        kernel (str): "linear", "rbf", or "poly". Kernel to use for MIDA. Defaults to "linear".
+        mu (float): Hyperparameter of the l2 penalty. Defaults to 1.0.
+        eta (float): Hyperparameter of the label dependence. Defaults to 1.0.
+        augmentation (bool): Whether using covariates as augment features. Defaults to False.
+        kernel_params (dict or None): Parameters for the kernel. Defaults to None.
+
+    References:
+        [1] Yan, K., Kou, L. and Zhang, D., 2018. Learning domain-invariant subspace using domain features and
+            independence maximization. IEEE transactions on cybernetics, 48(1), pp.288-299.
+    """
+
+    def __init__(
+        self, n_components, kernel="linear", lambda_=1.0, mu=1.0, eta=1.0, augmentation=False, kernel_params=None,
+    ):
+        self.n_components = n_components
+        self.kernel = kernel
+        self.mu = mu
+        self.eta = eta
+        self.augmentation = augmentation
+        if kernel_params is None:
+            self.kernel_params = {}
+        else:
+            self.kernel_params = kernel_params
+        self._label_binarizer = LabelBinarizer(pos_label=1, neg_label=-1)
+        self._centerer = KernelCenterer()
+        self.x_fit = None
+
+    def _get_kernel(self, x, y=None):
+        if self.kernel in ["linear", "rbf", "poly"]:
+            params = self.kernel_params or {}
+        else:
+            raise ValueError("Pre-computed kernel not supported")
+        return pairwise_kernels(x, y, metric=self.kernel, filter_params=True, **params)
+
+    def fit(self, x, y=None, covariates=None):
+        """
+        Args:
+            x : array-like. Input data, shape (n_samples, n_features)
+            y : array-like. Labels, shape (nl_samples,)
+            covariates : array-like. Domain co-variates, shape (n_samples, n_co-variates)
+
+        Note:
+            Unsupervised MIDA is performed if y is None.
+            Semi-supervised MIDA is performed is y is not None.
+        """
+        if self.augmentation and type(covariates) == np.ndarray:
+            x = np.concatenate((x, covariates), axis=1)
+
+        # Kernel matrix
+        kernel_x = self._get_kernel(x)
+        kernel_x[np.isnan(kernel_x)] = 0
+
+        # Solve the optimization problem
+        self._fit(kernel_x, y, covariates)
+        self.x_fit = x
+
+        return self
+
+    def _fit(self, kernel_x, y, covariates=None):
+        """solve MIDA
+
+        Args:
+            kernel_x: array-like, kernel matrix of input data x, shape (n_samples, n_samples)
+            y: array-like. Labels, shape (nl_samples,)
+            covariates: array-like. Domain co-variates, shape (n_samples, n_covariates)
+
+        Returns:
+            self
+        """
+        n_samples = kernel_x.shape[0]
+        # Identity (unit) matrix
+        unit_mat = np.eye(n_samples)
+        # Centering matrix
+        ctr_mat = unit_mat - 1.0 / n_samples * np.ones((n_samples, n_samples))
+
+        kernel_x = self._centerer.fit_transform(kernel_x)
+        if type(covariates) == np.ndarray:
+            kernel_c = np.dot(covariates, covariates.T)
+        else:
+            kernel_c = np.zeros((n_samples, n_samples))
+
+        if y is not None:
+            n_labeled = y.shape[0]
+            if n_labeled > n_samples:
+                raise ValueError("Number of labels exceeds number of samples")
+            y_mat_ = self._label_binarizer.fit_transform(y)
+            y_mat = np.zeros((n_samples, y_mat_.shape[1]))
+            y_mat[:n_labeled, :] = y_mat_
+            ker_y = np.dot(y_mat, y_mat.T)
+            obj = multi_dot(
+                [
+                    kernel_x,
+                    self.mu * ctr_mat
+                    + self.eta * multi_dot([ctr_mat, ker_y, ctr_mat])
+                    - multi_dot([ctr_mat, kernel_c, ctr_mat]),
+                    kernel_x.T,
+                ]
+            )
+        else:
+            obj = multi_dot([kernel_x, self.mu * ctr_mat - multi_dot([ctr_mat, kernel_c, ctr_mat]), kernel_x.T])
+
+        eig_values, eig_vectors = linalg.eigh(obj, subset_by_index=[n_samples - self.n_components, n_samples - 1])
+        idx_sorted = eig_values.argsort()[::-1]
+
+        self.eig_values_ = eig_values[idx_sorted]
+        self.U = eig_vectors[:, idx_sorted]
+        self.U = np.asarray(self.U, dtype=np.float)
+
+        return self
+
+    def fit_transform(self, x, y=None, covariates=None):
+        """
+        Args:
+            x : array-like, shape (n_samples, n_features)
+            y : array-like, shape (n_samples,)
+            covariates : array-like, shape (n_samples, n_covariates)
+
+        Returns:
+            x_transformed : array-like, shape (n_samples, n_components)
+        """
+        self.fit(x, y, covariates)
+
+        return self.transform(x, covariates)
+
+    def transform(self, x, covariates=None):
+        """
+        Args:
+            x : array-like, shape (n_samples, n_features)
+            covariates : array-like, augmentation features, shape (n_samples, n_covariates)
+        Returns:
+            x_transformed : array-like, shape (n_samples, n_components)
+        """
+        check_is_fitted(self, "x_fit")
+        if type(covariates) == np.ndarray and self.augmentation:
+            x = np.concatenate((x, covariates), axis=1)
+        kernel_x = self._centerer.transform(
+            pairwise_kernels(x, self.x_fit, metric=self.kernel, filter_params=True, **self.kernel_params)
+        )
+
+        return np.dot(kernel_x, self.U)
diff --git a/kale/interpret/model_weights.py b/kale/interpret/model_weights.py
index 906f53d7d..c26b18be1 100644
--- a/kale/interpret/model_weights.py
+++ b/kale/interpret/model_weights.py
@@ -1,5 +1,5 @@
 # =============================================================================
-# Author: Shuo Zhou, szhou20@sheffield.ac.uk
+# Author: Shuo Zhou, shuo.zhou@sheffield.ac.uk
 #         Haiping Lu, h.lu@sheffield.ac.uk or hplu@ieee.org
 # =============================================================================
 
diff --git a/kale/interpret/visualize.py b/kale/interpret/visualize.py
index b6e0c855e..8dc638032 100644
--- a/kale/interpret/visualize.py
+++ b/kale/interpret/visualize.py
@@ -1,10 +1,11 @@
 # =============================================================================
-# Author: Shuo Zhou, szhou20@sheffield.ac.uk
+# Author: Shuo Zhou, shuo.zhou@sheffield.ac.uk
 #         Haiping Lu, h.lu@sheffield.ac.uk or hplu@ieee.org
 # =============================================================================
 
 import matplotlib.pyplot as plt
 import numpy as np
+import seaborn as sns
 
 
 def _none2dict(kwarg):
@@ -123,3 +124,57 @@ def plot_multi_images(
         plt.title(image_titles[i])
 
     return fig
+
+
+def distplot_1d(
+    data,
+    labels=None,
+    xlabel=None,
+    ylabel=None,
+    title=None,
+    figsize=None,
+    colors=None,
+    title_kwargs=None,
+    hist_kwargs=None,
+):
+    """Plot distribution of 1D data.
+
+    Args:
+        data (array-like or list): Data to plot.
+        labels (list, optional): List of labels for each data. Defaults to None.
+        xlabel (str, optional): Label for x-axis. Defaults to None.
+        ylabel (str, optional): Label for y-axis. Defaults to None.
+        title (str, optional): Title of the plot. Defaults to None.
+        figsize (tuple, optional): Figure size. Defaults to None.
+        colors (str, optional): Color of the line. Defaults to None.
+        title_kwargs (dict, optional): Key word arguments for title. Defaults to None.
+        hist_kwargs (dict, optional): Key word arguments for histogram. Defaults to None.
+
+    Returns:
+        [matplotlib.figure.Figure]: Figure to plot.
+    """
+    hist_kwargs = _none2dict(hist_kwargs)
+    title_kwargs = _none2dict(title_kwargs)
+
+    fig = plt.figure(figsize=figsize)
+    if colors is None:
+        colors = plt.get_cmap("Set1").colors
+
+    if type(data) != list:
+        data = [data]
+
+    if labels is None:
+        labels = np.range(len(data))
+
+    for i in range(len(data)):
+        sns.histplot(data[i], color=colors[i], label=labels[i], **hist_kwargs)
+
+    if title is not None:
+        plt.title(title, **title_kwargs)
+    if xlabel is not None:
+        plt.xlabel(xlabel)
+    if ylabel is not None:
+        plt.ylabel(ylabel)
+    plt.legend()
+
+    return fig
diff --git a/kale/pipeline/mpca_trainer.py b/kale/pipeline/mpca_trainer.py
index 8499aeeec..7913ecf84 100644
--- a/kale/pipeline/mpca_trainer.py
+++ b/kale/pipeline/mpca_trainer.py
@@ -1,5 +1,5 @@
 # =============================================================================
-# Author: Shuo Zhou, szhou20@sheffield.ac.uk
+# Author: Shuo Zhou, shuo.zhou@sheffield.ac.uk
 #         Haiping Lu, h.lu@sheffield.ac.uk or hplu@ieee.org
 # =============================================================================
 
@@ -26,7 +26,7 @@
 from sklearn.svm import LinearSVC, SVC
 from sklearn.utils.validation import check_is_fitted
 
-from ..embed.mpca import MPCA
+from ..embed.factorization import MPCA
 
 param_c_grids = list(np.logspace(-4, 2, 7))
 classifiers = {
diff --git a/kale/pipeline/multi_domain_adapter.py b/kale/pipeline/multi_domain_adapter.py
index 937c406f7..c9bfb8b14 100644
--- a/kale/pipeline/multi_domain_adapter.py
+++ b/kale/pipeline/multi_domain_adapter.py
@@ -1,11 +1,15 @@
 # =============================================================================
-# Author: Shuo Zhou, szhou20@sheffield.ac.uk/sz144@outlook.com
+# Author: Shuo Zhou, shuo.zhou@sheffield.ac.uk/sz144@outlook.com
 # =============================================================================
 """Multi-source domain adaptation pipelines
 """
 
 import torch
 import torch.nn as nn
+from sklearn.base import BaseEstimator, ClassifierMixin
+from sklearn.metrics.pairwise import pairwise_kernels
+from sklearn.preprocessing import LabelBinarizer
+from torch.linalg import multi_dot
 from torch.nn.functional import one_hot
 
 import kale.predict.losses as losses
@@ -333,3 +337,118 @@ def cls_discrepancy(self, x):
                 cls_disc += torch.mean(torch.abs(cls_disc_))
 
         return cls_disc * 2 / (n_domains * (n_domains - 1))
+
+
+class CoIRLS(BaseEstimator, ClassifierMixin):
+    """Covariate-Independence Regularized Least Squares (CoIRLS)
+
+    Args:
+        kernel (str, optional): {"linear", "rbf", "poly"}. Kernel to use. Defaults to "linear".
+        kernel_kwargs (dict or None, optional): Hyperparameter for the kernel. Defaults to None.
+        alpha (float, optional): Hyperparameter of the l2 (Ridge) penalty. Defaults to 1.0.
+        lambda_ (float, optional): Hyperparameter of the covariate dependence.  Defaults to 1.0.
+
+    Reference:
+        [1] Zhou, S., 2022. Interpretable Domain-Aware Learning for Neuroimage Classification (Doctoral dissertation,
+            University of Sheffield).
+        [2] Zhou, S., Li, W., Cox, C.R., & Lu, H. (2020). Side Information Dependence as a Regularizer for Analyzing
+            Human Brain Conditions across Cognitive Experiments. AAAI 2020, New York, USA.
+    """
+
+    def __init__(self, kernel="linear", kernel_kwargs=None, alpha=1.0, lambda_=1.0):
+        super().__init__()
+        self.kernel = kernel
+        self.model = None
+        self.alpha = alpha
+        self.lambda_ = lambda_
+        if kernel_kwargs is None:
+            self.kernel_kwargs = dict()
+        else:
+            self.kernel_kwargs = kernel_kwargs
+        self.device = torch.device("cuda:0" if torch.cuda.is_available() else "cpu")
+        self.losses = {"ovr": [], "pred": [], "code": [], "reg": []}
+        self.x = None
+        self._label_binarizer = LabelBinarizer(pos_label=1, neg_label=-1)
+        self.coef_ = None
+
+    def fit(self, x, y, covariates):
+        """fit a model with input data, labels and covariates
+
+        Args:
+            x (np.ndarray or tensor): shape (n_samples, n_features)
+            y (np.ndarray or tensor): shape (n_samples, )
+            covariates (np.ndarray or tensor): (n_samples, n_covariates)
+        """
+        self._label_binarizer.fit(y)
+        x = torch.as_tensor(x)
+        x = torch.cat([x, torch.ones(x.shape[0], 1)], 1)
+        y = torch.as_tensor(y)
+        krnl_x = torch.as_tensor(
+            pairwise_kernels(x.detach().numpy(), metric=self.kernel, filter_params=True, **self.kernel_kwargs),
+            dtype=torch.float,
+        )
+        n_samples = x.shape[0]
+        n_classes = torch.unique(y).shape[0]
+        n_labeled = y.shape[0]
+        unit_mat = torch.eye(n_samples)
+        ctr_mat = unit_mat - 1.0 / n_samples * torch.ones((n_samples, n_samples))
+        mat_j = torch.zeros((n_samples, n_samples))
+        mat_j[:n_labeled, :n_labeled] = torch.eye(n_labeled)
+
+        covariates = torch.as_tensor(covariates, dtype=torch.float)
+        krnl_cov = torch.mm(covariates, covariates.T)
+
+        if n_classes == 2:
+            mat_y = torch.zeros((n_samples, 1))
+        else:
+            mat_y = torch.zeros((n_samples, n_classes))
+        mat_y[:n_labeled, :] = torch.as_tensor(self._label_binarizer.fit_transform(y))
+        mat_y = torch.as_tensor(mat_y)
+
+        mat_q = torch.mm(mat_j, krnl_x) + self.alpha * unit_mat
+        mat_q += self.lambda_ * multi_dot((ctr_mat, krnl_cov, ctr_mat, krnl_x))
+
+        self.coef_ = torch.linalg.solve(mat_q, mat_y)
+
+        self.x = x
+
+    def predict(self, x):
+        """Predict labels for data x
+
+        Args:
+            x (np.ndarray or tensor): Samples need prediction, shape (n_samples, n_features)
+
+        Returns:
+            y (np.ndarray): Predicted labels, shape (n_samples, )
+        """
+        out = self.decision_function(x)
+        if self._label_binarizer.y_type_ == "binary":
+            pred = self._label_binarizer.inverse_transform(torch.sign(out).view(-1))
+        else:
+            pred = self._label_binarizer.inverse_transform(out)
+
+        return pred
+
+    def decision_function(self, x):
+        """Compute decision scores for data x
+
+        Args:
+            x (np.ndarray or tensor): Samples need decision scores, shape (n_samples, n_features)
+
+        Returns:
+            scores (np.ndarray): Decision scores, shape (n_samples, )
+        """
+        x = torch.as_tensor(x)
+        x = torch.cat([x, torch.ones(x.shape[0], 1)], dim=1)
+        krnl_x = torch.as_tensor(
+            pairwise_kernels(
+                x.detach().numpy(),
+                self.x.detach().numpy(),
+                metric=self.kernel,
+                filter_params=True,
+                **self.kernel_kwargs,
+            ),
+            dtype=torch.float,
+        )
+
+        return torch.mm(krnl_x, self.coef_)
diff --git a/setup.py b/setup.py
index 55a174f3a..1cc13bb4f 100644
--- a/setup.py
+++ b/setup.py
@@ -15,6 +15,7 @@
     "ipykernel",
     "ipython",
     "matplotlib",
+    "nilearn",
     "numpy>=1.18.0",
     "pandas",
     "Pillow",
@@ -27,10 +28,11 @@
     "scikit-image>=0.16.2",
     "scikit-learn>=0.23.2",
     "scipy>=1.5.4",
+    "seaborn",
     "tensorly>=0.5.1",
-    "torch>=1.7.0",
+    "torch>=1.10.0",
     "torchsummary>=1.5.0",
-    "torchvision",  # >=0.8.1
+    "torchvision>=0.11.0",
     "yacs>=0.1.7",
 ]
 


Test patch:
diff --git a/.github/workflows/test.yml b/.github/workflows/test.yml
index 722f5a2ee..a4ebc1786 100644
--- a/.github/workflows/test.yml
+++ b/.github/workflows/test.yml
@@ -36,13 +36,13 @@ jobs:
         run: |
           conda install -c conda-forge rdkit=2020.09.1.0
           pip install pyparsing==2.4.7
-          pip install torch==1.7.1+cpu -f https://download.pytorch.org/whl/torch_stable.html
-          pip install torchvision==0.8.2+cpu -f https://download.pytorch.org/whl/torch_stable.html
-          pip install torch-scatter==2.0.7 -f https://pytorch-geometric.com/whl/torch-1.7.0+cpu.html
-          pip install torch-sparse==0.6.12 -f https://pytorch-geometric.com/whl/torch-1.7.0+cpu.html
-          pip install torch-cluster==1.5.9 -f https://pytorch-geometric.com/whl/torch-1.7.0+cpu.html
-          pip install torch-spline-conv==1.2.1 -f https://pytorch-geometric.com/whl/torch-1.7.0+cpu.html
-          pip install torch-geometric==1.7.0 --ignore-installed certifi
+          pip install torch==1.10.0+cpu -f https://download.pytorch.org/whl/torch_stable.html
+          pip install torchvision==0.11.0+cpu -f https://download.pytorch.org/whl/torch_stable.html
+          pip install torch-scatter==2.0.9 -f https://pytorch-geometric.com/whl/torch-1.10.0+cpu.html
+          pip install torch-sparse==0.6.12 -f https://pytorch-geometric.com/whl/torch-1.10.0+cpu.html
+          pip install torch-cluster==1.5.9 -f https://pytorch-geometric.com/whl/torch-1.10.0+cpu.html
+          pip install torch-spline-conv==1.2.1 -f https://pytorch-geometric.com/whl/torch-1.10.0+cpu.html
+          pip install torch-geometric
           pip install -e .[dev]
         shell: bash -l {0}
       - name: Run tests
diff --git a/tests/embed/test_mpca.py b/tests/embed/test_factorization.py
similarity index 73%
rename from tests/embed/test_mpca.py
rename to tests/embed/test_factorization.py
index 3c8b92263..e9785f4f3 100644
--- a/tests/embed/test_mpca.py
+++ b/tests/embed/test_factorization.py
@@ -4,9 +4,11 @@
 import pytest
 from numpy import testing
 from scipy.io import loadmat
+from sklearn.datasets import make_blobs
+from sklearn.preprocessing import OneHotEncoder
 from tensorly.tenalg import multi_mode_dot
 
-from kale.embed.mpca import MPCA
+from kale.embed.factorization import MIDA, MPCA
 from kale.utils.download import download_file_by_url
 
 N_COMPS = [1, 50, 100]
@@ -80,3 +82,27 @@ def test_mpca_against_baseline(gait, baseline_model):
         # check whether each eigen-vector column is equal to/opposite of corresponding baseline eigen-vector column
         # testing.assert_allclose(abs(mpca.proj_mats[i]), abs(baseline_proj_mats[i]))
         testing.assert_allclose(mpca.proj_mats[i] ** 2, baseline_proj_mats[i] ** 2, rtol=relative_tol)
+
+
+@pytest.mark.parametrize("kernel", ["linear", "rbf"])
+@pytest.mark.parametrize("augmentation", [True, False])
+def test_mida(kernel, augmentation):
+    np.random.seed(29118)
+    # Generate toy data
+    n_samples = 200
+
+    xs, ys = make_blobs(n_samples, n_features=3, centers=[[0, 0, 0], [0, 2, 1]], cluster_std=[0.3, 0.35])
+    xt, yt = make_blobs(n_samples, n_features=3, centers=[[2, -2, 2], [2, 0.2, -1]], cluster_std=[0.35, 0.4])
+    x = np.concatenate((xs, xt), axis=0)
+
+    covariates = np.zeros(n_samples * 2)
+    covariates[:n_samples] = 1
+
+    enc = OneHotEncoder(handle_unknown="ignore")
+    covariates_mat = enc.fit_transform(covariates.reshape(-1, 1)).toarray()
+    mida = MIDA(n_components=2, kernel=kernel, augmentation=augmentation)
+    x_transformed = mida.fit_transform(x, covariates=covariates_mat)
+    testing.assert_equal(x_transformed.shape, (n_samples * 2, 2))
+
+    x_transformed = mida.fit_transform(x, ys, covariates=covariates_mat)
+    testing.assert_equal(x_transformed.shape, (n_samples * 2, 2))
diff --git a/tests/embed/test_gripnet.py b/tests/embed/test_gripnet.py
index 6dbb3d33e..23a70ea1b 100644
--- a/tests/embed/test_gripnet.py
+++ b/tests/embed/test_gripnet.py
@@ -6,7 +6,7 @@
 from kale.embed.gripnet import TypicalGripNetEncoder
 from kale.utils.download import download_file_by_url
 
-pose_url = "https://github.com/pykale/data/raw/main/graphs/pose.pt"
+pose_url = "https://github.com/pykale/data/raw/main/graphs/pose_pyg_2.pt"
 
 
 @pytest.fixture(scope="module")
diff --git a/tests/pipeline/test_deep_dti.py b/tests/pipeline/test_deep_dti.py
index 83a0d648e..869f3a87e 100644
--- a/tests/pipeline/test_deep_dti.py
+++ b/tests/pipeline/test_deep_dti.py
@@ -29,8 +29,8 @@ def test_deep_data(download_path):
     save_parameters = {"seed": 2020, "batch_size": 256}
     model = DeepDTATrainer(drug_encoder, target_encoder, decoder, lr=0.001, ci_metric=True, **save_parameters).eval()
     trainer = pl.Trainer(max_epochs=1, gpus=0)
-    trainer.fit(model, train_dataloader=train_dataloader, val_dataloaders=valid_dataloader)
-    trainer.test(test_dataloaders=test_dataloader)
+    trainer.fit(model, train_dataloaders=train_dataloader, val_dataloaders=valid_dataloader)
+    trainer.test(dataloaders=test_dataloader)
     assert isinstance(model.drug_encoder, CNNEncoder)
     assert isinstance(model.target_encoder, CNNEncoder)
     assert isinstance(model.decoder, MLPDecoder)
diff --git a/tests/pipeline/test_multi_domain_adapter.py b/tests/pipeline/test_multi_domain_adapter.py
index a317b7dfd..12e6017e5 100644
--- a/tests/pipeline/test_multi_domain_adapter.py
+++ b/tests/pipeline/test_multi_domain_adapter.py
@@ -1,10 +1,14 @@
+import matplotlib.figure
 import pytest
 import torch
+from sklearn.metrics import accuracy_score
+from torch.nn.functional import one_hot
 
 from kale.embed.image_cnn import ResNet18Feature
+from kale.interpret.visualize import distplot_1d
 from kale.loaddata.image_access import ImageAccess
 from kale.loaddata.multi_domain import MultiDomainAdapDataset
-from kale.pipeline.multi_domain_adapter import create_ms_adapt_trainer
+from kale.pipeline.multi_domain_adapter import CoIRLS, create_ms_adapt_trainer
 from kale.predict.class_domain_nets import ClassNetSmallImage
 from tests.helpers.pipe_test_helper import ModelTestHelper
 
@@ -60,3 +64,45 @@ def test_multi_source(method, input_dimension, office_caltech_access, testing_cf
     )
     kwargs = {"limit_train_batches": 0.1, "limit_val_batches": 0.3, "limit_test_batches": 0.2}
     ModelTestHelper.test_model(model, train_params, **kwargs)
+
+
+@pytest.mark.parametrize("kernel", ["linear", "rbf"])
+def test_coirls(kernel, office_caltech_access):
+    dataset = MultiDomainAdapDataset(office_caltech_access)
+    dataset.prepare_data_loaders()
+    dataloader = dataset.get_domain_loaders(split="train", batch_size=100)
+    feature_network = ResNet18Feature()
+    x, y, z = next(iter(dataloader))
+    tgt_idx = torch.where(z == 0)
+    src_idx = torch.where(z != 0)
+
+    x_feat = feature_network(x)
+    z_one_hot = one_hot(z)
+    clf = CoIRLS(kernel=kernel, alpha=1.0)
+
+    x_train = torch.cat((x_feat[src_idx], x_feat[tgt_idx]))
+    y_train = y[src_idx]
+    z_train = torch.cat((z_one_hot[src_idx], z_one_hot[tgt_idx]))
+    clf.fit(x_train, y_train, z_train)
+    y_pred = clf.predict(x_feat[tgt_idx])
+
+    acc = accuracy_score(y[tgt_idx], y_pred)
+
+    assert 0 <= acc <= 1
+
+    if kernel == "linear":
+        scores = [clf.decision_function(x_feat[src_idx])[:, 0], clf.decision_function(x_feat[tgt_idx])[:, 0]]
+        domain_labels = ["source", "target"]
+        title = "Decision score distribution"
+        title_kwargs = {"fontsize": 14, "fontweight": "bold"}
+        hist_kwargs = {"alpha": 0.7}
+        fig = distplot_1d(
+            scores,
+            labels=domain_labels,
+            xlabel="Decision Scores",
+            ylabel="Count",
+            title=title,
+            title_kwargs=title_kwargs,
+            hist_kwargs=hist_kwargs,
+        )
+        assert type(fig) == matplotlib.figure.Figure