diff --git a/src/stan/mcmc/chainset.hpp b/src/stan/mcmc/chainset.hpp
index 9eb9369621..e380504ecc 100644
--- a/src/stan/mcmc/chainset.hpp
+++ b/src/stan/mcmc/chainset.hpp
@@ -299,10 +299,15 @@ class chainset {
    * @return parameter value at quantile
    */
   double quantile(const int index, const double prob) const {
-    // Ensure the probability is within [0, 1]
     Eigen::MatrixXd draws = samples(index);
     Eigen::Map<Eigen::VectorXd> map(draws.data(), draws.size());
-    return stan::math::quantile(map, prob);
+    double result;
+    try {
+      result = stan::math::quantile(map, prob);
+    } catch (const std::logic_error& e) {
+      return std::numeric_limits<double>::quiet_NaN();
+    }
+    return result;
   }
 
   /**
@@ -332,7 +337,14 @@ class chainset {
     Eigen::MatrixXd draws = samples(index);
     Eigen::Map<Eigen::VectorXd> map(draws.data(), draws.size());
     std::vector<double> probs_vec(probs.data(), probs.data() + probs.size());
-    std::vector<double> quantiles = stan::math::quantile(map, probs_vec);
+    std::vector<double> quantiles;
+    try {
+      quantiles = stan::math::quantile(map, probs_vec);
+    } catch (const std::logic_error& e) {
+      Eigen::VectorXd nans(probs.size());
+      nans.setConstant(std::numeric_limits<double>::quiet_NaN());
+      return nans;
+    }
     return Eigen::Map<Eigen::VectorXd>(quantiles.data(), quantiles.size());
   }
 
diff --git a/src/test/unit/mcmc/chainset_test.cpp b/src/test/unit/mcmc/chainset_test.cpp
index 01cb18da42..0c26069876 100644
--- a/src/test/unit/mcmc/chainset_test.cpp
+++ b/src/test/unit/mcmc/chainset_test.cpp
@@ -220,6 +220,6 @@ TEST_F(McmcChains, quantile_tests) {
   probs << 0.0, 0.01, 0.05, 0.95, 0.99, 1.0;
   Eigen::VectorXd stepsize_quantiles = datagen_chains.quantiles("stepsize__", probs);
   for (size_t i = 0; i < probs.size(); ++i) {
-    EXPECT_EQ(stepsize_quantiles(i), std::numeric_limits<double>::quiet_NaN());
+    EXPECT_TRUE(std::isnan(stepsize_quantiles(i)));
   }
 }