diff --git a/.Rbuildignore b/.Rbuildignore
index 1ff3218..b1b9b14 100644
--- a/.Rbuildignore
+++ b/.Rbuildignore
@@ -14,6 +14,7 @@
 ^vignettes/test$
 ^README\.Rmd$
 ^attic$
+^\.superpowers$
 ^doc$
 ^Meta$
 ^.vscode$
diff --git a/DESCRIPTION b/DESCRIPTION
index 6755492..e97a594 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,6 +1,6 @@
 Package: randomPlantedForest
 Title: Random Planted Forest: A Directly Interpretable Tree Ensemble
-Version: 0.3.0
+Version: 0.3.0.9000
 Authors@R: c(
     person(c("Joseph", "Theo"), "Meyer", role = "aut"),
     person("Munir", "Hiabu", role = "aut"),
@@ -24,6 +24,7 @@ Imports:
     Rcpp,
     utils
 Suggests:
+    bundle,
     knitr,
     mvtnorm,
     recipes,
diff --git a/NAMESPACE b/NAMESPACE
index dbf99cd..b62174b 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -1,5 +1,6 @@
 # Generated by roxygen2: do not edit by hand
 
+S3method(bundle::bundle,rpf)
 S3method(predict,rpf)
 S3method(print,rpf)
 S3method(print,rpf_forest)
@@ -16,6 +17,9 @@ export(predict_components)
 export(preprocess_predictors_predict)
 export(purify)
 export(rpf)
+export(rpf_is_valid)
+export(rpf_marshal)
+export(rpf_unmarshal)
 import(checkmate)
 importFrom(Rcpp,loadModule)
 importFrom(Rcpp,sourceCpp)
diff --git a/NEWS.md b/NEWS.md
index 5427b22..29a971f 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -1,3 +1,17 @@
+# randomPlantedForest 0.3.0.9000
+
+* New `rpf_marshal()` / `rpf_unmarshal()` serialize a fitted forest to a plain
+  R list and back, making `saveRDS()`-based storage of rpf models possible (#52).
+  Purified forests restore their purified state directly; training data is only
+  embedded with `include_data = TRUE` (required to `purify()` after restoring).
+  Blobs record the blob format version and the package version they were
+  created with; restoring under an older package than the one that saved the
+  model warns. Malformed or corrupt blobs error instead of reading out of bounds.
+* New `rpf_is_valid()` checks whether an rpf object's internal model is usable;
+  `predict()`, `purify()` and `predict_components()` now give an actionable
+  error for rpf objects restored via `readRDS()` without marshaling.
+* New `bundle::bundle()` method for rpf models wrapping the marshaling API.
+
 # randomPlantedForest 0.3.0
 
 ## Major changes (#61)
diff --git a/R/bundle.R b/R/bundle.R
new file mode 100644
index 0000000..23f8111
--- /dev/null
+++ b/R/bundle.R
@@ -0,0 +1,25 @@
+#' Bundle an rpf model
+#'
+#' Method for [bundle::bundle()] wrapping [rpf_marshal()]/[rpf_unmarshal()],
+#' so rpf models work with the standard tidymodels serialization workflow.
+#' Training data is not included; see [rpf_marshal()] for the implications.
+#'
+#' @param x An [rpf] object.
+#' @param ... Unused.
+#' @return An object of class `bundled_rpf` for `bundle()`, or the restored
+#'   [rpf] object for `unbundle()`.
+#' @examplesIf requireNamespace("bundle", quietly = TRUE)
+#' fit <- rpf(mpg ~ wt + cyl, data = mtcars, ntrees = 10)
+#' b <- bundle::bundle(fit)
+#' tmp <- tempfile(fileext = ".rds")
+#' saveRDS(b, tmp)
+#' restored <- bundle::unbundle(readRDS(tmp))
+#' predict(restored, mtcars)
+#' @exportS3Method bundle::bundle
+bundle.rpf <- function(x, ...) {
+  bundle::bundle_constr(
+    object = rpf_marshal(x),
+    situate = bundle::situate_constr(function(object) randomPlantedForest::rpf_unmarshal(object)),
+    desc_class = "rpf"
+  )
+}
diff --git a/R/marshal.R b/R/marshal.R
new file mode 100644
index 0000000..9f943a4
--- /dev/null
+++ b/R/marshal.R
@@ -0,0 +1,126 @@
+#' Serialize and restore Random Planted Forests
+#'
+#' `rpf_marshal()` converts an [rpf] object into a plain R list safe for
+#' [saveRDS()]; `rpf_unmarshal()` reverses it. The C++ forest behind an `rpf`
+#' object does not survive R serialization, so this explicit round-trip is
+#' required to store or transfer fitted models.
+#'
+#' Training data is only included with `include_data = TRUE`; without it, a
+#' restored forest can predict (including purified prediction if the forest
+#' was purified before marshaling) but can never be purified afterwards.
+#'
+#' @param x An object of class `rpf`.
+#' @param include_data `[FALSE]`: Store training data in the blob, enabling
+#'   [purify()] after restoring.
+#' @param blob An object of class `rpf_marshaled` created by `rpf_marshal()`.
+#' @return `rpf_marshal()` returns a list of class `rpf_marshaled`;
+#'   `rpf_unmarshal()` returns an object of class [rpf].
+#' @export
+#' @examples
+#' fit <- rpf(mpg ~ wt + cyl, data = mtcars, ntrees = 10)
+#' blob <- rpf_marshal(fit)
+#' tmp <- tempfile(fileext = ".rds")
+#' saveRDS(blob, tmp)
+#' restored <- rpf_unmarshal(readRDS(tmp))
+#' all.equal(predict(fit, mtcars), predict(restored, mtcars))
+rpf_marshal <- function(x, include_data = FALSE) {
+  checkmate::assert_class(x, "rpf")
+  checkmate::assert_flag(include_data)
+  check_rpf_alive(x)
+
+  purified <- x$fit$is_purified()
+  fit_state <- list(
+    version = 1L,
+    pkg_version = utils::packageVersion("randomPlantedForest"),
+    model = x$fit$get_model(),
+    bounds = x$fit$get_bounds(),
+    dims = x$fit$get_shape(),
+    purified = purified,
+    grid = if (purified) x$fit$get_grid_leaves(),
+    data = if (include_data) x$fit$get_data(),
+    # $forest (export_forest = TRUE) duplicates $model; rebuilt on unmarshal
+    had_forest = !is.null(x$forest)
+  )
+
+  out <- unclass(x)
+  out$fit <- NULL
+  out$forest <- NULL
+  out$fit_state <- fit_state
+  structure(out, class = "rpf_marshaled")
+}
+
+#' @rdname rpf_marshal
+#' @export
+rpf_unmarshal <- function(blob) {
+  checkmate::assert_class(blob, "rpf_marshaled")
+  state <- blob$fit_state
+  if (!identical(state$version, 1L)) {
+    stop("Unsupported rpf_marshaled version: ", state$version)
+  }
+  installed <- utils::packageVersion("randomPlantedForest")
+  if (!is.null(state$pkg_version) && state$pkg_version > installed) {
+    warning(
+      "This model was marshaled with randomPlantedForest ", state$pkg_version,
+      " but version ", installed, " is installed. Restoring may not be reliable.",
+      call. = FALSE
+    )
+  }
+
+  pars <- rpf_param_vector(blob$params, blob$mode)
+  fit <- if (blob$mode == "classification") {
+    methods::new(ClassificationRPF, blob$params$loss, pars)
+  } else {
+    methods::new(RandomPlantedForest, pars)
+  }
+  fit$set_shape(
+    state$dims$feature_size,
+    state$dims$value_size,
+    state$dims$sample_size,
+    state$bounds$lower,
+    state$bounds$upper
+  )
+  if (!is.null(state$data)) {
+    fit$set_training_data(as.matrix(state$data$Y), as.matrix(state$data$X))
+  }
+  fit$set_model(state$model)
+  if (state$purified) {
+    fit$set_grid_leaves(state$grid)
+  }
+
+  forest <- NULL
+  if (isTRUE(state$had_forest)) {
+    forest <- fit$get_model()
+    class(forest) <- "rpf_forest"
+  }
+
+  out <- blob
+  out$fit_state <- NULL
+  out <- unclass(out)
+  fields <- out[setdiff(names(out), c("fit", "blueprint", "forest"))]
+  do.call(new_rpf, c(list(fit = fit, blueprint = blob$blueprint, forest = forest), fields))
+}
+
+#' Check whether an rpf object's C++ forest is still alive
+#'
+#' The C++ forest does not survive [saveRDS()]; an `rpf` object restored via
+#' [readRDS()] without [rpf_marshal()]/[rpf_unmarshal()] is unusable.
+#' @param x An object of class `rpf`.
+#' @return `TRUE` if the underlying model can be used, `FALSE` otherwise.
+#' @export
+rpf_is_valid <- function(x) {
+  checkmate::assert_class(x, "rpf")
+  !inherits(try(x$fit$is_purified(), silent = TRUE), "try-error")
+}
+
+check_rpf_alive <- function(x) {
+  if (!rpf_is_valid(x)) {
+    stop(
+      "The C++ forest behind this rpf object is gone - most likely it was ",
+      "saved with saveRDS() and restored with readRDS().\n",
+      "Use blob <- rpf_marshal(x) before saving and rpf_unmarshal(blob) ",
+      "after loading. See ?rpf_marshal.",
+      call. = FALSE
+    )
+  }
+  invisible(x)
+}
diff --git a/R/predict_components.R b/R/predict_components.R
index 6dc9910..bf6be47 100644
--- a/R/predict_components.R
+++ b/R/predict_components.R
@@ -70,6 +70,7 @@
 #'
 predict_components <- function(object, new_data, max_interaction = NULL, predictors = NULL) {
   checkmate::assert_class(object, classes = "rpf")
+  check_rpf_alive(object)
 
   # if max_interaction is not provided, we use the max from the rpf model fit
   if (is.null(max_interaction)) {
diff --git a/R/predict_rpf.R b/R/predict_rpf.R
index f42c012..54f2988 100644
--- a/R/predict_rpf.R
+++ b/R/predict_rpf.R
@@ -30,6 +30,7 @@
 #' rpfit <- rpf(y = mtcars$mpg, x = mtcars[, c("cyl", "wt")])
 #' predict(rpfit, mtcars[, c("cyl", "wt")])
 predict.rpf <- function(object, new_data, type = ifelse(object$mode == "regression", "numeric", "prob"), ...) {
+  check_rpf_alive(object)
   # Enforces column order, type, column names, etc
   processed <- hardhat::forge(new_data, object$blueprint)
 
diff --git a/R/purify.R b/R/purify.R
index 5d92216..d752d5b 100644
--- a/R/purify.R
+++ b/R/purify.R
@@ -43,6 +43,7 @@ purify.default <- function(x, ...) {
 #' @importFrom utils capture.output
 purify.rpf <- function(x, ..., maxp_interaction = NULL, mode = 2L, nthreads = NULL) {
   checkmate::assert_class(x, "rpf")
+  check_rpf_alive(x)
   checkmate::assert_int(mode, lower = 1, upper = 2)
   if (!is.null(nthreads)) {
     checkmate::assert_int(nthreads, lower = 1)
@@ -66,5 +67,6 @@ purify.rpf <- function(x, ..., maxp_interaction = NULL, mode = 2L, nthreads = NU
 #' @rdname purify
 is_purified <- function(x) {
   checkmate::assert_class(x, "rpf")
+  check_rpf_alive(x)
   x$fit$is_purified()
 }
diff --git a/R/rpf.R b/R/rpf.R
index f1aec4b..d02d13a 100644
--- a/R/rpf.R
+++ b/R/rpf.R
@@ -371,26 +371,30 @@ rpf_bridge <- function(
   checkmate::assert_flag(delete_leaves)
   checkmate::assert_choice(split_structure, choices = c("res_trees", "cur_trees_2", "cur_trees_1", "leaves", "hist"))
 
-  fit <- rpf_impl(
-    Y = outcomes$outcomes,
-    X = predictors$predictors_matrix,
-    mode = outcomes$mode,
-    max_interaction = max_interaction,
+  params <- list(
+    loss = loss,
     ntrees = ntrees,
+    max_interaction = max_interaction,
     splits = splits,
     split_try = split_try,
     t_try = t_try,
     split_decay_rate = split_decay_rate,
     max_candidates = max_candidates,
     delete_leaves = delete_leaves,
+    split_structure = split_structure,
+    delta = delta,
+    epsilon = epsilon,
     deterministic = deterministic,
     nthreads = nthreads,
     purify = purify,
-    cv = cv,
-    loss = loss,
-    delta = delta,
-    epsilon = epsilon,
-    split_structure = split_structure
+    cv = cv
+  )
+
+  fit <- rpf_impl(
+    Y = outcomes$outcomes,
+    X = predictors$predictors_matrix,
+    mode = outcomes$mode,
+    params = params
   )
 
   # Optionally export a compact R list representation of the forest.
@@ -405,24 +409,7 @@ rpf_bridge <- function(
     blueprint = processed$blueprint,
     mode = outcomes$mode,
     factor_levels = predictors$factor_levels,
-    params = list(
-      loss = loss, # FIXME: Dedup, requires changes in tests and predict
-      ntrees = ntrees,
-      max_interaction = max_interaction,
-      splits = splits,
-      split_try = split_try,
-      t_try = t_try,
-      split_decay_rate = split_decay_rate,
-      max_candidates = max_candidates,
-      delete_leaves = delete_leaves,
-      split_structure = split_structure,
-      delta = delta,
-      epsilon = epsilon,
-      deterministic = deterministic,
-      nthreads = nthreads,
-      purify = purify,
-      cv = cv
-    ),
+    params = params,
     forest = forest
   )
 }
@@ -437,80 +424,49 @@ new_rpf <- function(fit, blueprint, ...) {
   )
 }
 
+# Assemble the positional parameter vector for the C++ constructors.
+# 13 elements for regression, 15 (+ delta, epsilon) for classification.
+rpf_param_vector <- function(params, mode) {
+  split_mode <- switch(
+    params$split_structure,
+    res_trees = 0L,
+    cur_trees_2 = 1L,
+    cur_trees_1 = 2L,
+    leaves = 3L,
+    hist = 4L
+  )
+  base <- c(
+    params$max_interaction,
+    params$ntrees,
+    params$splits,
+    params$split_try,
+    params$t_try,
+    params$purify,
+    params$deterministic,
+    params$nthreads,
+    params$cv,
+    params$split_decay_rate,
+    params$max_candidates,
+    params$delete_leaves,
+    split_mode
+  )
+  if (mode == "classification") {
+    base <- c(base, params$delta, params$epsilon)
+  }
+  base
+}
+
 # Main fitting function and interface to C++ implementation
-rpf_impl <- function(
-  Y,
-  X,
-  mode = c("regression", "classification"),
-  max_interaction = 1,
-  ntrees = 50,
-  splits = 30,
-  split_try = 10,
-  t_try = 0.4,
-  deterministic = FALSE,
-  nthreads = 1,
-  purify = FALSE,
-  cv = FALSE,
-  split_decay_rate = 0.1,
-  max_candidates = 50,
-  delete_leaves = TRUE,
-  loss = "L2",
-  delta = 0,
-  epsilon = 0.1,
-  split_structure = "leaves"
-) {
+rpf_impl <- function(Y, X, mode = c("regression", "classification"), params) {
   # Final input validation, should be superfluous
   checkmate::assert_matrix(X, mode = "numeric", any.missing = FALSE)
   mode <- match.arg(mode)
-  split_structure <- match.arg(split_structure, c("leaves", "res_trees", "cur_trees_2", "cur_trees_1", "hist"))
-  # map split_structure string to numeric mode for C++
-  split_mode <- switch(split_structure, res_trees = 0L, cur_trees_2 = 1L, cur_trees_1 = 2L, leaves = 3L, hist = 4L)
+  pars <- rpf_param_vector(params, mode)
 
   if (mode == "classification") {
-    fit <- new(
-      ClassificationRPF,
-      Y,
-      X,
-      loss,
-      c(
-        max_interaction,
-        ntrees,
-        splits,
-        split_try,
-        t_try,
-        purify,
-        deterministic,
-        nthreads,
-        cv,
-        split_decay_rate,
-        max_candidates,
-        delete_leaves,
-        split_mode,
-        delta,
-        epsilon
-      )
-    )
-  } else if (mode == "regression") {
-    fit <- new(
-      RandomPlantedForest,
-      Y,
-      X,
-      c(
-        max_interaction,
-        ntrees,
-        splits,
-        split_try,
-        t_try,
-        purify,
-        deterministic,
-        nthreads,
-        cv,
-        split_decay_rate,
-        max_candidates,
-        delete_leaves,
-        split_mode
-      )
-    )
+    fit <- new(ClassificationRPF, Y, X, params$loss, pars)
+  } else {
+    fit <- new(RandomPlantedForest, Y, X, pars)
   }
 
   fit
diff --git a/_pkgdown.yml b/_pkgdown.yml
index d642cc8..790f387 100644
--- a/_pkgdown.yml
+++ b/_pkgdown.yml
@@ -30,4 +30,7 @@ reference:
   - starts_with("print.")
   - starts_with("str.")
   - preprocess_predictors_predict
+  - rpf_marshal
+  - rpf_is_valid
+  - bundle.rpf
 
diff --git a/man/bundle.rpf.Rd b/man/bundle.rpf.Rd
new file mode 100644
index 0000000..4d8f634
--- /dev/null
+++ b/man/bundle.rpf.Rd
@@ -0,0 +1,32 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/bundle.R
+\name{bundle.rpf}
+\alias{bundle.rpf}
+\title{Bundle an rpf model}
+\usage{
+\method{bundle}{rpf}(x, ...)
+}
+\arguments{
+\item{x}{An \link{rpf} object.}
+
+\item{...}{Unused.}
+}
+\value{
+An object of class \code{bundled_rpf} for \code{bundle()}, or the restored
+\link{rpf} object for \code{unbundle()}.
+}
+\description{
+Method for \code{\link[bundle:bundle]{bundle::bundle()}} wrapping \code{\link[=rpf_marshal]{rpf_marshal()}}/\code{\link[=rpf_unmarshal]{rpf_unmarshal()}},
+so rpf models work with the standard tidymodels serialization workflow.
+Training data is not included; see \code{\link[=rpf_marshal]{rpf_marshal()}} for the implications.
+}
+\examples{
+\dontshow{if (requireNamespace("bundle", quietly = TRUE)) withAutoprint(\{ # examplesIf}
+fit <- rpf(mpg ~ wt + cyl, data = mtcars, ntrees = 10)
+b <- bundle::bundle(fit)
+tmp <- tempfile(fileext = ".rds")
+saveRDS(b, tmp)
+restored <- bundle::unbundle(readRDS(tmp))
+predict(restored, mtcars)
+\dontshow{\}) # examplesIf}
+}
diff --git a/man/rpf_is_valid.Rd b/man/rpf_is_valid.Rd
new file mode 100644
index 0000000..284eb7d
--- /dev/null
+++ b/man/rpf_is_valid.Rd
@@ -0,0 +1,18 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/marshal.R
+\name{rpf_is_valid}
+\alias{rpf_is_valid}
+\title{Check whether an rpf object's C++ forest is still alive}
+\usage{
+rpf_is_valid(x)
+}
+\arguments{
+\item{x}{An object of class \code{rpf}.}
+}
+\value{
+\code{TRUE} if the underlying model can be used, \code{FALSE} otherwise.
+}
+\description{
+The C++ forest does not survive \code{\link[=saveRDS]{saveRDS()}}; an \code{rpf} object restored via
+\code{\link[=readRDS]{readRDS()}} without \code{\link[=rpf_marshal]{rpf_marshal()}}/\code{\link[=rpf_unmarshal]{rpf_unmarshal()}} is unusable.
+}
diff --git a/man/rpf_marshal.Rd b/man/rpf_marshal.Rd
new file mode 100644
index 0000000..f7bae4c
--- /dev/null
+++ b/man/rpf_marshal.Rd
@@ -0,0 +1,42 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/marshal.R
+\name{rpf_marshal}
+\alias{rpf_marshal}
+\alias{rpf_unmarshal}
+\title{Serialize and restore Random Planted Forests}
+\usage{
+rpf_marshal(x, include_data = FALSE)
+
+rpf_unmarshal(blob)
+}
+\arguments{
+\item{x}{An object of class \code{rpf}.}
+
+\item{include_data}{\verb{[FALSE]}: Store training data in the blob, enabling
+\code{\link[=purify]{purify()}} after restoring.}
+
+\item{blob}{An object of class \code{rpf_marshaled} created by \code{rpf_marshal()}.}
+}
+\value{
+\code{rpf_marshal()} returns a list of class \code{rpf_marshaled};
+\code{rpf_unmarshal()} returns an object of class \link{rpf}.
+}
+\description{
+\code{rpf_marshal()} converts an \link{rpf} object into a plain R list safe for
+\code{\link[=saveRDS]{saveRDS()}}; \code{rpf_unmarshal()} reverses it. The C++ forest behind an \code{rpf}
+object does not survive R serialization, so this explicit round-trip is
+required to store or transfer fitted models.
+}
+\details{
+Training data is only included with \code{include_data = TRUE}; without it, a
+restored forest can predict (including purified prediction if the forest
+was purified before marshaling) but can never be purified afterwards.
+}
+\examples{
+fit <- rpf(mpg ~ wt + cyl, data = mtcars, ntrees = 10)
+blob <- rpf_marshal(fit)
+tmp <- tempfile(fileext = ".rds")
+saveRDS(blob, tmp)
+restored <- rpf_unmarshal(readRDS(tmp))
+all.equal(predict(fit, mtcars), predict(restored, mtcars))
+}
diff --git a/src/include/cpf.hpp b/src/include/cpf.hpp
index f993f95..900b68a 100644
--- a/src/include/cpf.hpp
+++ b/src/include/cpf.hpp
@@ -10,9 +10,17 @@ class ClassificationRPF : public RandomPlantedForest
   using RandomPlantedForest::calcOptimalSplit;
   ClassificationRPF(const NumericMatrix &samples_Y, const NumericMatrix &samples_X,
                     const String loss = "L2", const NumericVector parameters = {1, 50, 30, 10, 0.4, 0, 0, 0, 0, 0.1, 0, 0.1, 50,1});
+  // Params-only constructor: parses configuration and loss but loads no data
+  // and does not fit. Used by rpf_unmarshal() to rebuild a serialized forest.
+  ClassificationRPF(const String loss, const NumericVector parameters);
   void set_parameters(StringVector keys, NumericVector values);
   ~ClassificationRPF(){};
 
+protected:
+  // Maps a loss-name string to `loss` and `calcLoss` (verbatim move of the
+  // former inline if/else chain).
+  void set_loss_function(const String &loss);
+
 private:
   double delta;
   double epsilon;
diff --git a/src/include/helper.hpp b/src/include/helper.hpp
index c80cbfd..3821451 100644
--- a/src/include/helper.hpp
+++ b/src/include/helper.hpp
@@ -83,6 +83,9 @@ namespace utils
         throw std::invalid_argument("Index out of range.");
       return entries[index];
     }
+    // Flat storage access for serialization (column-major, first index fastest)
+    std::vector<T> &flat() { return entries; }
+    const std::vector<T> &flat() const { return entries; }
     std::vector<int> dims;
     size_t n_entries = 1;
 
diff --git a/src/include/rpf.hpp b/src/include/rpf.hpp
index ff28a93..22c2d9d 100644
--- a/src/include/rpf.hpp
+++ b/src/include/rpf.hpp
@@ -35,6 +35,25 @@ class RandomPlantedForest
   RandomPlantedForest(const NumericMatrix &samples_Y, const NumericMatrix &samples_X,
                       const NumericVector parameters = {1, 50, 30, 10, 0.4, 0, 0, 0, 0, 0.1, 50, 1, 3});
   RandomPlantedForest(){};
+  // Params-only constructor: parses configuration but loads no data and does
+  // not fit. Used by rpf_unmarshal() to rebuild a serialized forest.
+  RandomPlantedForest(const NumericVector parameters);
+  // Restore shape metadata without training data (serialization path).
+  void set_shape(int feature_size_in, int value_size_in, int sample_size_in,
+                 const NumericVector lower, const NumericVector upper);
+  // Store training data without fitting and without recomputing bounds.
+  void set_training_data(const NumericMatrix &samples_Y, const NumericMatrix &samples_X);
+  // Export training data (for rpf_marshal(include_data = TRUE)).
+  List get_data();
+  // Export bounds/shape metadata (for rpf_marshal()).
+  List get_bounds();
+  List get_shape();
+  // Restore tree structure from get_model()'s output (serialization path).
+  void set_model(List &model);
+  // Export per-tree GridLeaves (purified leaf grid) for serialization.
+  List get_grid_leaves();
+  // Restore per-tree GridLeaves from get_grid_leaves()'s output; sets purified = true.
+  void set_grid_leaves(List &grid);
   // Load or replace data without fitting; computes bounds and resets state.
   void set_data(const NumericMatrix &samples_Y, const NumericMatrix &samples_X);
   // Predict for a matrix or a single vector. `components = {0}` means the full
@@ -63,6 +82,8 @@ class RandomPlantedForest
   bool is_purified();
   
 protected:
+  // Parse the flat parameter vector shared by both constructors.
+  void parse_parameters(const std::vector<double> &pars);
   // Internal per-family worker (grid-based mode 1)
   void purify_3_family(TreeFamily &curr_family, int maxp_interaction);
   // Internal per-family worker for fast exact purifier (mode 2)
diff --git a/src/lib/cpf.cpp b/src/lib/cpf.cpp
index 0bebd2a..f73d6ad 100644
--- a/src/lib/cpf.cpp
+++ b/src/lib/cpf.cpp
@@ -34,22 +34,10 @@ using namespace rpf_utils;
 
 // loss moved to lib/losses_*.cpp
 
-// constructor with parameters split_try, t_try, purify_forest, deterministic, nthreads
-ClassificationRPF::ClassificationRPF(const NumericMatrix &samples_Y, const NumericMatrix &samples_X,
-                                     const String loss, const NumericVector parameters)
-    : RandomPlantedForest( 
-        samples_Y, 
-        samples_X, 
-        // pass first 13 parameters to base (includes split_structure)
-        parameters.size() >= 13 ? parameters[Rcpp::Range(0, 12)] : parameters[Rcpp::Range(0, 11)] 
-      )
+// Maps a loss-name string to `loss` and `calcLoss` (verbatim move of the
+// former inline if/else chain from the data-fitting constructor).
+void ClassificationRPF::set_loss_function(const String &loss)
 {
-
-  // Ensure correct Rcpp RNG state
-  Rcpp::RNGScope scope;
-
-  // initialize class members
-  std::vector<double> pars = to_std_vec(parameters);
   if (loss == "L1")
   {
     this->loss = LossType::L1;
@@ -106,6 +94,25 @@ ClassificationRPF::ClassificationRPF(const NumericMatrix &samples_Y, const Numer
     this->loss = LossType::L2;
     this->calcLoss = &ClassificationRPF::L2_loss;
   }
+}
+
+// constructor with parameters split_try, t_try, purify_forest, deterministic, nthreads
+ClassificationRPF::ClassificationRPF(const NumericMatrix &samples_Y, const NumericMatrix &samples_X,
+                                     const String loss, const NumericVector parameters)
+    : RandomPlantedForest(
+        samples_Y,
+        samples_X,
+        // pass first 13 parameters to base (includes split_structure)
+        parameters.size() >= 13 ? parameters[Rcpp::Range(0, 12)] : parameters[Rcpp::Range(0, 11)]
+      )
+{
+
+  // Ensure correct Rcpp RNG state
+  Rcpp::RNGScope scope;
+
+  // initialize class members
+  std::vector<double> pars = to_std_vec(parameters);
+  set_loss_function(loss);
   if (pars.size() != 15)
   {
     Rcout << "Wrong number of parameters - set to default." << std::endl;
@@ -147,6 +154,21 @@ ClassificationRPF::ClassificationRPF(const NumericMatrix &samples_Y, const Numer
   this->set_data(samples_Y, samples_X);
 }
 
+// Params-only constructor: parses configuration and loss but loads no data
+// and does not fit. Used by rpf_unmarshal() to rebuild a serialized forest.
+ClassificationRPF::ClassificationRPF(const String loss, const NumericVector parameters)
+    : RandomPlantedForest(
+          parameters.size() >= 13 ? NumericVector(parameters[Rcpp::Range(0, 12)])
+                                  : NumericVector(parameters[Rcpp::Range(0, 11)]))
+{
+  std::vector<double> pars = to_std_vec(parameters);
+  if (pars.size() != 15)
+    Rcpp::stop("ClassificationRPF requires 15 parameters, got %d", (int)pars.size());
+  set_loss_function(loss);
+  this->delta = pars[13];
+  this->epsilon = pars[14];
+}
+
 // Mode 1: cur_trees_2 (classification variant)
 Split ClassificationRPF::calcOptimalSplit_curTrees2(const std::vector<std::vector<double>> &Y, const std::vector<std::vector<double>> &X,
                                            std::vector<SplitCandidate> &possible_splits, TreeFamily &curr_family, std::vector<std::vector<double>> &weights)
diff --git a/src/lib/purify.cpp b/src/lib/purify.cpp
index db0a4da..f36a8cd 100644
--- a/src/lib/purify.cpp
+++ b/src/lib/purify.cpp
@@ -383,6 +383,10 @@ void RandomPlantedForest::purify_fast_exact_family(TreeFamily &curr_family, int
 // Unified purifier entry: mode 1 = grid path, mode 2 = fast exact path
 void RandomPlantedForest::purify(int maxp_interaction, int nthreads_param, int mode)
 {
+  if (sample_size == 0 || X.empty())
+    Rcpp::stop("Cannot purify: no training data available. If this forest was "
+               "restored with rpf_unmarshal(), marshal it with include_data = TRUE.");
+
   // Determine threads: if user provided >0, use it; otherwise default to
   // min(object-configured nthreads, hardware concurrency)
   unsigned int threads_to_use = 0;
diff --git a/src/lib/rpf.cpp b/src/lib/rpf.cpp
index 858d9dc..3c8453d 100644
--- a/src/lib/rpf.cpp
+++ b/src/lib/rpf.cpp
@@ -64,6 +64,24 @@ void RandomPlantedForest::L2_loss(Split &split)
   }
 }
 
+void RandomPlantedForest::parse_parameters(const std::vector<double> &pars)
+{
+  this->max_interaction = pars[0];
+  this->n_trees = pars[1];
+  this->n_splits = pars[2];
+  this->split_try = pars[3];
+  this->t_try = pars[4];
+  this->purify_forest = pars[5];
+  this->deterministic = pars[6];
+  this->nthreads = pars[7];
+  this->cross_validate = pars[8];
+  this->split_decay_rate_ = pars[9];
+  this->max_candidates_   = static_cast<size_t>(pars[10]);
+  this->delete_leaves   = (pars[11] != 0);
+  // map: 0=res_trees, 1=cur_trees_2, 2=cur_trees_1, 3=leaves, 4=hist
+  this->split_structure_mode_ = (pars.size() >= 13) ? static_cast<int>(pars[12]) : 3;
+}
+
 // constructor (parsing includes split_structure)
 RandomPlantedForest::RandomPlantedForest(const NumericMatrix &samples_Y, const NumericMatrix &samples_X,
                                          const NumericVector parameters)
@@ -76,24 +94,174 @@ RandomPlantedForest::RandomPlantedForest(const NumericMatrix &samples_Y, const N
   }
   else
   {
-    this->max_interaction = pars[0];
-    this->n_trees = pars[1];
-    this->n_splits = pars[2];
-    this->split_try = pars[3];
-    this->t_try = pars[4];
-    this->purify_forest = pars[5];
-    this->deterministic = pars[6];
-    this->nthreads = pars[7];
-    this->cross_validate = pars[8];
-    this->split_decay_rate_ = pars[9];
-    this->max_candidates_   = static_cast<size_t>(pars[10]);
-    this->delete_leaves   = (pars[11] != 0);
-    // map: 0=res_trees, 1=cur_trees_2, 2=cur_trees_1, 3=leaves, 4=hist
-    this->split_structure_mode_ = (pars.size() >= 13) ? static_cast<int>(pars[12]) : 3;
+    parse_parameters(pars);
   }
   this->set_data(samples_Y, samples_X);
 }
 
+// Params-only constructor: parses configuration but loads no data and does
+// not fit. Used by rpf_unmarshal() to rebuild a serialized forest.
+RandomPlantedForest::RandomPlantedForest(const NumericVector parameters)
+{
+  std::vector<double> pars = to_std_vec(parameters);
+  if (pars.size() != 12 && pars.size() != 13)
+    Rcpp::stop("RandomPlantedForest requires 12 or 13 parameters, got %d", pars.size());
+  parse_parameters(pars);
+}
+
+void RandomPlantedForest::set_shape(int feature_size_in, int value_size_in, int sample_size_in,
+                                    const NumericVector lower, const NumericVector upper)
+{
+  this->feature_size = feature_size_in;
+  this->value_size = value_size_in;
+  this->sample_size = sample_size_in;
+  this->lower_bounds = to_std_vec(lower);
+  this->upper_bounds = to_std_vec(upper);
+  this->n_leaves = std::vector<int>(feature_size, 1);
+}
+
+void RandomPlantedForest::set_training_data(const NumericMatrix &samples_Y, const NumericMatrix &samples_X)
+{
+  // Restore path only: bounds are NOT recomputed (they come from the blob via
+  // set_shape), but shapes must be consistent with the restored forest so the
+  // purify path cannot read out of bounds.
+  if (samples_X.ncol() != feature_size)
+    Rcpp::stop("Corrupt training data: X has %d columns, expected %d.",
+               samples_X.ncol(), feature_size);
+  if (samples_X.nrow() != samples_Y.nrow())
+    Rcpp::stop("Corrupt training data: X has %d rows but Y has %d.",
+               samples_X.nrow(), samples_Y.nrow());
+  this->Y = to_std_vec(samples_Y);
+  this->X = to_std_vec(samples_X);
+  this->sample_size = X.size();
+}
+
+List RandomPlantedForest::get_data()
+{
+  return List::create(Named("X") = from_std_vec(X), Named("Y") = from_std_vec(Y));
+}
+
+List RandomPlantedForest::get_bounds()
+{
+  return List::create(Named("lower") = from_std_vec(lower_bounds),
+                      Named("upper") = from_std_vec(upper_bounds));
+}
+
+List RandomPlantedForest::get_shape()
+{
+  return List::create(Named("feature_size") = feature_size,
+                      Named("value_size") = (int)value_size,
+                      Named("sample_size") = (int)sample_size);
+}
+
+void RandomPlantedForest::set_model(List &model)
+{
+  size_t n_families = model.size();
+  tree_families = std::vector<TreeFamily>(n_families);
+  for (size_t i = 0; i < n_families; ++i) {
+    List family = model[i];
+    List variables = family["variables"];
+    List values = family["values"];
+    List intervals = family["intervals"];
+    size_t n_trees_fam = variables.size();
+    for (size_t j = 0; j < n_trees_fam; ++j) {
+      IntegerVector tree_variables = variables[j];
+      std::set<int> dims(tree_variables.begin(), tree_variables.end());
+      List tree_values = values[j];
+      List tree_intervals = intervals[j];
+      size_t n_leaves = tree_values.size();
+      std::vector<Leaf> leaves(n_leaves);
+      for (size_t k = 0; k < n_leaves; ++k) {
+        // get_model() builds leaf_values via push_back() onto a default-constructed
+        // NumericMatrix, which yields a plain vector without matrix dims - read it
+        // back as a vector rather than casting to NumericMatrix.
+        leaves[k].value = as<std::vector<double>>(tree_values[k]);
+        NumericMatrix leaf_intervals = tree_intervals[k];
+        if (leaf_intervals.nrow() < 2 || leaf_intervals.ncol() < feature_size)
+          Rcpp::stop("Corrupt model data: leaf interval matrix has dimensions %dx%d, expected 2x%d.",
+                     leaf_intervals.nrow(), leaf_intervals.ncol(), feature_size);
+        std::vector<Interval> ivs(feature_size);
+        for (int l = 0; l < feature_size; ++l)
+          ivs[l] = Interval{leaf_intervals(0, l), leaf_intervals(1, l)};
+        leaves[k].intervals = ivs;
+      }
+      tree_families[i].insert(
+          std::make_pair(dims, std::make_shared<DecisionTree>(DecisionTree(dims, leaves))));
+    }
+  }
+  purified = false;
+}
+
+List RandomPlantedForest::get_grid_leaves()
+{
+  List families;
+  for (auto &family : tree_families) {
+    List trees;
+    List fam_lim_list;
+    bool lim_captured = false;
+    for (auto &tree : family) {
+      if (!lim_captured) {
+        // lim_list is identical for every tree in a family; store once
+        List ll;
+        for (auto &v : tree.second->GridLeaves.lim_list) ll.push_back(from_std_vec(v));
+        fam_lim_list = ll;
+        lim_captured = true;
+      }
+      auto &vals = tree.second->GridLeaves.values;
+      const auto &flat = vals.flat();
+      NumericMatrix vmat(flat.size(), value_size);
+      for (size_t e = 0; e < flat.size(); ++e)
+        for (size_t p = 0; p < value_size && p < flat[e].size(); ++p)
+          vmat(e, p) = flat[e][p];
+      trees.push_back(List::create(
+          Named("variables") = from_std_set(tree.first),
+          Named("dims") = from_std_vec(vals.dims),
+          Named("values") = vmat));
+    }
+    families.push_back(List::create(Named("lim_list") = fam_lim_list,
+                                    Named("trees") = trees));
+  }
+  return families;
+}
+
+void RandomPlantedForest::set_grid_leaves(List &grid)
+{
+  if ((size_t)grid.size() != tree_families.size())
+    Rcpp::stop("Grid data does not match the number of tree families.");
+  for (size_t i = 0; i < tree_families.size(); ++i) {
+    List family = grid[i];
+    List fam_lim_list = family["lim_list"];
+    std::vector<std::vector<double>> lim_list;
+    for (int l = 0; l < fam_lim_list.size(); ++l)
+      lim_list.push_back(to_std_vec(NumericVector(fam_lim_list[l])));
+    List trees = family["trees"];
+    for (int j = 0; j < trees.size(); ++j) {
+      List tr = trees[j];
+      IntegerVector vars = tr["variables"];
+      std::set<int> dims_set(vars.begin(), vars.end());
+      auto it = tree_families[i].find(dims_set);
+      if (it == tree_families[i].end())
+        Rcpp::stop("Grid data references a tree not present in the forest.");
+      std::vector<int> mdims = to_std_vec(IntegerVector(tr["dims"]));
+      for (int d : mdims)
+        if (d <= 0)
+          Rcpp::stop("Corrupt grid data: non-positive grid dimension %d.", d);
+      NumericMatrix vmat = tr["values"];
+      utils::Matrix<std::vector<double>> values(mdims, std::vector<double>(value_size, 0));
+      auto &flat = values.flat();
+      if ((size_t)vmat.nrow() != flat.size() || (size_t)vmat.ncol() != value_size)
+        Rcpp::stop("Corrupt grid data: values matrix has dimensions %dx%d, expected %dx%d.",
+                   vmat.nrow(), vmat.ncol(), (int)flat.size(), (int)value_size);
+      for (size_t e = 0; e < flat.size(); ++e)
+        for (size_t p = 0; p < value_size; ++p)
+          flat[e][p] = vmat(e, p);
+      it->second->GridLeaves.values = values;
+      it->second->GridLeaves.lim_list = lim_list;
+    }
+  }
+  purified = true;
+}
+
 // --------------- calcOptimalSplit per mode ---------------
 
 // Mode 3: leaves implementation moved to lib/splits_leaves.cpp
diff --git a/src/randomPlantedForest.cpp b/src/randomPlantedForest.cpp
index 4bb7b5a..772ba85 100644
--- a/src/randomPlantedForest.cpp
+++ b/src/randomPlantedForest.cpp
@@ -9,7 +9,14 @@ RCPP_MODULE(mod_rpf)
 
   class_<RandomPlantedForest>("RandomPlantedForest")
       .constructor<const NumericMatrix, const NumericMatrix, const NumericVector>()
+      .constructor<const NumericVector>()
       .method("set_data", &RandomPlantedForest::set_data)
+      .method("set_shape", &RandomPlantedForest::set_shape)
+      .method("set_training_data", &RandomPlantedForest::set_training_data)
+      .method("get_data", &RandomPlantedForest::get_data)
+      .method("get_bounds", &RandomPlantedForest::get_bounds)
+      .method("get_shape", &RandomPlantedForest::get_shape)
+      .method("set_model", &RandomPlantedForest::set_model)
       .method("get_parameters", &RandomPlantedForest::get_parameters)
       .method("cross_validation", &RandomPlantedForest::cross_validation)
       .method("predict_matrix", &RandomPlantedForest::predict_matrix)
@@ -19,10 +26,13 @@ RCPP_MODULE(mod_rpf)
       .method("print", &RandomPlantedForest::print)
       .method("set_parameters", &RandomPlantedForest::set_parameters)
       .method("get_model", &RandomPlantedForest::get_model)
+      .method("get_grid_leaves", &RandomPlantedForest::get_grid_leaves)
+      .method("set_grid_leaves", &RandomPlantedForest::set_grid_leaves)
       .method("is_purified", &RandomPlantedForest::is_purified);
 
   class_<ClassificationRPF>("ClassificationRPF")
       .derives<RandomPlantedForest>("RandomPlantedForest")
       .constructor<const NumericMatrix, const NumericMatrix, const String, const NumericVector>()
+      .constructor<const String, const NumericVector>()
       .method("set_parameters", &ClassificationRPF::set_parameters);
 }
diff --git a/tests/testthat/test-bundle.R b/tests/testthat/test-bundle.R
new file mode 100644
index 0000000..236aaa6
--- /dev/null
+++ b/tests/testthat/test-bundle.R
@@ -0,0 +1,9 @@
+test_that("bundle round-trip", {
+  skip_if_not_installed("bundle")
+  fit <- rpf(mpg ~ wt + cyl, data = mtcars, ntrees = 5)
+  b <- bundle::bundle(fit)
+  tmp <- tempfile(fileext = ".rds")
+  saveRDS(b, tmp)
+  restored <- bundle::unbundle(readRDS(tmp))
+  expect_identical(predict(restored, mtcars), predict(fit, mtcars))
+})
diff --git a/tests/testthat/test-marshal.R b/tests/testthat/test-marshal.R
new file mode 100644
index 0000000..55c5cf9
--- /dev/null
+++ b/tests/testthat/test-marshal.R
@@ -0,0 +1,181 @@
+test_that("params-only constructors create empty, unfitted objects", {
+  # 13-element parameter vector: max_interaction, ntrees, splits, split_try,
+  # t_try, purify, deterministic, nthreads, cv, split_decay_rate,
+  # max_candidates, delete_leaves, split_mode
+  pars <- c(1, 50, 30, 10, 0.4, 0, 0, 1, 0, 0.1, 50, 1, 3)
+  fit <- new(RandomPlantedForest, pars)
+  expect_length(fit$get_model(), 0)
+  expect_false(fit$is_purified())
+
+  cfit <- new(ClassificationRPF, "L2", c(pars, 0, 0.1))
+  expect_length(cfit$get_model(), 0)
+
+  expect_error(new(ClassificationRPF, "L2", pars), "15 parameters")
+})
+
+test_that("get_data returns training data from a fitted forest", {
+  X <- as.matrix(mtcars[, c("wt", "cyl")])
+  Y <- as.matrix(mtcars$mpg)
+  fit <- new(RandomPlantedForest, Y, X, c(1, 10, 30, 10, 0.4, 0, 0, 1, 0, 0.1, 50, 1, 3))
+  dat <- fit$get_data()
+  expect_equal(unname(as.matrix(dat$X)), unname(X))
+  expect_equal(unname(as.matrix(dat$Y)), unname(Y))
+})
+
+expect_roundtrip <- function(fit, new_data) {
+  p_before <- predict(fit, new_data)
+  blob <- rpf_marshal(fit)
+  tmp <- tempfile(fileext = ".rds")
+  saveRDS(blob, tmp)
+  restored <- rpf_unmarshal(readRDS(tmp))
+  expect_s3_class(restored, "rpf")
+  expect_identical(predict(restored, new_data), p_before)
+  restored
+}
+
+test_that("regression round-trip, no data", {
+  fit <- rpf(mpg ~ wt + cyl, data = mtcars, ntrees = 10)
+  expect_roundtrip(fit, mtcars)
+})
+
+test_that("binary classification round-trip, no data", {
+  # NOTE: `am` is pre-converted to a factor column rather than using
+  # `as.factor(am)` inline in the formula - hardhat::mold() does not capture
+  # an inline outcome transform in the blueprint ptype, which breaks
+  # `predict()` independently of marshaling (pre-existing hardhat/formula
+  # limitation, not specific to this task).
+  dat <- transform(mtcars, am = factor(am))
+  fit <- rpf(am ~ wt + hp, data = dat, ntrees = 10, loss = "logit")
+  restored <- expect_roundtrip(fit, dat)
+  expect_identical(
+    predict(restored, dat, type = "prob"),
+    predict(fit, dat, type = "prob")
+  )
+})
+
+test_that("multiclass round-trip, no data", {
+  fit <- rpf(Species ~ ., data = iris, ntrees = 10, loss = "L2")
+  restored <- expect_roundtrip(fit, iris)
+  expect_identical(
+    predict(restored, iris, type = "prob"),
+    predict(fit, iris, type = "prob")
+  )
+})
+
+test_that("purified regression round-trip without data", {
+  fit <- rpf(mpg ~ wt + cyl + hp, data = mtcars, ntrees = 10, max_interaction = 2)
+  purify(fit)
+  comp_before <- predict_components(fit, mtcars)
+  restored <- expect_roundtrip(fit, mtcars)
+  expect_true(is_purified(restored))
+  expect_identical(predict_components(restored, mtcars), comp_before)
+})
+
+test_that("purified multiclass round-trip without data", {
+  fit <- rpf(Species ~ ., data = iris, ntrees = 10, loss = "L2")
+  purify(fit)
+  comp_before <- predict_components(fit, iris)
+  restored <- expect_roundtrip(fit, iris)
+  expect_true(is_purified(restored))
+  expect_identical(predict_components(restored, iris), comp_before)
+})
+
+test_that("marshaled blob contains no external pointers", {
+  fit <- rpf(mpg ~ wt + cyl, data = mtcars, ntrees = 5)
+  blob <- rpf_marshal(fit)
+  has_extptr <- function(x) {
+    if (typeof(x) == "externalptr") {
+      return(TRUE)
+    }
+    # unclass so classed list-likes (e.g. package_version, whose `[[` returns
+    # another package_version) recurse over their raw contents and terminate
+    x <- unclass(x)
+    if (is.list(x)) {
+      return(any(vapply(x, has_extptr, logical(1))))
+    }
+    FALSE
+  }
+  expect_false(has_extptr(unclass(blob)))
+})
+
+test_that("include_data = TRUE allows purify after restore", {
+  fit <- rpf(mpg ~ wt + cyl, data = mtcars, ntrees = 10, max_interaction = 2)
+
+  restored <- rpf_unmarshal(rpf_marshal(fit, include_data = TRUE))
+  purify(restored)
+
+  # restored shares the identical forest and training data with fit,
+  # so purifying each independently must give identical predictions
+  purify(fit)
+  expect_identical(predict(restored, mtcars), predict(fit, mtcars))
+})
+
+test_that("purify on a data-free restored forest errors informatively", {
+  fit <- rpf(mpg ~ wt + cyl, data = mtcars, ntrees = 5)
+  restored <- rpf_unmarshal(rpf_marshal(fit))
+  expect_error(purify(restored), "include_data")
+})
+
+test_that("restored-without-marshal rpf gives actionable error", {
+  fit <- rpf(mpg ~ wt + cyl, data = mtcars, ntrees = 5)
+  tmp <- tempfile(fileext = ".rds")
+  saveRDS(fit, tmp)
+  zombie <- readRDS(tmp)
+  expect_false(rpf_is_valid(zombie))
+  expect_error(predict(zombie, mtcars), "rpf_marshal")
+  expect_error(purify(zombie), "rpf_marshal")
+  expect_error(predict_components(zombie, mtcars), "rpf_marshal")
+  expect_error(rpf_marshal(zombie), "readRDS")
+  expect_true(rpf_is_valid(fit))
+})
+
+test_that("rpf_unmarshal errors on malformed blob missing fit_state", {
+  expect_error(rpf_unmarshal(structure(list(), class = "rpf_marshaled")), "Unsupported")
+})
+
+test_that("corrupt blobs error instead of reading out of bounds", {
+  fit <- rpf(mpg ~ wt + cyl, data = mtcars, ntrees = 5)
+
+  # Interval matrix narrower than feature_size
+  blob <- rpf_marshal(fit)
+  blob$fit_state$model[[1]]$intervals[[1]][[1]] <- matrix(0, nrow = 2, ncol = 1)
+  expect_error(rpf_unmarshal(blob), "interval matrix")
+
+  # Purified grid: values matrix with wrong dimensions
+  pfit <- rpf(mpg ~ wt + cyl, data = mtcars, ntrees = 5)
+  purify(pfit)
+  # trees[[1]] is the null tree whose 1x1 values matrix is legitimately tiny;
+  # corrupt a real component tree instead
+  pblob <- rpf_marshal(pfit)
+  pblob$fit_state$grid[[1]]$trees[[2]]$values <- matrix(0, nrow = 1, ncol = 1)
+  expect_error(rpf_unmarshal(pblob), "values matrix")
+
+  # Purified grid: non-positive grid dimension
+  pblob2 <- rpf_marshal(pfit)
+  pblob2$fit_state$grid[[1]]$trees[[2]]$dims <- -1L
+  expect_error(rpf_unmarshal(pblob2), "non-positive")
+
+  # Training data with wrong number of feature columns
+  dblob <- rpf_marshal(fit, include_data = TRUE)
+  dblob$fit_state$data$X <- dblob$fit_state$data$X[, 1, drop = FALSE]
+  expect_error(rpf_unmarshal(dblob), "columns")
+})
+
+test_that("exported forest is rebuilt instead of duplicated in the blob", {
+  fit <- rpf(mpg ~ wt + cyl, data = mtcars, ntrees = 5, export_forest = TRUE)
+  blob <- rpf_marshal(fit)
+  expect_null(blob$forest)
+  restored <- rpf_unmarshal(blob)
+  expect_s3_class(restored$forest, "rpf_forest")
+  expect_identical(restored$forest, fit$forest)
+})
+
+test_that("blob records package version and warns when saved with a newer one", {
+  fit <- rpf(mpg ~ wt + cyl, data = mtcars, ntrees = 5)
+  blob <- rpf_marshal(fit)
+  expect_identical(blob$fit_state$pkg_version, utils::packageVersion("randomPlantedForest"))
+
+  blob$fit_state$pkg_version <- package_version("999.0.0")
+  expect_warning(restored <- rpf_unmarshal(blob), "999.0.0")
+  expect_identical(predict(restored, mtcars), predict(fit, mtcars))
+})