diff --git a/.Rbuildignore b/.Rbuildignore
index a71d2d5..f586e90 100644
--- a/.Rbuildignore
+++ b/.Rbuildignore
@@ -34,3 +34,7 @@ build-pkg*
 .valgrindrc
 ^doc$
 ^Meta$
+^_pkgdown\.yml$
+^docs$
+^pkgdown$
+^\.github$
diff --git a/.github/.gitignore b/.github/.gitignore
new file mode 100644
index 0000000..2d19fc7
--- /dev/null
+++ b/.github/.gitignore
@@ -0,0 +1 @@
+*.html
diff --git a/.github/workflows/pkgdown.yaml b/.github/workflows/pkgdown.yaml
new file mode 100644
index 0000000..ed7650c
--- /dev/null
+++ b/.github/workflows/pkgdown.yaml
@@ -0,0 +1,48 @@
+# Workflow derived from https://github.com/r-lib/actions/tree/v2/examples
+# Need help debugging build failures? Start at https://github.com/r-lib/actions#where-to-find-help
+on:
+  push:
+    branches: [main, master]
+  pull_request:
+    branches: [main, master]
+  release:
+    types: [published]
+  workflow_dispatch:
+
+name: pkgdown
+
+jobs:
+  pkgdown:
+    runs-on: ubuntu-latest
+    # Only restrict concurrency for non-PR jobs
+    concurrency:
+      group: pkgdown-${{ github.event_name != 'pull_request' || github.run_id }}
+    env:
+      GITHUB_PAT: ${{ secrets.GITHUB_TOKEN }}
+    permissions:
+      contents: write
+    steps:
+      - uses: actions/checkout@v3
+
+      - uses: r-lib/actions/setup-pandoc@v2
+
+      - uses: r-lib/actions/setup-r@v2
+        with:
+          use-public-rspm: true
+
+      - uses: r-lib/actions/setup-r-dependencies@v2
+        with:
+          extra-packages: any::pkgdown, local::.
+          needs: website
+
+      - name: Build site
+        run: pkgdown::build_site_github_pages(new_process = FALSE, install = FALSE)
+        shell: Rscript {0}
+
+      - name: Deploy to GitHub pages 🚀
+        if: github.event_name != 'pull_request'
+        uses: JamesIves/github-pages-deploy-action@v4.4.1
+        with:
+          clean: false
+          branch: gh-pages
+          folder: docs
diff --git a/.gitignore b/.gitignore
index 3b63ec0..c95bebc 100644
--- a/.gitignore
+++ b/.gitignore
@@ -73,3 +73,4 @@ src/*.dll
 inst/doc
 /doc/
 /Meta/
+docs
diff --git a/DESCRIPTION b/DESCRIPTION
index a800cf9..e2c167e 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -3,9 +3,10 @@ Type: Package
 Title: Network Diffusion Algorithms
 Version: 0.2.1
 Date: 2018-04-20
-Authors@R: person("Simon", "Dirmeier",
+Authors@R: c(person("Simon", "Dirmeier",
                   email = "simon.dirmeier@gmx.de",
-                  role = c("aut", "cre"))
+                  role = c("aut")),
+             person("Yinchun", "Su", email = "me@yinchun.su", role = c("ctb", "cre")))
 Maintainer: Simon Dirmeier <simon.dirmeier@gmx.de>
 Description: Implementation of network diffusion algorithms such as
   heat diffusion or Markov random walks. Network diffusion algorithms generally
@@ -20,6 +21,7 @@ License: GPL (>=3)
 Depends: R (>= 3.4)
 LazyData: TRUE
 LinkingTo: Rcpp, RcppEigen
+biocViews: sparseMatrixStats
 Imports:
     Rcpp,
     igraph,
@@ -34,7 +36,7 @@ Suggests:
     lintr,
     Matrix
 VignetteBuilder: knitr
-RoxygenNote: 7.2.3
+RoxygenNote: 7.3.0
 Encoding: UTF-8
 Config/testthat/edition: 3
 SystemRequirements: C++17
diff --git a/R/diffusr-package.R b/R/diffusr-package.R
index 8c9d501..eeb69b5 100644
--- a/R/diffusr-package.R
+++ b/R/diffusr-package.R
@@ -25,9 +25,6 @@
 #' @author Simon Dirmeier \email{simon.dirmeier@@gmx.de}
 #' @name diffusr-package
 #'
-#' @docType package
-#' @keywords package
-#'
 #' @references
 #'
 #' Tong, H., Faloutsos, C., & Pan, J. Y. (2006),
diff --git a/R/heat_diffusion.R b/R/heat_diffusion.R
index fd3c571..e31c13b 100644
--- a/R/heat_diffusion.R
+++ b/R/heat_diffusion.R
@@ -21,15 +21,17 @@
 #' Graph diffusion using a heat diffusion process on a Laplacian matrix.
 #'
 #' @description An amount of starting heat gets distribution using the
-#' Laplacian matrix of a graph. Every iteration (or time interval) \code{t}
-#'heat streams from the starting nodes into surrounding nodes.
+#'   Laplacian matrix of a graph. Every iteration (or time interval) \eqn{t}
+#'   heat streams from the starting nodes into surrounding nodes.
 #'
 #' @export
 #'
-#' @param h0   an \code{n x p}-dimensional numeric non-negative vector/matrix
-#'  of starting temperatures
-#' @param graph  an (\code{n x n})-dimensional numeric non-negative adjacence
-#'  matrix representing the graph
+#' @param h0   an \eqn{n \times p}-dimensional numeric non-negative
+#'   \code{\link[base]{matrix}} (or
+#'   \code{\link[Matrix:dgCMatrix-class]{dgCMatrix}},
+#'   \code{\link[base]{vector}}) of starting temperatures
+#' @param graph  an (\eqn{n \times n})-dimensional numeric non-negative
+#'   adjacence matrix representing the graph
 #' @param t  time point when heat is measured
 #' @param ...  additional parameters
 #' @return  returns the heat on every node as numeric vector
@@ -41,7 +43,9 @@
 #' @importFrom Rcpp sourceCpp
 #'
 #' @references
-#' \url{https://en.wikipedia.org/wiki/Laplacian_matrix} \cr
+#'
+#' \url{https://en.wikipedia.org/wiki/Laplacian_matrix}
+#'
 #' \url{https://en.wikipedia.org/wiki/Heat_equation}
 #'
 #' @examples
diff --git a/R/is.dgCMatrix.R b/R/is.dgCMatrix.R
index 859b15d..d9e6a44 100644
--- a/R/is.dgCMatrix.R
+++ b/R/is.dgCMatrix.R
@@ -1,4 +1,4 @@
-#' Check if a matrix is `dgCMatrix` class
+#' Check if a matrix is \code{dgCMatrix} class
 #'
 #' @param mat The matrix to check
 #'
diff --git a/R/mat_util.R b/R/mat_util.R
index 92191f7..54aff7f 100644
--- a/R/mat_util.R
+++ b/R/mat_util.R
@@ -23,11 +23,10 @@
 #' @export
 #'
 #' @param obj  \code{\link[base]{matrix}} (or
-#'   \code{\link[Matrix:dgCMatrix-class]{dgCMatrix}}, \link[base]{vector}) that
-#'   is stochstically normalized
+#'   \code{\link[Matrix:dgCMatrix-class]{dgCMatrix}},
+#'   \code{\link[base]{vector}}) that is stochstically normalized
 #' @param ...  additional params
-#' @return  returns the normalized \code{\link[base]{matrix}} (or
-#'   \code{\link[Matrix:dgCMatrix-class]{dgCMatrix}}, \link[base]{vector})
+#' @return  returns the normalized matrix/vector)
 #'
 #' @importFrom checkmate assert check_matrix test_numeric test_atomic_vector
 #'
diff --git a/R/mrw.R b/R/mrw.R
index 92eb7a1..d539716 100644
--- a/R/mrw.R
+++ b/R/mrw.R
@@ -21,21 +21,21 @@
 #' Graph diffusion using a Markov random walk
 #'
 #' @description
-#' A Markov Random Walk takes an inital distribution \code{p0}
+#' A Markov Random Walk takes an inital distribution \eqn{p_0}
 #' and calculates the stationary distribution of that.
-#' The diffusion process is regulated by a restart probability \code{r} which
+#' The diffusion process is regulated by a restart probability \eqn{r} which
 #' controls how often the MRW jumps back to the initial values.
 #'
-#' @param p0  an \code{n x p}-dimensional numeric non-negative vector/matrix
+#' @param p0  an \eqn{n \times p}-dimensional numeric non-negative vector/matrix
 #'  representing the starting distribution of the Markov chain
 #'  (does not need to sum to one).
 #'
-#' @param graph  an (\code{n x n})-dimensional numeric non-negative adjacence
+#' @param graph  an \eqn{n \times p}-dimensional numeric non-negative adjacence
 #'   \code{\link[base]{matrix}} (or
 #'   \code{\link[Matrix:dgCMatrix-class]{dgCMatrix}}) representing the graph
 #'
-#' @param r  a scalar between (0, 1). restart probability if a Markov random
-#' walk with restart is desired
+#' @param r  a scalar between \eqn{(0, 1)}. restart probability if a Markov
+#'   random walk with restart is desired
 #'
 #' @param thresh  threshold for breaking the iterative computation of the
 #'  stationary distribution. If the absolute difference of the distribution at
@@ -53,7 +53,8 @@
 #' @param correct.for.hubs if \code{TRUE} multiplies a correction factor to the
 #'  nodes, such that the random walk gets not biased to nodes with high
 #'  degree. In that case the original input matrix will be normalized as:
-#'  \deqn{ P(j | i) =  1 /degree(i) *  min(1, degree(j)/degree(j))}
+#'  \deqn{ P(j | i) = \dfrac{1}{\text{degree}(i)} \times
+#'    \min \left(1, \dfrac{\text{degree}(j)}{\text{degree}(i)}\right)}
 #'  \emph{Note that this will not consider edge weights.}
 #'
 #' @param allow.ergodic Allow multiple components in a graph.
@@ -62,9 +63,9 @@
 #'
 #' @return  returns a list with the following elements
 #'  \itemize{
-#'   \item p.inf  the stationary distribution as numeric vector
-#'   \item transition.matrix the column normalized transition matrix used for
-#'         the random walk
+#'   \item \code{p.inf}  the stationary distribution as numeric vector
+#'   \item \code{transition.matrix} the column normalized transition matrix used
+#'         for the random walk
 #'  }
 #'
 #' @references
diff --git a/R/neighbors.R b/R/neighbors.R
index b715cff..4963521 100644
--- a/R/neighbors.R
+++ b/R/neighbors.R
@@ -23,7 +23,7 @@
 #' @description
 #' For every node in a set of nodes the graph gets traversed along the node's
 #' shortest paths to its neighbors. Nearest neighbors are added until a maximum
-#' depth of \code{k} is reached. For settings where there are more than \code{k}
+#' depth of \eqn{k} is reached. For settings where there are more than \eqn{k}
 #' neighbors having the same distance, all neighbors are returned.
 #'
 #' @export
@@ -33,8 +33,8 @@
 #'
 #' @param graph an (\eqn{n \times n})-dimensional numeric non-negative adjacence
 #'   \code{\link[base]{matrix}} (or
-#'   \code{\link[Matrix:dgCMatrix-class]{dgCMatrix}}, \link[base]{vector})
-#'   representing the graph
+#'   \code{\link[Matrix:dgCMatrix-class]{dgCMatrix}},
+#'   \code{\link[base]{vector}}) representing the graph
 #'
 #' @param k the depth of the nearest neighbor search, e.g. the depth of the
 #'   graph traversal
diff --git a/README.md b/README.md
index d4929ac..0b8905b 100644
--- a/README.md
+++ b/README.md
@@ -20,13 +20,15 @@ Before installation, we recommended you install Intel oneAPI Math Kernel Library
 Windows users can download oneMKL from [Intel's website](https://www.intel.com/content/www/us/en/developer/tools/oneapi/onemkl-download.html) and install it in the default directory. The default directory is: `C:\Program Files (x86)\Intel\oneAPI`.
 
 Debian users can download oneMKL using apt in the Debian non-free repo:
-```{bash}
+
+``` bash
 # Install oneMKL version 2020.4.304-4
 sudo apt install intel-mkl-full
 ```
 
 Or using the Intel repo:
-```{bash}
+
+``` bash
 # Set up the repository and signed the entry
 wget -O- https://apt.repos.intel.com/intel-gpg-keys/GPG-PUB-KEY-INTEL-SW-PRODUCTS.PUB \
 | gpg --dearmor | sudo tee /usr/share/keyrings/oneapi-archive-keyring.gpg > /dev/null
@@ -38,7 +40,8 @@ sudo apt install intel-oneapi-mkl
 ```
 
 Fedora users can download oneMKL by using dnf:
-```{bash}
+
+``` bash
 # Create dnf repository file
 tee > /tmp/oneAPI.repo << EOF
 [oneAPI]
@@ -57,12 +60,14 @@ sudo dnf install intel-oneapi-mkl
 ## Installation
 
 Install `diffusr` using:
-```{r}
+
+``` r
 install.packages("diffusr")
 ```
 
 Alternatively use the latest version from github:
-```{r}
+
+``` r
 devtools::install_github("randef1ned/diffusr", build_vignettes = TRUE)
 ```
 
diff --git a/_pkgdown.yml b/_pkgdown.yml
new file mode 100644
index 0000000..a6167be
--- /dev/null
+++ b/_pkgdown.yml
@@ -0,0 +1,4 @@
+url: https://diffusr.yinchun.su
+template:
+  bootstrap: 5
+
diff --git a/man/diffusr-package.Rd b/man/diffusr-package.Rd
index 6635e49..e516216 100644
--- a/man/diffusr-package.Rd
+++ b/man/diffusr-package.Rd
@@ -1,6 +1,5 @@
 % Generated by roxygen2: do not edit by hand
 % Please edit documentation in R/diffusr-package.R
-\docType{package}
 \name{diffusr-package}
 \alias{diffusr-package}
 \title{diffusr}
@@ -32,4 +31,3 @@ across pathways and protein complexes.
 \author{
 Simon Dirmeier \email{simon.dirmeier@gmx.de}
 }
-\keyword{package}
diff --git a/man/heat.diffusion.Rd b/man/heat.diffusion.Rd
index c9f86f2..9b9f997 100644
--- a/man/heat.diffusion.Rd
+++ b/man/heat.diffusion.Rd
@@ -7,11 +7,13 @@
 heat.diffusion(h0, graph, t = 0.5, ...)
 }
 \arguments{
-\item{h0}{an \code{n x p}-dimensional numeric non-negative vector/matrix
-of starting temperatures}
+\item{h0}{an \eqn{n \times p}-dimensional numeric non-negative
+\code{\link[base]{matrix}} (or
+\code{\link[Matrix:dgCMatrix-class]{dgCMatrix}},
+\code{\link[base]{vector}}) of starting temperatures}
 
-\item{graph}{an (\code{n x n})-dimensional numeric non-negative adjacence
-matrix representing the graph}
+\item{graph}{an (\eqn{n \times n})-dimensional numeric non-negative
+adjacence matrix representing the graph}
 
 \item{t}{time point when heat is measured}
 
@@ -22,8 +24,8 @@ returns the heat on every node as numeric vector
 }
 \description{
 An amount of starting heat gets distribution using the
-Laplacian matrix of a graph. Every iteration (or time interval) \code{t}
-heat streams from the starting nodes into surrounding nodes.
+  Laplacian matrix of a graph. Every iteration (or time interval) \eqn{t}
+  heat streams from the starting nodes into surrounding nodes.
 }
 \examples{
 # count of nodes
@@ -36,6 +38,7 @@ graph <- matrix(abs(rnorm(n*n)), n, n)
 heat <- heat.diffusion(h0, graph)
 }
 \references{
-\url{https://en.wikipedia.org/wiki/Laplacian_matrix} \cr
+\url{https://en.wikipedia.org/wiki/Laplacian_matrix}
+
 \url{https://en.wikipedia.org/wiki/Heat_equation}
 }
diff --git a/man/is.dgCMatrix.Rd b/man/is.dgCMatrix.Rd
index 5cd188e..dd879de 100644
--- a/man/is.dgCMatrix.Rd
+++ b/man/is.dgCMatrix.Rd
@@ -2,7 +2,7 @@
 % Please edit documentation in R/is.dgCMatrix.R
 \name{is.dgCMatrix}
 \alias{is.dgCMatrix}
-\title{Check if a matrix is `dgCMatrix` class}
+\title{Check if a matrix is \code{dgCMatrix} class}
 \usage{
 is.dgCMatrix(mat)
 }
@@ -13,5 +13,5 @@ is.dgCMatrix(mat)
 A boolean variable
 }
 \description{
-Check if a matrix is `dgCMatrix` class
+Check if a matrix is \code{dgCMatrix} class
 }
diff --git a/man/nearest.neighbors.Rd b/man/nearest.neighbors.Rd
index 8e943df..8f8e304 100644
--- a/man/nearest.neighbors.Rd
+++ b/man/nearest.neighbors.Rd
@@ -12,8 +12,8 @@ for which the algorithm is applied iteratively}
 
 \item{graph}{an (\eqn{n \times n})-dimensional numeric non-negative adjacence
 \code{\link[base]{matrix}} (or
-\code{\link[Matrix:dgCMatrix-class]{dgCMatrix}}, \link[base]{vector})
-representing the graph}
+\code{\link[Matrix:dgCMatrix-class]{dgCMatrix}},
+\code{\link[base]{vector}}) representing the graph}
 
 \item{k}{the depth of the nearest neighbor search, e.g. the depth of the
 graph traversal}
@@ -26,7 +26,7 @@ returns the kNN nodes as list of integer vectors of node indexes
 \description{
 For every node in a set of nodes the graph gets traversed along the node's
 shortest paths to its neighbors. Nearest neighbors are added until a maximum
-depth of \code{k} is reached. For settings where there are more than \code{k}
+depth of \eqn{k} is reached. For settings where there are more than \eqn{k}
 neighbors having the same distance, all neighbors are returned.
 }
 \examples{
diff --git a/man/normalize.stochastic.Rd b/man/normalize.stochastic.Rd
index 8b49422..46e0155 100644
--- a/man/normalize.stochastic.Rd
+++ b/man/normalize.stochastic.Rd
@@ -8,14 +8,13 @@ normalize.stochastic(obj, ...)
 }
 \arguments{
 \item{obj}{\code{\link[base]{matrix}} (or
-\code{\link[Matrix:dgCMatrix-class]{dgCMatrix}}, \link[base]{vector}) that
-is stochstically normalized}
+\code{\link[Matrix:dgCMatrix-class]{dgCMatrix}},
+\code{\link[base]{vector}}) that is stochstically normalized}
 
 \item{...}{additional params}
 }
 \value{
-returns the normalized \code{\link[base]{matrix}} (or
-  \code{\link[Matrix:dgCMatrix-class]{dgCMatrix}}, \link[base]{vector})
+returns the normalized matrix/vector)
 }
 \description{
 Create a stochastically normalized matrix/vector
diff --git a/man/random.walk.Rd b/man/random.walk.Rd
index d64d790..4220e4c 100644
--- a/man/random.walk.Rd
+++ b/man/random.walk.Rd
@@ -17,16 +17,16 @@ random.walk(
 )
 }
 \arguments{
-\item{p0}{an \code{n x p}-dimensional numeric non-negative vector/matrix
+\item{p0}{an \eqn{n \times p}-dimensional numeric non-negative vector/matrix
 representing the starting distribution of the Markov chain
 (does not need to sum to one).}
 
-\item{graph}{an (\code{n x n})-dimensional numeric non-negative adjacence
+\item{graph}{an \eqn{n \times p}-dimensional numeric non-negative adjacence
 \code{\link[base]{matrix}} (or
 \code{\link[Matrix:dgCMatrix-class]{dgCMatrix}}) representing the graph}
 
-\item{r}{a scalar between (0, 1). restart probability if a Markov random
-walk with restart is desired}
+\item{r}{a scalar between \eqn{(0, 1)}. restart probability if a Markov
+random walk with restart is desired}
 
 \item{niter}{maximal number of iterations for computation of the
 Markov chain. If \code{thresh} is not reached, then \code{niter} is used as
@@ -44,7 +44,8 @@ computed solving the analytical solution or rather iteratively}
 \item{correct.for.hubs}{if \code{TRUE} multiplies a correction factor to the
 nodes, such that the random walk gets not biased to nodes with high
 degree. In that case the original input matrix will be normalized as:
-\deqn{ P(j | i) =  1 /degree(i) *  min(1, degree(j)/degree(j))}
+\deqn{ P(j | i) = \dfrac{1}{\text{degree}(i)} \times
+  \min \left(1, \dfrac{\text{degree}(j)}{\text{degree}(i)}\right)}
 \emph{Note that this will not consider edge weights.}}
 
 \item{allow.ergodic}{Allow multiple components in a graph.}
@@ -54,15 +55,15 @@ degree. In that case the original input matrix will be normalized as:
 \value{
 returns a list with the following elements
  \itemize{
-  \item p.inf  the stationary distribution as numeric vector
-  \item transition.matrix the column normalized transition matrix used for
-        the random walk
+  \item \code{p.inf}  the stationary distribution as numeric vector
+  \item \code{transition.matrix} the column normalized transition matrix used
+        for the random walk
  }
 }
 \description{
-A Markov Random Walk takes an inital distribution \code{p0}
+A Markov Random Walk takes an inital distribution \eqn{p_0}
 and calculates the stationary distribution of that.
-The diffusion process is regulated by a restart probability \code{r} which
+The diffusion process is regulated by a restart probability \eqn{r} which
 controls how often the MRW jumps back to the initial values.
 }
 \examples{