diff --git a/DESCRIPTION b/DESCRIPTION
index 7d0563f..8abaa30 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -45,6 +45,7 @@ Suggests:
     knitr,
     nnet,
     rmarkdown,
+    scales,
     spelling,
     testthat,
     tidyr
diff --git a/inst/WORDLIST b/inst/WORDLIST
index 279b54b..e75b5a4 100644
--- a/inst/WORDLIST
+++ b/inst/WORDLIST
@@ -14,6 +14,7 @@ Fienberg
 frac
 fulltime
 ggplot
+ggeffects
 gradeschool
 GSS
 https
diff --git a/vignettes/fig/wlf-alt-plot-1.png b/vignettes/fig/wlf-alt-plot-1.png
index 99e5e76..c35f974 100644
Binary files a/vignettes/fig/wlf-alt-plot-1.png and b/vignettes/fig/wlf-alt-plot-1.png differ
diff --git a/vignettes/fig/wlf-effect-plot-1-1.png b/vignettes/fig/wlf-effect-plot-1-1.png
index 2ce6dd3..961ab5c 100644
Binary files a/vignettes/fig/wlf-effect-plot-1-1.png and b/vignettes/fig/wlf-effect-plot-1-1.png differ
diff --git a/vignettes/fig/wlf-effect-plot-2-1.png b/vignettes/fig/wlf-effect-plot-2-1.png
index 0a4b008..d58bbd1 100644
Binary files a/vignettes/fig/wlf-effect-plot-2-1.png and b/vignettes/fig/wlf-effect-plot-2-1.png differ
diff --git a/vignettes/fig/wlf-plot-1.png b/vignettes/fig/wlf-plot-1.png
index e8acfa1..04a23b9 100644
Binary files a/vignettes/fig/wlf-plot-1.png and b/vignettes/fig/wlf-plot-1.png differ
diff --git a/vignettes/nestedLogit.Rmd b/vignettes/nestedLogit.Rmd
index 7d1d775..ca7668b 100644
--- a/vignettes/nestedLogit.Rmd
+++ b/vignettes/nestedLogit.Rmd
@@ -459,10 +459,15 @@ various levels of the other predictors, and to compose these into a single figur
 #| fig.cap = "**plot method**: Predicted probabilities and 95 percent pointwise confidence envelopes of not working, working part-time, and working full-time",
 #| fig.show = "hold"
 op <- par(mfcol=c(1, 2), mar=c(4, 4, 3, 1) + 0.1)
-plot(wlf.nested, "hincome", list(children="absent"), # left panel
-     xlab="Husband's Income", legend.location="top")
-plot(wlf.nested, "hincome", list(children="present"), # right panel
-     xlab="Husband's Income", legend=FALSE)
+col <- scales::hue_pal()(3)                  # ggplot discrete colors
+plot(wlf.nested, "hincome",                  # left panel
+     other = list(children="absent"),
+     xlab = "Husband's Income", 
+     legend.location="top", col = col)
+plot(wlf.nested, "hincome",                  # right panel
+     other = list(children="present"),
+     xlab = "Husband's Income", 
+     legend=FALSE, col = col)
 par(op)
 ```
 
@@ -491,7 +496,8 @@ An alternative style of effect plot shows the fit as a stacked-area graph (in th
 #| fig.height = 4,
 #| fig.cap = "Stacked-area predictor effect plots for the nested-logit model fit to the `Womenlf` data"
 plot(predictorEffects(wlf.nested),
-     axes=list(y=list(style="stacked")))
+     axes=list(y=list(style="stacked")),
+     lines=list(col=scales::hue_pal()(3)))
 ```
 
 The computation and display of effect graphs are highly customizable. For details, see the documentation for the **effects** package.
@@ -553,10 +559,15 @@ Here's a graph of the alternative model:
 #| fig.cap = "Predicted probabilities for the alternative nested-dichotomies model",
 #| fig.show = "hold"
 op <- par(mfcol=c(1, 2), mar=c(4, 4, 3, 1) + 0.1)
-plot(wlf.nested.alt, "hincome", list(children="absent"), # left panel
-     xlab="Husband's Income", legend.location="top")
-plot(wlf.nested.alt, "hincome", list(children="present"), # right panel
-     xlab="Husband's Income", legend=FALSE)
+col <- scales::hue_pal()(3)
+plot(wlf.nested.alt, "hincome",           # left panel
+     others = list(children="absent"), 
+     xlab="Husband's Income", 
+     legend.location="top", col = col)
+plot(wlf.nested.alt, "hincome",           # right panel
+     others = list(children="present"),
+     xlab="Husband's Income", 
+     legend=FALSE, col = col)
 par(op)
 ```
 Compare this to the previous graph for the original specification.