diff --git a/stan/math/fwd/prob/std_normal_log_qf.hpp b/stan/math/fwd/prob/std_normal_log_qf.hpp
index 2dc0a855cf7..52f53dc107f 100644
--- a/stan/math/fwd/prob/std_normal_log_qf.hpp
+++ b/stan/math/fwd/prob/std_normal_log_qf.hpp
@@ -16,15 +16,7 @@ namespace math {
 template <typename T>
 inline fvar<T> std_normal_log_qf(const fvar<T>& p) {
   const T xv = std_normal_log_qf(p.val_);
-  int p_sign = 1;
-  auto p_d = p.d_;
-  if (p.d_ < 0) {
-    p_sign = -1;
-    p_d *= -1;
-  }
-  return fvar<T>(
-      xv,
-      p_sign * exp(p.val_ + log(p_d) - NEG_LOG_SQRT_TWO_PI + 0.5 * square(xv)));
+  return fvar<T>(xv, p.d_ * exp(p.val_ - std_normal_lpdf(xv)));
 }
 }  // namespace math
 }  // namespace stan
diff --git a/test/unit/math/mix/prob/std_normal_log_qf_test.cpp b/test/unit/math/mix/prob/std_normal_log_qf_test.cpp
index 1ae44dd9c27..8a2908da869 100644
--- a/test/unit/math/mix/prob/std_normal_log_qf_test.cpp
+++ b/test/unit/math/mix/prob/std_normal_log_qf_test.cpp
@@ -46,3 +46,14 @@ TEST_F(AgradRev, mathMixMatFunLog_stdNormalLogQfVarmat) {
   }
   expect_ad_vector_matvar(f, A);
 }
+
+TEST_F(AgradRev, GradientStabilityStdNormalLogQf) {
+  auto f = [](const auto& y) {
+    return stan::math::sum(stan::math::std_normal_log_qf(y));
+  };
+
+  Eigen::VectorXd y1(2);
+  y1 << -10, -2;
+
+  stan::test::expect_ad(f, y1);
+}