diff --git a/DESCRIPTION b/DESCRIPTION
index ac13e0d..f1f4ac4 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,6 +1,6 @@
 Package: ripserr
 Title: Calculate Persistent Homology with Ripser-Based Engines
-Version: 1.0.0
+Version: 1.0.0.0003
 Authors@R:
     c(person(given = "Raoul R.",
              family = "Wadhwa",
@@ -79,7 +79,9 @@ Depends: R (>= 3.5.0)
 Imports: 
     Rcpp (>= 1.0),
     stats (>= 3.0),
-    utils (>= 3.0)
+    utils (>= 3.0),
+    phutil,
+    lifecycle
 Suggests: 
     covr (>= 3.5),
     knitr (>= 1.29),
diff --git a/NAMESPACE b/NAMESPACE
index 5e6e424..652e518 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -2,6 +2,7 @@
 
 S3method(cubical,array)
 S3method(cubical,matrix)
+S3method(cubical,numeric)
 S3method(head,PHom)
 S3method(print,PHom)
 S3method(tail,PHom)
@@ -17,6 +18,7 @@ export(as.PHom)
 export(cubical)
 export(cubical.array)
 export(cubical.matrix)
+export(cubical.numeric)
 export(is.PHom)
 export(vietoris_rips)
 export(vietoris_rips.data.frame)
diff --git a/NEWS.md b/NEWS.md
index ab6ec3b..9a98984 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -1,17 +1,31 @@
 # next version
 
-## permission for deaths before births and superlevel set filtrations of rasters
+### output class
 
-A logical argument `sublevel` has been added to `cubical` that, when `FALSE`, will pre- and post-transform raster data in order to obtain superlevel set persistent homology.
-Enabling this, an assertion that all `birth < death` has been removed from checks of persistence data.
+The new `return_class` parameter allows the user to specify whether to output persistence data in the legacy `'PHom'` class or the `'persistence'` class from {phutil}.
+It defaults to `'PHom'` but will switch to `'persistence'` when the `'PHom'` class is deprecated in a future version.
+
+## Vietoris-Rips PH
 
-## sliding window embeddings of multivariable time series (breaking change)
+### sliding window embeddings of multivariable time series (breaking change)
 
 Previously only univariable time series could be passed to `vietoris_rips()` via the sliding window embedding (used for quasi-attractor detection).
 An additional unexported embedding function has been written to handle multivariable time series. (It underperforms the original function on univariable time series, which therefore continue to rely on the original.)
 
 Furthermore, whereas `data_dim` previously defaulted to `2`, it now defaults to the number of observations per time unit of a time series as recovered by `tsp()`. (The behavior for unclassed numeric vectors remains unchanged.)
 
+## cubical PH
+
+### functionality for 1-dimensional arrays
+
+`cubical()` can now handle 1-dimensional arrays (for which no dedicated source code exists) by treating them as 2-dimensional (with an expanse of 1 in the second dimension).
+This enables the new method `cubical.numeric()` to accept vectors.
+
+### deaths before births and superlevel set filtrations
+
+A logical argument `sublevel` has been added to `cubical` that, when `FALSE`, will pre- and post-transform raster data in order to obtain superlevel set persistent homology.
+Enabling this, an assertion that all `birth < death` has been removed from checks of persistence data.
+
 # ripserr 1.0.0
 
 This major version replaces an outdated version of the Ripser C++ library with its current version.
diff --git a/R/PHom.R b/R/PHom.R
index 3b9172a..aced7b9 100644
--- a/R/PHom.R
+++ b/R/PHom.R
@@ -11,6 +11,13 @@ new_PHom <- function(x = data.frame(dimension = integer(0),
                      dim_col = 1,
                      b_col = 2,
                      d_col = 3) {
+  lifecycle::deprecate_soft(
+    "1.1.0",
+    I("'PHom' class"),
+    with = I("'persistence' from the {phutil} package"),
+    id = "PHom"
+  )
+  
   # parameter validity checked in PHom helper
   
   # construct df for PHom object
@@ -118,6 +125,13 @@ valid_colval <- function(df, val, val_name) {
 #' # print feature details to confirm accuracy
 #' print.data.frame(df_phom)
 PHom <- function(x, dim_col = 1, birth_col = 2, death_col = 3) {
+  lifecycle::deprecate_soft(
+    "1.1.0",
+    I("'PHom' class"),
+    with = I("'persistence' from the {phutil} package"),
+    id = "PHom"
+  )
+  
   ## basic parameter checks (column nums/names are valid, etc.)
   if (!is.data.frame(x)) {
     x <- as.data.frame(x)
@@ -163,6 +177,13 @@ PHom <- function(x, dim_col = 1, birth_col = 2, death_col = 3) {
 #' # print feature details to confirm accuracy
 #' print.data.frame(df_phom)
 as.PHom <- function(x, dim_col = 1, birth_col = 2, death_col = 3) {
+  lifecycle::deprecate_soft(
+    "1.1.0",
+    I("'PHom' class"),
+    with = I("'persistence' from the {phutil} package"),
+    id = "PHom"
+  )
+  
   x <- as.data.frame(x)
   
   return(PHom(x))
@@ -191,6 +212,13 @@ as.PHom <- function(x, dim_col = 1, birth_col = 2, death_col = 3) {
 #' # confirm that persistence data is NOT valid
 #' is.PHom(df)
 is.PHom <- function(x) {
+  lifecycle::deprecate_soft(
+    "1.1.0",
+    I("'PHom' class"),
+    with = I("'persistence' from the {phutil} package"),
+    id = "PHom"
+  )
+  
   # use validate to implement checks
   return(
     validate_PHom(x = x, error = FALSE)
@@ -243,7 +271,7 @@ print.PHom <- function(x, ...) {
                         "; max = ",
                         signif(max(x$death[is.finite(x$death)], na.rm = TRUE),
                                digits = 5),
-                        "."),
+                        ".\n"),
                  "")
   
   cat(paste(ans1, ans2, ans3, sep = "\n\n"))
diff --git a/R/cubical.R b/R/cubical.R
index 82d0dae..9ee63d3 100644
--- a/R/cubical.R
+++ b/R/cubical.R
@@ -21,22 +21,28 @@
 #' @param ... other relevant parameters
 #' @rdname cubical
 #' @export cubical
-#' @return `PHom` object
+#' @return `"PHom"` or [`"persistence"`][phutil::as_persistence] object
 #' @examples 
 #' 
+#' # 1-dim example
+#' dataset <- rnorm(1500)
+#' dim(dataset) <- 1500
+#' cubical_hom1 <- cubical(dataset)
+#' 
 #' # 2-dim example
 #' dataset <- rnorm(10 ^ 2)
 #' dim(dataset) <- rep(10, 2)
-#' cubical_hom2 <- cubical(dataset)
+#' ( cubical_hom2 <- cubical(dataset) )
 #' 
 #' # 3-dim example
 #' dataset <- rnorm(8 ^ 3)
 #' dim(dataset) <- rep(8, 3)
-#' cubical_hom3 <- cubical(dataset)
+#' ( cubical_hom3 <- cubical(dataset) )
 #' 
 #' # 4-dim example
 #' dataset <- rnorm(5 ^ 4)
 #' dim(dataset) <- rep(5, 4)
+#' ( cubical_hom4 <- cubical(dataset) )
 #' 
 #' # sublevel versus superlevel
 #' cubical(volcano)
@@ -56,12 +62,17 @@ cubical <- function(dataset, ...) {
 #'   see Kaji et al. (2020) <https://arxiv.org/abs/2005.12692> for details
 #' @param sublevel logical; whether to take the sublevel set filtration or else
 #'   the superlevel set filtration
+#' @param return_class class of output object; either `"PHom"` (default; legacy)
+#'   or `"persistence"` (from the
+#'   **[phutil](https://cran.r-project.org/package=phutil)** package)
 #' @export cubical.array
 #' @export
 cubical.array <- function(
     dataset,
-    threshold = 9999, method = "lj",
+    threshold = 9999,
+    method = "lj",
     sublevel = TRUE,
+    return_class = c("PHom", "persistence"),
     ...
 ) {
   # do this before checks since it modifies `dataset`
@@ -74,6 +85,8 @@ cubical.array <- function(
                       method = method)
   validate_arr_cub(dataset)
   
+  # if dataset is 1-dimensional, treat it as 2-dimensional
+  if (length(dim(dataset)) == 1L) dim(dataset) <- c(dim(dataset), 1L)
   
   # transform method parameter for C++ function
   method_int <- switch(method,
@@ -124,7 +137,24 @@ cubical.array <- function(
   ans$dimension <- as.integer(ans$dimension)
   
   # convert data frame to a PHom object
-  ans <- new_PHom(ans)
+  ans <- switch(
+    match.arg(return_class),
+    PHom = {
+      lifecycle::deprecate_soft(
+        "1.1.0",
+        I("'PHom' class"),
+        with = I("'persistence' from the {phutil} package"),
+        id = "PHom"
+      )
+      new_PHom(ans)
+    },
+    persistence = phutil::as_persistence(
+      ans,
+      engine = "ripserr::cubical",
+      filtration = "cubical",
+      parameters = list(threshold = threshold, method = method)
+    )
+  )
   
   # return
   return(ans)
@@ -146,4 +176,17 @@ cubical.matrix <- function(dataset, ...) {
   
   # return
   return(ans)
-}
\ No newline at end of file
+}
+
+#' @rdname cubical
+#' @export cubical.numeric
+#' @export
+cubical.numeric <- function(dataset, ...) {
+  # convert the numeric vector to a 1-dimensional array
+  dataset <- as.array(dataset, dim = 1L)
+  
+  # calculate persistent homology using cubical.array
+  ans <- cubical.array(dataset, ...)
+  
+  return(ans)
+}
diff --git a/R/utility.R b/R/utility.R
index 30026cb..539e077 100644
--- a/R/utility.R
+++ b/R/utility.R
@@ -122,7 +122,7 @@ validate_arr_cub <- function(dataset) {
   error_class(dataset, "dataset", "array")
   
   # dataset should have either 2, 3, or 4 dimensions (only ones supported)
-  if (!(length(dim(dataset)) %in% c(2, 3, 4))) {
+  if (!(length(dim(dataset)) %in% seq(4))) {
     stop(paste("dataset parameter must have either 2, 3, or 4 dimensions,",
                "passed argument has", length(dim(dataset)), "dimensions"))
   }
@@ -139,7 +139,11 @@ validate_arr_cub <- function(dataset) {
   }
   
   # make sure dataset is not too large
-  if (length(dim(dataset)) == 2) {
+  if (length(dim(dataset)) == 1L) {
+    if (dim(dataset) > 2000) {
+      stop(paste("Max size for dim 1 = 2000; passed size =", dim(dataset)))
+    }
+  } else if (length(dim(dataset)) == 2) {
     if (dim(dataset)[1] > 2000 |
         dim(dataset)[2] > 1000) {
       stop(paste("Max size for dim 2 = 2000 x 1000; passed size =",
diff --git a/R/vietoris_rips.R b/R/vietoris_rips.R
index 5f9d79a..0d0dbbb 100644
--- a/R/vietoris_rips.R
+++ b/R/vietoris_rips.R
@@ -34,7 +34,7 @@
 #' @param ... other relevant parameters
 #' @rdname vietoris_rips
 #' @export vietoris_rips
-#' @return `PHom` object
+#' @return `"PHom"` or [`"persistence"`][phutil::as_persistence] object
 #' @examples
 #'
 #' # create a 2-d point cloud of a circle (100 points)
@@ -43,7 +43,7 @@
 #' pt.cloud <- cbind(cos(rand.angle), sin(rand.angle))
 #'
 #' # calculate persistent homology (num.pts by 3 numeric matrix)
-#' pers.hom <- vietoris_rips(pt.cloud)
+#' ( pers.hom <- vietoris_rips(pt.cloud) )
 #' 
 #' # sliding window persistent homology
 #' ( ld.phom <- vietoris_rips(ldeaths, data_dim = 12) )
@@ -83,13 +83,16 @@ vietoris_rips.data.frame <- function(dataset, ...) {
   return(ans)
 }
 
+#' @rdname vietoris_rips
 #' @param max_dim maximum dimension of persistent homology features to be
 #'   calculated
 #' @param dim deprecated; passed to `max_dim` or ignored if `max_dim` is
 #'   specified
 #' @param threshold maximum simplicial complex diameter to explore
 #' @param p prime field in which to calculate persistent homology
-#' @rdname vietoris_rips
+#' @param return_class class of output object; either `"PHom"` (default; legacy)
+#'   or `"persistence"` (from the
+#'   **[phutil](https://cran.r-project.org/package=phutil)** package)
 #' @export vietoris_rips.matrix
 #' @export
 vietoris_rips.matrix <- function(
@@ -98,11 +101,31 @@ vietoris_rips.matrix <- function(
     threshold = -1,
     p = 2L,
     dim = NULL,
+    return_class = c("PHom", "persistence"),
     ...
 ) {
   
   # shortcut for special case (only 1 row should return empty PHom)
-  if (nrow(dataset) == 1L) return(new_PHom())
+  if (nrow(dataset) == 1L) {
+    return(switch(
+      match.arg(return_class),
+      PHom = {
+        lifecycle::deprecate_soft(
+          "1.1.0",
+          I("'PHom' class"),
+          with = I("'persistence' from the {phutil} package"),
+          id = "PHom"
+        )
+        new_PHom()
+      },
+      persistence = phutil::as_persistence(
+        matrix(NA_real_, nrow = 0L, ncol = 3L),
+        engine = "ripserr::vietoris_rips",
+        filtration = "Vietoris-Rips",
+        parameters = list(max_dim = max_dim, threshold = threshold, p = p)
+      )
+    ))
+  }
   
   # handle `dim` if passed
   if (! is.null(dim)) {
@@ -132,7 +155,24 @@ vietoris_rips.matrix <- function(
   ans <- ripser_cpp_dist(dataset, max_dim, threshold, 1., p)
   
   # coerce to 'PHom' class
-  ans <- new_PHom(ripser_ans_to_df(ans))
+  ans <- switch(
+    match.arg(return_class),
+    PHom = {
+      lifecycle::deprecate_soft(
+        "1.1.0",
+        I("'PHom' class"),
+        with = I("'persistence' from the {phutil} package"),
+        id = "PHom"
+      )
+      new_PHom(ripser_ans_to_df(ans))
+    },
+    persistence = phutil::as_persistence(
+      ans,
+      engine = "ripserr::vietoris_rips",
+      filtration = "Vietoris-Rips",
+      parameters = list(max_dim = max_dim, threshold = threshold, p = p)
+    )
+  )
   
   # return
   return(ans)
@@ -147,6 +187,7 @@ vietoris_rips.dist <- function(
     threshold = -1,
     p = 2L,
     dim = NULL,
+    return_class = c("PHom", "persistence"),
     ...
 ) {
   
@@ -178,12 +219,30 @@ vietoris_rips.dist <- function(
   ans <- ripser_cpp_dist(dataset, max_dim, threshold, 1., p)
   
   # coerce to 'PHom' class
-  ans <- new_PHom(ripser_ans_to_df(ans))
+  ans <- switch(
+    match.arg(return_class),
+    PHom = {
+      lifecycle::deprecate_soft(
+        "1.1.0",
+        I("'PHom' class"),
+        with = I("'persistence' from the {phutil} package"),
+        id = "PHom"
+      )
+      new_PHom(ripser_ans_to_df(ans))
+    },
+    persistence = phutil::as_persistence(
+      ans,
+      engine = "ripserr::vietoris_rips",
+      filtration = "Vietoris-Rips",
+      parameters = list(max_dim = max_dim, threshold = threshold, p = p)
+    )
+  )
   
   # return
   return(ans)
 }
 
+#' @rdname vietoris_rips
 #' @aliases vietoris_rips.numeric vietoris_rips.ts
 #' @importFrom stats tsp
 #' @param data_dim desired end data dimension (for `"ts"`, defaults to obs/time
@@ -191,7 +250,6 @@ vietoris_rips.dist <- function(
 #' @param dim_lag time series lag factor between dimensions
 #' @param sample_lag time series lag factor between samples (rows)
 #' @param method currently only allows `"qa"` (quasi-attractor method)
-#' @rdname vietoris_rips
 #' @export vietoris_rips.numeric
 #' @export
 vietoris_rips.numeric <- function(
diff --git a/man/cubical.Rd b/man/cubical.Rd
index b2f1ea9..3387961 100644
--- a/man/cubical.Rd
+++ b/man/cubical.Rd
@@ -4,13 +4,23 @@
 \alias{cubical}
 \alias{cubical.array}
 \alias{cubical.matrix}
+\alias{cubical.numeric}
 \title{Calculating Persistent Homology via a Cubical Complex}
 \usage{
 cubical(dataset, ...)
 
-\method{cubical}{array}(dataset, threshold = 9999, method = "lj", sublevel = TRUE, ...)
+\method{cubical}{array}(
+  dataset,
+  threshold = 9999,
+  method = "lj",
+  sublevel = TRUE,
+  return_class = c("PHom", "persistence"),
+  ...
+)
 
 \method{cubical}{matrix}(dataset, ...)
+
+\method{cubical}{numeric}(dataset, ...)
 }
 \arguments{
 \item{dataset}{object on which to calculate persistent homology}
@@ -24,9 +34,13 @@ see Kaji et al. (2020) \url{https://arxiv.org/abs/2005.12692} for details}
 
 \item{sublevel}{logical; whether to take the sublevel set filtration or else
 the superlevel set filtration}
+
+\item{return_class}{class of output object; either \code{"PHom"} (default; legacy)
+or \code{"persistence"} (from the
+\strong{\href{https://cran.r-project.org/package=phutil}{phutil}} package)}
 }
 \value{
-\code{PHom} object
+\code{"PHom"} or \code{\link[phutil:persistence]{"persistence"}} object
 }
 \description{
 This function is an R wrapper for the CubicalRipser C++ library to calculate
@@ -50,19 +64,25 @@ persistence diagram in which \code{death} precedes \code{birth}.
 }
 \examples{
 
+# 1-dim example
+dataset <- rnorm(1500)
+dim(dataset) <- 1500
+cubical_hom1 <- cubical(dataset)
+
 # 2-dim example
 dataset <- rnorm(10 ^ 2)
 dim(dataset) <- rep(10, 2)
-cubical_hom2 <- cubical(dataset)
+( cubical_hom2 <- cubical(dataset) )
 
 # 3-dim example
 dataset <- rnorm(8 ^ 3)
 dim(dataset) <- rep(8, 3)
-cubical_hom3 <- cubical(dataset)
+( cubical_hom3 <- cubical(dataset) )
 
 # 4-dim example
 dataset <- rnorm(5 ^ 4)
 dim(dataset) <- rep(5, 4)
+( cubical_hom4 <- cubical(dataset) )
 
 # sublevel versus superlevel
 cubical(volcano)
diff --git a/man/vietoris_rips.Rd b/man/vietoris_rips.Rd
index 0b3ac4a..664d7d3 100644
--- a/man/vietoris_rips.Rd
+++ b/man/vietoris_rips.Rd
@@ -15,9 +15,25 @@ vietoris_rips(dataset, ...)
 
 \method{vietoris_rips}{data.frame}(dataset, ...)
 
-\method{vietoris_rips}{matrix}(dataset, max_dim = 1L, threshold = -1, p = 2L, dim = NULL, ...)
+\method{vietoris_rips}{matrix}(
+  dataset,
+  max_dim = 1L,
+  threshold = -1,
+  p = 2L,
+  dim = NULL,
+  return_class = c("PHom", "persistence"),
+  ...
+)
 
-\method{vietoris_rips}{dist}(dataset, max_dim = 1L, threshold = -1, p = 2L, dim = NULL, ...)
+\method{vietoris_rips}{dist}(
+  dataset,
+  max_dim = 1L,
+  threshold = -1,
+  p = 2L,
+  dim = NULL,
+  return_class = c("PHom", "persistence"),
+  ...
+)
 
 \method{vietoris_rips}{numeric}(
   dataset,
@@ -49,6 +65,10 @@ calculated}
 \item{dim}{deprecated; passed to \code{max_dim} or ignored if \code{max_dim} is
 specified}
 
+\item{return_class}{class of output object; either \code{"PHom"} (default; legacy)
+or \code{"persistence"} (from the
+\strong{\href{https://cran.r-project.org/package=phutil}{phutil}} package)}
+
 \item{data_dim}{desired end data dimension (for \code{"ts"}, defaults to obs/time
 if > 1)}
 
@@ -59,7 +79,7 @@ if > 1)}
 \item{method}{currently only allows \code{"qa"} (quasi-attractor method)}
 }
 \value{
-\code{PHom} object
+\code{"PHom"} or \code{\link[phutil:persistence]{"persistence"}} object
 }
 \description{
 This function is an R wrapper for the Ripser C++ library to
@@ -97,7 +117,7 @@ rand.angle <- runif(num.pts, 0, 2*pi)
 pt.cloud <- cbind(cos(rand.angle), sin(rand.angle))
 
 # calculate persistent homology (num.pts by 3 numeric matrix)
-pers.hom <- vietoris_rips(pt.cloud)
+( pers.hom <- vietoris_rips(pt.cloud) )
 
 # sliding window persistent homology
 ( ld.phom <- vietoris_rips(ldeaths, data_dim = 12) )
diff --git a/tests/testthat/test-cubical-2dim.R b/tests/testthat/test-cubical-2dim.R
index 72cb044..b175d90 100644
--- a/tests/testthat/test-cubical-2dim.R
+++ b/tests/testthat/test-cubical-2dim.R
@@ -54,4 +54,14 @@ test_that("2-dim cubical returns same values as validated tests", {
   # check means of births and deaths to ensure close enough
   expect_equal(mean(test_output$birth), mean(output_data$birth), tolerance = 0.025)
   expect_equal(mean(test_output$death), mean(output_data$death), tolerance = 0.025)
-})
\ No newline at end of file
+})
+
+test_that("specified class is returned", {
+  # calculate 'PHom' object
+  expect_s3_class(cubical(test_data, return_class = "PHom"),
+                  "PHom")
+  
+  # calculate 'persistence' object
+  expect_s3_class(cubical(test_data, return_class = "persistence"),
+                  "persistence")
+})
diff --git a/tests/testthat/test-cubical-3dim.R b/tests/testthat/test-cubical-3dim.R
index cd75677..7ab33b7 100644
--- a/tests/testthat/test-cubical-3dim.R
+++ b/tests/testthat/test-cubical-3dim.R
@@ -1,14 +1,14 @@
 context("cubical 3-dim")
 library("ripserr")
 
+# reproducibility
+set.seed(42)
+
+# create data
+test_data <- rnorm(10 ^ 3)
+dim(test_data) <- rep(10, 3)
+
 test_that("basic 3-dim cubical works", {
-  # reproducibility
-  set.seed(42)
-  
-  # create data
-  test_data <- rnorm(10 ^ 3)
-  dim(test_data) <- rep(10, 3)
-  
   # create cubical complex
   cub_comp <- ripserr::cubical(test_data)
   
@@ -54,4 +54,14 @@ test_that("3-dim calculation returns same values as validated tests", {
   # check means of births and deaths to ensure close enough
   expect_equal(mean(test_output$birth), mean(output_data$birth))
   expect_equal(mean(test_output$death), mean(output_data$death))
-})
\ No newline at end of file
+})
+
+test_that("specified class is returned", {
+  # calculate 'PHom' object
+  expect_s3_class(cubical(test_data, return_class = "PHom"),
+                  "PHom")
+  
+  # calculate 'persistence' object
+  expect_s3_class(cubical(test_data, return_class = "persistence"),
+                  "persistence")
+})
diff --git a/tests/testthat/test-cubical-4dim.R b/tests/testthat/test-cubical-4dim.R
index 568d368..58ecfee 100644
--- a/tests/testthat/test-cubical-4dim.R
+++ b/tests/testthat/test-cubical-4dim.R
@@ -1,14 +1,14 @@
 context("cubical 4-dim")
 library("ripserr")
 
+# reproducibility
+set.seed(42)
+
+# create dataset to be used
+tester <- rnorm(5 ^ 4)
+dim(tester) <- rep(5, 4)
+
 test_that("basic 4-dim cubical works", {
-  # reproducibility
-  set.seed(42)
-  
-  # create dataset to be used
-  tester <- rnorm(5 ^ 4)
-  dim(tester) <- rep(5, 4)
-  
   # calculate cubical complex for 4-dim voxel data
   cub_comp <- ripserr::cubical(tester)
   
@@ -53,4 +53,14 @@ test_that("4-dim calculation returns same values as validated tests", {
   # check means of births and deaths to ensure close enough
   expect_equal(mean(test_output$birth), mean(output_data$birth), tolerance = 0.025)
   expect_equal(mean(test_output$death), mean(output_data$death), tolerance = 0.025)
-})
\ No newline at end of file
+})
+
+test_that("specified class is returned", {
+  # calculate 'PHom' object
+  expect_s3_class(cubical(tester, return_class = "PHom"),
+                  "PHom")
+  
+  # calculate 'persistence' object
+  expect_s3_class(cubical(tester, return_class = "persistence"),
+                  "persistence")
+})
diff --git a/tests/testthat/test-vr.R b/tests/testthat/test-vr.R
index e920ab4..525d96b 100644
--- a/tests/testthat/test-vr.R
+++ b/tests/testthat/test-vr.R
@@ -69,3 +69,21 @@ test_that("consistency across generic methods for time series", {
   # compare persistent homology across classes
   expect_equal(num_phom, ts_phom)
 })
+
+test_that("specified class is returned", {
+  # calculate 'PHom' object for each class
+  expect_s3_class(vietoris_rips(circle_mat, return_class = "PHom"),
+                  "PHom")
+  expect_s3_class(vietoris_rips(circle_df, return_class = "PHom"),
+                  "PHom")
+  expect_s3_class(vietoris_rips(circle_dist, return_class = "PHom"),
+                  "PHom")
+  
+  # calculate 'persistence' object for each class
+  expect_s3_class(vietoris_rips(circle_mat, return_class = "persistence"),
+                  "persistence")
+  expect_s3_class(vietoris_rips(circle_df, return_class = "persistence"),
+                  "persistence")
+  expect_s3_class(vietoris_rips(circle_dist, return_class = "persistence"),
+                  "persistence")
+})