diff --git a/NEWS.md b/NEWS.md
index 02d50a85f..d70bf061a 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -88,15 +88,20 @@ mu-referencing style to run the optimization.
   intervals without bands on each of the percentiles; You can also
   choose not to match the secondary bands limits with `levels` but use
   your own `ci=0.99` for instance
-  
+
 - A new function was introduced `meanProbs()` which calculates the
   mean and expected quantiles under either the normal or t
   distribution
-    
+
 - When calculating the intervals for `rxode2` simulated objects you
   can also use `mean=TRUE` to use the mean for the first level of
   confidence using `meanProbs()`
 
+- When plotting the `confint` derived intervals from an `rxode2`
+  simulation, you can now subset based on a simulated value like
+  `plot(ci, Cc)` which will only plot the variable `Cc` that you
+  summarized even if you also summarized `eff` (for instance).
+
 ## Internal new features
 
 - Add `as.model()` for list expressions, which implies `model(ui) <-
diff --git a/R/plot.R b/R/plot.R
index f2d486a8c..54512602b 100644
--- a/R/plot.R
+++ b/R/plot.R
@@ -358,8 +358,26 @@ plot.rxSolve <- function(x, y, ..., log = "", xlab = "Time", ylab = "") {
 #' @export
 plot.rxSolveConfint1 <- function(x, y, ..., xlab = "Time", ylab = "", log = "") {
   .data <- NULL
+  .y <- as.character(substitute(y))
+  .call0 <- match.call()[-(1:2)]
+  .call <- as.list(.call0)
+  .w <- which(names(.call) %in% c("x", "y", "log", "xlab", "ylab"))
+  if (length(.w) > 0) {
+    .call <- .call[-.w]
+  }
+  .cmts <- c(
+    as.character(substitute(y)),
+    names(vapply(as.character(.call), `c`, character(1), USE.NAMES=FALSE)),
+    as.character(unlist(.call))
+  )
+  .cmts <- .cmts[.cmts != ""]
+  .cmts <- unique(.cmts)
   .lvl <- attr(class(x), ".rx")$lvl
   .parm <- attr(class(x), ".rx")$parm
+  if (length(.cmts) > 0) {
+    .parm <- intersect(.parm, .cmts)
+    x  <- x[x$trt %in% .parm,]
+  }
   .by <- attr(class(x), ".rx")$by
   .aes <- aes(.data$time, .data$eff)
   .facet <- NULL
@@ -421,9 +439,27 @@ plot.rxSolveConfint1 <- function(x, y, ..., xlab = "Time", ylab = "", log = "")
 #' @export
 plot.rxSolveConfint2 <- function(x, y, ..., xlab = "Time", ylab = "", log = "") {
   .data <- NULL
+  .y <- as.character(substitute(y))
+  .call0 <- match.call()[-(1:2)]
+  .call <- as.list(.call0)
+  .w <- which(names(.call) %in% c("x", "y", "log", "xlab", "ylab"))
+  if (length(.w) > 0) {
+    .call <- .call[-.w]
+  }
+  .cmts <- c(
+    as.character(substitute(y)),
+    names(vapply(as.character(.call), `c`, character(1), USE.NAMES=FALSE)),
+    as.character(unlist(.call))
+  )
+  .cmts <- .cmts[.cmts != ""]
+  .cmts <- unique(.cmts)
+  .parm <- attr(class(x), ".rx")$parm
+  if (length(.cmts) > 0) {
+    .parm <- intersect(.parm, .cmts)
+    x  <- x[x$trt %in% .parm,]
+  }
   .lvl <- attr(class(x), ".rx")$lvl
   .ci  <- attr(class(x), ".rx")$ci
-  .parm <- attr(class(x), ".rx")$parm
   .by <- attr(class(x), ".rx")$by
   .aes <- aes(.data$time, .data$p50,
     color = .data$Percentile,
diff --git a/tests/testthat/test-zzzz-plot.R b/tests/testthat/test-zzzz-plot.R
index 4085b0648..9823a3ff2 100644
--- a/tests/testthat/test-zzzz-plot.R
+++ b/tests/testthat/test-zzzz-plot.R
@@ -153,6 +153,8 @@ rxTest({
 
       ci1.C2 <- confint(sim, "C2", ci=0.99)
 
+      ci1.C3 <- confint(sim,  ci=0.99)
+
       ci1.C2.e <- confint(sim, "C2", by="extra")
 
       ci1.C2.eff <- confint(sim, c("C2", "eff"))
@@ -160,11 +162,11 @@ rxTest({
       ci1.C2.eff.e <- confint(sim, c("C2", "eff"), by="extra")
 
       sim2 <- rxSolve(m2, ev, omega = omega, nSub = 2500, keep="extra")
-      
+
       sim2R <- rxSolve(m2, evR, omega = omega, nSub = 2500)
 
       ci2.C2 <- confint(sim2, "C2")
-      
+
       ci2.C2.e <- confint(sim2, "C2", by="extra")
 
       ci2.C2.eff <- confint(sim2, c("C2", "eff"))
@@ -203,6 +205,8 @@ rxTest({
           vdiffr::expect_doppelganger(paste0("plot-", .xgxtxt, "all-log-yx-r"), suppressWarnings(sR %>% plot(log = "yx")))
 
           vdiffr::expect_doppelganger(paste0("plot-ci1c2", .xgxtxt), suppressWarnings(ci1.C2 %>% plot()))
+          vdiffr::expect_doppelganger(paste0("plot-ci1c3", .xgxtxt, "-centr"), suppressWarnings(ci1.C3 %>% plot("centr")))
+          vdiffr::expect_doppelganger(paste0("plot-ci1c3", .xgxtxt, "-full"), suppressWarnings(ci1.C3 %>% plot()))
           vdiffr::expect_doppelganger(paste0("plot-ci1c2", .xgxtxt, "log-x"), suppressWarnings(ci1.C2 %>% plot(log = "x")))
           vdiffr::expect_doppelganger(paste0("plot-ci1c2", .xgxtxt, "log-y"), suppressWarnings(ci1.C2 %>% plot(log = "y")))
           vdiffr::expect_doppelganger(paste0("plot-ci1c2", .xgxtxt, "log-xy"), suppressWarnings(ci1.C2 %>% plot(log = "xy")))