diff --git a/DESCRIPTION b/DESCRIPTION
index aa46d92..e3ea3a3 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,7 +1,7 @@
 Package: miaTime
 Type: Package
 Title: Microbiome Time Series Analysis
-Version: 0.99.3
+Version: 0.99.4
 Authors@R:
      c(person(given = "Leo", family = "Lahti", role = c("aut"),
               email = "leo.lahti@iki.fi",
diff --git a/NAMESPACE b/NAMESPACE
index af83c00..b3a95c3 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -2,15 +2,19 @@
 
 export(addBaselineDivergence)
 export(addBimodality)
+export(addShortTermChange)
 export(addStepwiseDivergence)
 export(getBaselineDivergence)
 export(getBimodality)
+export(getShortTermChange)
 export(getStepwiseDivergence)
 exportMethods(addBaselineDivergence)
 exportMethods(addBimodality)
+exportMethods(addShortTermChange)
 exportMethods(addStepwiseDivergence)
 exportMethods(getBaselineDivergence)
 exportMethods(getBimodality)
+exportMethods(getShortTermChange)
 exportMethods(getStepwiseDivergence)
 import(S4Vectors)
 import(SingleCellExperiment)
@@ -24,7 +28,10 @@ importFrom(dplyr,arrange)
 importFrom(dplyr,group_by)
 importFrom(dplyr,lag)
 importFrom(dplyr,mutate)
+importFrom(dplyr,select)
+importFrom(dplyr,summarise)
 importFrom(dplyr,summarize)
+importFrom(dplyr,sym)
 importFrom(dplyr,ungroup)
 importFrom(tidyr,pivot_wider)
 importFrom(tidyr,unnest)
diff --git a/NEWS b/NEWS
index 4ba2d9f..238506d 100755
--- a/NEWS
+++ b/NEWS
@@ -1,3 +1,6 @@
+Changes in version 0.99.4
++ Added functions for calculating short term change metrics
+
 Changes in version 0.99.3
 + Handle repeated samples in divergence functions by calculating average
 
@@ -25,4 +28,3 @@ Changes in version 0.1.5
 + Changed the function name getTimeDivergence (deprecated) into getStepwiseDivergence
 + Added internal utility functions
 + Added SilvermanAGutData
-
diff --git a/R/AllGenerics.R b/R/AllGenerics.R
index ed83829..3e4f15e 100644
--- a/R/AllGenerics.R
+++ b/R/AllGenerics.R
@@ -1,22 +1,22 @@
 # All generic methods are listed here
 
-#' @rdname addBaselineDivergence
+#' @rdname getBaselineDivergence
 #' @export
 setGeneric("getBaselineDivergence", signature = "x", function(x, ...)
     standardGeneric("getBaselineDivergence"))
 
-#' @rdname addBaselineDivergence
+#' @rdname getBaselineDivergence
 #' @export
 setGeneric("addBaselineDivergence", signature = "x", function(x, ...)
     standardGeneric("addBaselineDivergence"))
 
-#' @rdname addStepwiseDivergence
+#' @rdname getStepwiseDivergence
 #' @export
 #'
 setGeneric("getStepwiseDivergence", signature = c("x"), function(x, ...)
     standardGeneric("getStepwiseDivergence"))
 
-#' @rdname addStepwiseDivergence
+#' @rdname getStepwiseDivergence
 #' @export
 setGeneric("addStepwiseDivergence", signature = "x", function(x, ...)
     standardGeneric("addStepwiseDivergence"))
@@ -30,3 +30,13 @@ setGeneric("getBimodality", signature = "x", function(x, ...)
 #' @export
 setGeneric("addBimodality", signature = "x", function(x, ...)
     standardGeneric("addBimodality"))
+
+#' @rdname getShortTermChange
+#' @export
+setGeneric("addShortTermChange", signature = "x", function(x, ...)
+    standardGeneric("addShortTermChange"))
+
+#' @rdname getShortTermChange
+#' @export
+setGeneric("getShortTermChange", signature = "x", function(x, ...)
+    standardGeneric("getShortTermChange"))
diff --git a/R/getBaselineDivergence.R b/R/getBaselineDivergence.R
index a3c3fa6..be648b5 100644
--- a/R/getBaselineDivergence.R
+++ b/R/getBaselineDivergence.R
@@ -1,4 +1,6 @@
-#' @name addBaselineDivergence
+#' @name
+#' getBaselineDivergence
+#'
 #' @export
 #'
 #' @title
@@ -63,7 +65,6 @@
 #'
 #' @examples
 #' library(miaTime)
-#' library(mia)
 #'
 #' data(hitchip1006)
 #' tse <- transformAssay(hitchip1006, method = "relabundance")
@@ -94,7 +95,7 @@
 #'
 NULL
 
-#' @rdname addBaselineDivergence
+#' @rdname getBaselineDivergence
 #' @export
 setMethod("getBaselineDivergence", signature = c(x = "SummarizedExperiment"),
     function(
@@ -155,7 +156,7 @@ setMethod("getBaselineDivergence", signature = c(x = "SummarizedExperiment"),
     }
 )
 
-#' @rdname addBaselineDivergence
+#' @rdname getBaselineDivergence
 #' @export
 setMethod("addBaselineDivergence", signature = c(x = "SummarizedExperiment"),
     function(
@@ -356,6 +357,10 @@ setMethod("addBaselineDivergence", signature = c(x = "SummarizedExperiment"),
 # This function get time difference between a sample and its reference sample
 #' @importFrom dplyr group_by mutate
 .get_time_difference <- function(x, time.col, reference, ...){
+    # This following line is to suppress "no visible binding for" messages
+    # in cmdcheck
+    temp_sample <- ":=" <- time <- .data <- NULL
+
     # Get timepoints
     time_point <- x[[time.col]]
     # Get reference time points
diff --git a/R/getBimodality.R b/R/getBimodality.R
index 39873aa..17b4c8d 100644
--- a/R/getBimodality.R
+++ b/R/getBimodality.R
@@ -58,7 +58,6 @@
 #'
 #' @examples
 #' library(miaTime)
-#' library(mia)
 #'
 #' data(SilvermanAGutData)
 #' tse <- SilvermanAGutData
diff --git a/R/getShortTermChange.R b/R/getShortTermChange.R
new file mode 100644
index 0000000..3feb2b1
--- /dev/null
+++ b/R/getShortTermChange.R
@@ -0,0 +1,167 @@
+#' @name
+#' getShortTermChange
+#'
+#' @export
+#'
+#' @title
+#' Short term changes in abundance
+#'
+#' @description
+#' Calculates short term changes in abundance of taxa using temporal
+#' abundance data.
+#'
+#' @details
+#' These functions can be utilized to calculate growth metrics for short term
+#' change. In specific, the functions calculate the metrics with the following
+#' equations:
+#'
+#' \deqn{time\_diff = time_{t} - time_{t-1}}
+#'
+#' \deqn{abundance\_diff = abundance_{t} - abundance_{t-1}}
+#'
+#' \deqn{growth\_rate = abundance\_diff - abundance_{t-1}}
+#'
+#' \deqn{rate\_of\_change = abundance\_diff - time\_diff}
+#'
+#' @references
+#' Ji, B.W., et al. (2020) Macroecological dynamics of gut microbiota.
+#' Nat Microbiol 5, 768–775 . doi: https://doi.org/10.1038/s41564-020-0685-1
+#'
+#' @return
+#' \code{getShortTermChange} returns \code{DataFrame} object containing
+#' the short term change in abundance over time for a microbe.
+#' \code{addShortTermChange}, on the other hand, returns a
+#' \code{\link[SummarizedExperiment:SummarizedExperiment]{SummarizedExperiment}}
+#' object with these results in its \code{metadata}.
+#'
+#' @inheritParams addBaselineDivergence
+#'
+#' @param name \code{Character scalar}. Specifies a name for storing
+#' short term results. (Default: \code{"short_term_change"})
+#'
+#' @param ... additional arguments.
+#' \itemize{
+#'   \item \code{time.interval}: \code{Integer scalar}. Indicates the increment
+#'   between time steps. By default, the function compares each sample to the
+#'   previous one. If you need to take every second, every third, or so, time
+#'   step, then increase this accordingly. (Default: \code{1L})
+#' }
+#
+#' @examples
+#' library(miaTime)
+#'
+#' # Load time series data
+#' data(minimalgut)
+#' tse <- minimalgut
+#'
+#' # Get relative abundances
+#' tse <- transformAssay(tse, method = "relabundance")
+#' # Calculate short term changes
+#' df <- getShortTermChange(
+#'   tse, assay.type = "relabundance", time.col = "Time.hr",
+#'   group = "StudyIdentifier")
+#'
+#' # Calculate the logarithm of the ratio described in Ji, B.W., et al. (2020)
+#' tse <- transformAssay(
+#'     tse, assay.type = "relabundance", method = "log10", pseudocount = TRUE)
+#' df <- getShortTermChange(
+#'   tse, assay.type = "log10", time.col = "Time.hr", group = "StudyIdentifier")
+#'
+#' # Select certain bacteria that have highest growth rate
+#' select <- df[["growth_rate"]] |> abs() |> order(decreasing = FALSE)
+#' select <- df[select, "FeatureID"] |> unique() |> head(3)
+#' df <- df[ which(df[["FeatureID"]] %in% select), ]
+#'
+#' # Plot results
+#' library(ggplot2)
+#' p <- ggplot(df, aes(x = Time.hr, y = growth_rate, colour = FeatureID)) +
+#'     geom_smooth(level = 0.5) +
+#'     facet_grid(. ~ StudyIdentifier, scales = "free") +
+#'     scale_y_log10()
+#' p
+#'
+#' @seealso
+#' \code{\link[miaTime:getStepwiseDivergence]{getStepwiseDivergence()}}
+NULL
+
+#' @rdname getShortTermChange
+#' @export
+setMethod("addShortTermChange", signature = c(x = "SummarizedExperiment"),
+    function(x, name = "short_term_change", ...){
+        temp <- .check_input(name, list("character scalar"))
+        x <- .check_and_get_altExp(x, ...)
+        args <- c(list(x = x), list(...)[!names(list(...)) %in% c("altexp")])
+        # Calculate short term change
+        res <- do.call(getShortTermChange, args)
+        # Add to metadata
+        x <- .add_values_to_metadata(x, name, res, ...)
+        return(x)
+    }
+)
+
+#' @rdname getShortTermChange
+#' @export
+setMethod("getShortTermChange", signature = c(x = "SummarizedExperiment"),
+    function(x, time.col, assay.type = "counts", group = NULL, ...){
+        ############################## Input check #############################
+        x <- .check_and_get_altExp(x, ...)
+        temp <- .check_input(
+            time.col, list("character scalar"), colnames(colData(x)))
+        if( !is.numeric(x[[time.col]]) ){
+            stop("'time.col' must specify numeric column from colData(x)",
+                call. = FALSE)
+        }
+        .check_assay_present(assay.type, x)
+        temp <- .check_input(
+            group, list(NULL, "character scalar"), colnames(colData(x)))
+        ########################### Input check end ############################
+        # Get data in long format
+        df <- meltSE(x, assay.type, add.col = c(time.col, group))
+        # Calculate metrics
+        res <- .calculate_growth_metrics(df, assay.type, time.col, group, ...)
+        return(res)
+    }
+)
+
+################################ HELP FUNCTIONS ################################
+
+# This function calculates the growth metrics from the data.frame.
+#' @importFrom dplyr arrange group_by sym summarise mutate select
+.calculate_growth_metrics <- function(
+        df, assay.type, time.col, group, time.interval = 1L, ...) {
+    # This following line is to suppress "no visible binding for" messages
+    # in cmdcheck
+    FeatureID <- .data <- value <- abundance_diff <- time_diff <- NULL
+    #
+    temp <- .check_input(time.interval, list("numeric scalar"))
+    # If there are replicated samples, give warning that average is calculated
+    if( anyDuplicated(df[, c("FeatureID", group, time.col)]) ){
+        message("The dataset contains replicated samples. The average ",
+                "abundance is calculated for each time point.")
+    }
+    # Sort data based on time
+    df <- df |> arrange( !!sym(time.col) )
+    # Group based on features and time points. If group was specified, take that
+    # also into account
+    if( !is.null(group) ){
+        df <- df |> group_by(!!sym(group), FeatureID, !!sym(time.col))
+    } else{
+        df <- df |> group_by(FeatureID, !!sym(time.col))
+    }
+    df <- df |>
+        # Summarize duplicated samples by taking an average
+        summarise(value = mean(.data[[assay.type]], na.rm = TRUE)) |>
+        # For each feature in a sample group, calculate growth metrics
+        mutate(
+            time_diff = !!sym(time.col) -
+                lag(!!sym(time.col), n = time.interval),
+            abundance_diff = value - lag(value, n = time.interval),
+            growth_rate = abundance_diff / lag(value, n = time.interval),
+            rate_of_change = abundance_diff / time_diff
+        ) |>
+        # Remove value column that includes average abundances
+        select(-value) |>
+        # Convert to DataFrame
+        DataFrame()
+    return(df)
+}
diff --git a/R/getStepwiseDivergence.R b/R/getStepwiseDivergence.R
index d2bef69..cd735d5 100644
--- a/R/getStepwiseDivergence.R
+++ b/R/getStepwiseDivergence.R
@@ -1,4 +1,6 @@
-#' @name addStepwiseDivergence
+#' @name
+#' getStepwiseDivergence
+#'
 #' @export
 #'
 #' @title
@@ -50,7 +52,7 @@
 #'
 NULL
 
-#' @rdname addStepwiseDivergence
+#' @rdname getStepwiseDivergence
 #' @export
 setMethod("getStepwiseDivergence", signature = c(x = "ANY"),
     function(
@@ -107,7 +109,7 @@ setMethod("getStepwiseDivergence", signature = c(x = "ANY"),
     }
 )
 
-#' @rdname addStepwiseDivergence
+#' @rdname getStepwiseDivergence
 #' @export
 setMethod("addStepwiseDivergence", signature = c(x = "SummarizedExperiment"),
     function(
diff --git a/R/miaTime.R b/R/miaTime.R
index 8a554bd..5e1b62f 100755
--- a/R/miaTime.R
+++ b/R/miaTime.R
@@ -3,7 +3,8 @@
 #' \code{miaTime} implements time series methods in \link[mia]{mia} ecosystem.
 #' @name miaTime-package
 #' @aliases miaTime
-#' @seealso \link[TreeSummarizedExperiment:TreeSummarizedExperiment-class]{TreeSummarizedExperiment} class
+#' @seealso
+#' \link[TreeSummarizedExperiment:TreeSummarizedExperiment-class]{TreeSummarizedExperiment} class
 "_PACKAGE"
 NULL
 
diff --git a/R/utils.R b/R/utils.R
index c9ebfc3..db7d017 100644
--- a/R/utils.R
+++ b/R/utils.R
@@ -176,3 +176,4 @@
 .check_assay_present <- mia:::.check_assay_present
 .add_values_to_colData <- mia:::.add_values_to_colData
 .check_and_get_altExp <- mia:::.check_and_get_altExp
+.add_values_to_metadata <- mia:::.add_values_to_metadata
diff --git a/man/addBaselineDivergence.Rd b/man/getBaselineDivergence.Rd
similarity index 99%
rename from man/addBaselineDivergence.Rd
rename to man/getBaselineDivergence.Rd
index d8f487f..83835f4 100644
--- a/man/addBaselineDivergence.Rd
+++ b/man/getBaselineDivergence.Rd
@@ -87,7 +87,6 @@ group (named list per group).
 }
 \examples{
 library(miaTime)
-library(mia)
 
 data(hitchip1006)
 tse <- transformAssay(hitchip1006, method = "relabundance")
diff --git a/man/getBimodality.Rd b/man/getBimodality.Rd
index 5bb53f3..d2730a2 100644
--- a/man/getBimodality.Rd
+++ b/man/getBimodality.Rd
@@ -69,7 +69,6 @@ bimodality and abundance to determine conditionally rare taxa (CRT).
 }
 \examples{
 library(miaTime)
-library(mia)
 
 data(SilvermanAGutData)
 tse <- SilvermanAGutData
diff --git a/man/getShortTermChange.Rd b/man/getShortTermChange.Rd
new file mode 100644
index 0000000..3ea241e
--- /dev/null
+++ b/man/getShortTermChange.Rd
@@ -0,0 +1,108 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/AllGenerics.R, R/getShortTermChange.R
+\name{addShortTermChange}
+\alias{addShortTermChange}
+\alias{getShortTermChange}
+\alias{addShortTermChange,SummarizedExperiment-method}
+\alias{getShortTermChange,SummarizedExperiment-method}
+\title{Short term changes in abundance}
+\usage{
+addShortTermChange(x, ...)
+
+getShortTermChange(x, ...)
+
+\S4method{addShortTermChange}{SummarizedExperiment}(x, name = "short_term_change", ...)
+
+\S4method{getShortTermChange}{SummarizedExperiment}(x, time.col, assay.type = "counts", group = NULL, ...)
+}
+\arguments{
+\item{x}{A
+\code{\link[SummarizedExperiment:SummarizedExperiment-class]{SummarizedExperiment}}
+object.}
+
+\item{...}{additional arguments.
+\itemize{
+\item \code{time.interval}: \code{Integer scalar}. Indicates the increment
+between time steps. By default, the function compares each sample to the
+previous one. If you need to take every second, every third, or so, time
+step, then increase this accordingly. (Default: \code{1L})
+}}
+
+\item{name}{\code{Character scalar}. Specifies a name for storing
+short term results. (Default: \code{"short_term_change"})}
+
+\item{time.col}{\code{Character scalar}. Specifies a name of the column from
+\code{colData} that identifies the sampling time points for the samples.}
+
+\item{assay.type}{\code{Character scalar}. Specifies which assay values are
+used in the dissimilarity estimation. (Default: \code{"counts"})}
+
+\item{group}{\code{Character scalar}. Specifies a name of the column from
+\code{colData} that identifies the grouping of the samples.
+(Default: \code{NULL})}
+}
+\value{
+\code{getShortTermChange} returns \code{DataFrame} object containing
+the short term change in abundance over time for a microbe.
+\code{addShortTermChange}, on the other hand, returns a
+\code{\link[SummarizedExperiment:SummarizedExperiment]{SummarizedExperiment}}
+object with these results in its \code{metadata}.
+}
+\description{
+Calculates short term changes in abundance of taxa using temporal
+abundance data.
+}
+\details{
+These functions can be utilized to calculate growth metrics for short term
+change. In specific, the functions calculate the metrics with the following
+equations:
+
+\deqn{time\_diff = time_{t} - time_{t-1}}
+
+\deqn{abundance\_diff = abundance_{t} - abundance_{t-1}}
+
+\deqn{growth\_rate = abundance\_diff - abundance_{t-1}}
+
+\deqn{rate\_of\_change = abundance\_diff - time\_diff}
+}
+\examples{
+library(miaTime)
+
+# Load time series data
+data(minimalgut)
+tse <- minimalgut
+
+# Get relative abundances
+tse <- transformAssay(tse, method = "relabundance")
+# Calculate short term changes
+df <- getShortTermChange(
+  tse, assay.type = "relabundance", time.col = "Time.hr",
+  group = "StudyIdentifier")
+
+# Calculate the logarithm of the ratio described in Ji, B.W., et al. (2020)
+tse <- transformAssay(
+    tse, assay.type = "relabundance", method = "log10", pseudocount = TRUE)
+df <- getShortTermChange(
+  tse, assay.type = "log10", time.col = "Time.hr", group = "StudyIdentifier")
+
+# Select certain bacteria that have highest growth rate
+select <- df[["growth_rate"]] |> abs() |> order(decreasing = FALSE)
+select <- df[select, "FeatureID"] |> unique() |> head(3)
+df <- df[ which(df[["FeatureID"]] \%in\% select), ]
+
+# Plot results
+library(ggplot2)
+p <- ggplot(df, aes(x = Time.hr, y = growth_rate, colour = FeatureID)) +
+    geom_smooth(level = 0.5) +
+    facet_grid(. ~ StudyIdentifier, scales = "free") +
+    scale_y_log10()
+p
+
+}
+\references{
+Ji, B.W., et al. (2020) Macroecological dynamics of gut microbiota.
+Nat Microbiol 5, 768–775 . doi: https://doi.org/10.1038/s41564-020-0685-1
+}
+\seealso{
+\code{\link[miaTime:getStepwiseDivergence]{getStepwiseDivergence()}}
+}
diff --git a/man/addStepwiseDivergence.Rd b/man/getStepwiseDivergence.Rd
similarity index 100%
rename from man/addStepwiseDivergence.Rd
rename to man/getStepwiseDivergence.Rd
diff --git a/tests/testthat/test-getShortTermChange.R b/tests/testthat/test-getShortTermChange.R
new file mode 100644
index 0000000..ce9672c
--- /dev/null
+++ b/tests/testthat/test-getShortTermChange.R
@@ -0,0 +1,213 @@
+
+test_that("getShortTermChange errors", {
+    tse <- makeTSE(nrow = 100, ncol = 20)
+    assayNames(tse) <- "counts"
+    tse[["Time"]] <- sample(seq(1, 5), 20, replace = TRUE)
+
+    df <- getShortTermChange(tse, time.col = "Test", assay.type = "counts") |>
+        expect_error()
+    df <- getShortTermChange(tse, time.col = "group", assay.type = "counts") |>
+        expect_error()
+    df <- getShortTermChange(tse, time.col = "group", assay.type = "test") |>
+        expect_error()
+    df <- getShortTermChange(tse, time.col = "Time", group = "counts") |>
+        expect_error()
+    df <- getShortTermChange(tse, time.col = "Time", time.interval = TRUE) |>
+        expect_error()
+    tse <- addShortTermChange(tse, time.col = "Time", name = TRUE) |>
+        expect_error()
+})
+
+# Test with time interval 1 and no grouping
+test_that("getShortTermChange calculates correctly without group", {
+    tse <- makeTSE(nrow = 100, ncol = 20)
+    assayNames(tse) <- "counts"
+    tse[["Time"]] <- sample(seq(1, 5), 20, replace = TRUE)
+
+    df <- getShortTermChange(tse, time.col = "Time", assay.type = "counts") |>
+        expect_message()
+
+    # Check the structure of the returned object
+    expect_s4_class(df, "DataFrame")
+    expect_true("time_diff" %in% colnames(df))
+    expect_true("abundance_diff" %in% colnames(df))
+    expect_true("growth_rate" %in% colnames(df))
+    expect_true("rate_of_change" %in% colnames(df))
+    expect_true("Time" %in% colnames(df))
+    expect_true("FeatureID" %in% colnames(df))
+
+    # Test 10 times with random combinations
+    for( i in seq_len(10) ){
+        # Get random feature and timepoint to test
+        feat <- sample(rownames(tse), 1)
+        time_point <- sample(tse[["Time"]], 1)
+        # Get index of previous sample
+        temp_df <- colData(tse)[ order(tse[["Time"]]), ]
+        time_points <- unique(temp_df[["Time"]])
+        prev_time <- which(time_points == time_point)-1
+        ref <- rep(NA_real_, 4)
+        # Calculate reference values manually if previous sample exists
+        if( prev_time > 0 ){
+            prev_time <- time_points[ prev_time ]
+            prev_sam <- rownames(temp_df[temp_df[["Time"]] == prev_time, ])
+            sam <- rownames(temp_df[temp_df[["Time"]] == time_point, ])
+            # Take certain samples
+            tse_curr <- tse[, sam]
+            tse_prev <- tse[, prev_sam]
+            # Get values
+            time_diff <- time_point - prev_time
+            val1 <- mean(assay(tse_curr, "counts")[feat, ], na.rm = TRUE)
+            val2 <- mean(assay(tse_prev, "counts")[feat, ], na.rm = TRUE)
+            ref <- c(
+                time_diff,
+                val1 - val2,
+                (val1 - val2) / val2,
+                (val1 - val2) / time_diff
+            )
+        }
+        names(ref) <- c(
+            "time_diff", "abundance_diff", "growth_rate", "rate_of_change")
+        res <- df[df[["FeatureID"]] == feat & df[["Time"]] == time_point, ]
+        res <- unlist(res[, names(ref)])
+        #
+        expect_equal(res, ref)
+    }
+})
+
+# Test with time interval 2 and no grouping
+test_that("getShortTermChange calculates correctly with time interval 2", {
+    tse <- makeTSE(nrow = 1000, ncol = 200)
+    assayNames(tse) <- "counts"
+    tse[["Time"]] <- sample(seq(1, 5), 200, replace = TRUE)
+    tse <- transformAssay(tse, method = "relabundance")
+    df <- getShortTermChange(
+        tse, time.col = "Time", assay.type = "relabundance",
+        time.interval = 2) |>
+        expect_message()
+
+    # Check the structure of the returned object
+    expect_s4_class(df, "DataFrame")
+    expect_true("time_diff" %in% colnames(df))
+    expect_true("abundance_diff" %in% colnames(df))
+    expect_true("growth_rate" %in% colnames(df))
+    expect_true("rate_of_change" %in% colnames(df))
+    expect_true("Time" %in% colnames(df))
+    expect_true("FeatureID" %in% colnames(df))
+
+    # Test 10 times with random combinations
+    for( i in seq_len(10) ){
+        # Get random feature and timepoint to test
+        feat <- sample(rownames(tse), 1)
+        time_point <- sample(tse[["Time"]], 1)
+        # Get index of previous sample
+        temp_df <- colData(tse)[ order(tse[["Time"]]), ]
+        time_points <- unique(temp_df[["Time"]])
+        prev_time <- which(time_points == time_point)-2
+        ref <- rep(NA_real_, 4)
+        # Calculate reference values manually if previous sample exists
+        if( prev_time > 0 ){
+            prev_time <- time_points[ prev_time ]
+            prev_sam <- rownames(temp_df[temp_df[["Time"]] == prev_time, ])
+            sam <- rownames(temp_df[temp_df[["Time"]] == time_point, ])
+            # Take certain samples
+            tse_curr <- tse[, sam]
+            tse_prev <- tse[, prev_sam]
+            # Get values
+            time_diff <- time_point - prev_time
+            val1 <- mean(assay(tse_curr, "relabundance")[feat, ], na.rm = TRUE)
+            val2 <- mean(assay(tse_prev, "relabundance")[feat, ], na.rm = TRUE)
+            ref <- c(
+                time_diff,
+                val1 - val2,
+                (val1 - val2) / val2,
+                (val1 - val2) / time_diff
+            )
+        }
+        names(ref) <- c(
+            "time_diff", "abundance_diff", "growth_rate", "rate_of_change")
+        res <- df[df[["FeatureID"]] == feat & df[["Time"]] == time_point, ]
+        res <- unlist(res[, names(ref)])
+        #
+        expect_equal(res, ref)
+    }
+})
+
+# Test with grouping
+test_that("getShortTermChange calculates correctly with group", {
+    tse <- makeTSE(nrow = 1000, ncol = 200)
+    assayNames(tse) <- "counts"
+    tse[["Time"]] <- sample(seq(1, 5), 200, replace = TRUE)
+    # Remove duplicated samples
+    tse <- tse[ , !duplicated(colData(tse)[, c("group", "Time")])]
+
+    df <- getShortTermChange(
+        tse, time.col = "Time", assay.type = "counts", group = "group") |>
+        expect_no_message()
+
+    # Check the structure of the returned object
+    expect_s4_class(df, "DataFrame")
+    expect_true("time_diff" %in% colnames(df))
+    expect_true("abundance_diff" %in% colnames(df))
+    expect_true("growth_rate" %in% colnames(df))
+    expect_true("rate_of_change" %in% colnames(df))
+    expect_true("Time" %in% colnames(df))
+    expect_true("FeatureID" %in% colnames(df))
+
+    # Test 10 times with random combinations
+    for( i in seq_len(10) ){
+        # Get random feature and timepoint to test
+        feat <- sample(rownames(tse), 1)
+        time_point <- sample(tse[["Time"]], 1)
+        group <- sample(tse[["group"]], 1)
+        # Get only specific group
+        temp_df <- colData(tse)[ tse[["group"]] == group, ]
+        # Get index of previous sample
+        temp_df <- temp_df[ order(temp_df[["Time"]]), ]
+        time_points <- unique(temp_df[["Time"]])
+        prev_time <- which(time_points == time_point)-1
+        ref <- rep(NA_real_, 4)
+        # Calculate reference values manually if previous sample exists
+        if( prev_time > 0 ){
+            prev_time <- time_points[ prev_time ]
+            prev_sam <- rownames(temp_df[temp_df[["Time"]] == prev_time, ])
+            sam <- rownames(temp_df[temp_df[["Time"]] == time_point, ])
+            # Take certain samples
+            tse_curr <- tse[, sam]
+            tse_prev <- tse[, prev_sam]
+            # Get values
+            time_diff <- time_point - prev_time
+            val1 <- mean(assay(tse_curr, "counts")[feat, ], na.rm = TRUE)
+            val2 <- mean(assay(tse_prev, "counts")[feat, ], na.rm = TRUE)
+            ref <- c(
+                time_diff,
+                val1 - val2,
+                (val1 - val2) / val2,
+                (val1 - val2) / time_diff
+            )
+        }
+        names(ref) <- c(
+            "time_diff", "abundance_diff", "growth_rate", "rate_of_change")
+        select <- df[["FeatureID"]] == feat & df[["Time"]] == time_point &
+            df[["group"]] == group
+        res <- df[select, ]
+        res <- unlist(res[, names(ref)])
+        #
+        expect_equal(res, ref)
+    }
+})
+
+# Test that getShortTimeChange and addShortTimeChange are equal
+test_that("get* and add* are equal", {
+    tse <- makeTSE(nrow = 1000, ncol = 200)
+    assayNames(tse) <- "counts"
+    tse[["Time"]] <- sample(seq(1, 5), 200, replace = TRUE)
+    # Remove duplicated samples
+    tse <- tse[ , !duplicated(colData(tse)[, c("Time")])]
+
+    df <- getShortTermChange(tse, time.col = "Time") |>
+        expect_no_message()
+    tse <- addShortTermChange(tse, time.col = "Time", name = "test") |>
+        expect_no_message()
+
+    expect_equal(metadata(tse)[["test"]], df)
+})