diff --git a/DESCRIPTION b/DESCRIPTION
index 3de10d6..a483e4b 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,7 +1,7 @@
 Package: BVAR
 Type: Package
 Title: Hierarchical Bayesian Vector Autoregression
-Version: 0.1.6.9000
+Version: 0.1.6
 Date: 2019-08-28
 Authors@R: c(person("Nikolas",  "Kuschnig", role = c("aut", "cre"), email = "nikolas.kuschnig@wu.ac.at"),
     person("Lukas", "Vashold", role = "aut"),
@@ -15,7 +15,7 @@ URL: https://github.com/nk027/bvar
 BugReports: https://github.com/nk027/bvar/issues
 Depends: R (>= 3.3.0)
 Imports: MASS, stats, graphics, utils
-Suggests: vars, coda
+Suggests: vars
 License: GPL-3
 Encoding: UTF-8
 LazyData: true
diff --git a/NAMESPACE b/NAMESPACE
index 5c9af7b..a4e802e 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -41,9 +41,8 @@ export(bv_metropolis)
 export(bv_mh)
 export(bv_minnesota)
 export(bv_mn)
+export(bv_plot)
 export(bv_plot_density)
-export(bv_plot_fcast)
-export(bv_plot_irf)
 export(bv_plot_trace)
 export(bv_priors)
 export(bv_psi)
@@ -60,13 +59,17 @@ importFrom(graphics,par)
 importFrom(graphics,plot)
 importFrom(graphics,polygon)
 importFrom(stats,arima)
+importFrom(stats,coef)
 importFrom(stats,density)
 importFrom(stats,dgamma)
+importFrom(stats,fitted)
 importFrom(stats,optim)
 importFrom(stats,predict)
 importFrom(stats,quantile)
+importFrom(stats,residuals)
 importFrom(stats,rnorm)
 importFrom(stats,runif)
 importFrom(stats,ts.plot)
+importFrom(stats,vcov)
 importFrom(utils,setTxtProgressBar)
 importFrom(utils,txtProgressBar)
diff --git a/R/80_plot.R b/R/80_plot.R
index 754ebfe..dc24ee4 100644
--- a/R/80_plot.R
+++ b/R/80_plot.R
@@ -34,6 +34,9 @@
 #'
 #' # Compare with second chain
 #' plot(x, chains = y)
+#'
+#' # Only plot the marginal likelihood's density
+#' plot(x, "dens", "ml")
 #' }
 plot.bvar <- function(
   x,
diff --git a/R/91_fitted.R b/R/91_fitted.R
index a91406c..2551b31 100644
--- a/R/91_fitted.R
+++ b/R/91_fitted.R
@@ -7,9 +7,6 @@
 #' @param conf_bands Numeric vector of desired confidence bands to apply.
 #' E.g. for bands at 5\%, 10\%, 90\% and 95\% set this to \code{c(0.05, 0.1)}.
 #' Note that the median, i.e. 0.5 is always included.
-#' @param n_thin Integer scalar. Every \emph{n_thin}'th draw in \emph{object}
-#' is used for forecasting, others are dropped. Defaults to \code{100L}, to
-#' prevent memory overflow.
 #'
 #' @param x Object of class \code{bvar_fitted} / \code{bvar_resid}.
 #' @param digits Integer scalar. Fed to \code{\link[base]{round}} and applied to
@@ -26,6 +23,8 @@
 #'
 #' @export
 #'
+#' @importFrom stats fitted residuals
+#'
 #' @examples
 #' \donttest{
 #' data <- matrix(rnorm(200), ncol = 2)
diff --git a/R/92_coef.R b/R/92_coef.R
index df85517..ce4c78c 100644
--- a/R/92_coef.R
+++ b/R/92_coef.R
@@ -20,6 +20,8 @@
 #'
 #' @export
 #'
+#' @importFrom stats coef vcov
+#'
 #' @examples
 #' \donttest{
 #' data <- matrix(rnorm(200), ncol = 2)
diff --git a/R/deprecated.R b/R/deprecated.R
index 9cff6e8..1f140d4 100644
--- a/R/deprecated.R
+++ b/R/deprecated.R
@@ -105,16 +105,14 @@ bv_plot_irf <- function(
   variables = NULL,
   vars_impulse = NULL,
   vars_response = NULL,
-  orientation = c("vertical", "horizontal"),
   mar = c(2, 2, 2, 0.5),
   ...) {
 
   .Deprecated("plot.bvar_irf")
-  plot_fcast(
+  plot_irf(
     x, conf_bands,
     variables = variables,
     vars_impulse = vars_impulse, vars_response = vars_response,
-    orientation = orientation,
     mar = mar, ... = ...)
 }
 
diff --git a/R/methods.R b/R/methods.R
deleted file mode 100644
index 512786a..0000000
--- a/R/methods.R
+++ /dev/null
@@ -1,10 +0,0 @@
-# summary.bvar <- function(x) {
-#
-# }
-#
-# logLik.bvar <- function(x) {}
-#
-# Phi.bvar <- function(x) {}
-#
-# Psi.bvar <- function(x) {}
-#
diff --git a/man/bv_plot_trace.Rd b/man/bv_plot_trace.Rd
index d554ce2..dd07e9a 100644
--- a/man/bv_plot_trace.Rd
+++ b/man/bv_plot_trace.Rd
@@ -20,6 +20,7 @@ marginal likelihood.}
 hyperparameter must be available and priors of the objects should match.}
 }
 \description{
+\strong{Deprecated. Use \code{\link{plot.bvar}} instead.}
 Diagnostic plots of the trace / density of a single hyperparameter.
 A parameter may be plotted across multiple iterations of \code{\link{bvar}}
 via the ellipsis parameter. Given that the settings for \code{\link{bvar}}
diff --git a/man/coef.bvar.Rd b/man/coef.bvar.Rd
index d98bc2a..8002052 100644
--- a/man/coef.bvar.Rd
+++ b/man/coef.bvar.Rd
@@ -26,8 +26,8 @@ Note that the median, i.e. 0.5 is always included.}
 
 \item{x}{Object of class \code{bvar_coefs} / \code{bvar_vcovs}.}
 
-\item{digits}{Integer scalar. Fed to \code{\link[base]{round}} and applied to
-numeric outputs (i.e. the quantiles).}
+\item{digits}{Integer scalar. Fed to \code{\link[base]{round}} and applied
+to numeric outputs (i.e. the quantiles).}
 }
 \value{
 Returns a numeric array of class \code{bvar_coefs} /
diff --git a/man/fitted.bvar.Rd b/man/fitted.bvar.Rd
index c16dc0b..c62b8b0 100644
--- a/man/fitted.bvar.Rd
+++ b/man/fitted.bvar.Rd
@@ -8,9 +8,9 @@
 \alias{print.bvar_resid}
 \title{Fitted and residual methods for Bayesian VARs}
 \usage{
-\method{fitted}{bvar}(object, conf_bands = 0.5, n_thin = 100L, ...)
+\method{fitted}{bvar}(object, conf_bands = 0.5, ...)
 
-\method{residuals}{bvar}(object, conf_bands = 0.5, n_thin = 100L, ...)
+\method{residuals}{bvar}(object, conf_bands = 0.5, ...)
 
 \method{plot}{bvar_resid}(x, vars = NULL, mar = c(2, 2, 2, 0.5), ...)
 
@@ -25,10 +25,6 @@
 E.g. for bands at 5\%, 10\%, 90\% and 95\% set this to \code{c(0.05, 0.1)}.
 Note that the median, i.e. 0.5 is always included.}
 
-\item{n_thin}{Integer scalar. Every \emph{n_thin}'th draw in \emph{object}
-is used for forecasting, others are dropped. Defaults to \code{100L}, to
-prevent memory overflow.}
-
 \item{...}{Other graphical parameters for \code{\link[graphics]{par}}.}
 
 \item{x}{Object of class \code{bvar_fitted} / \code{bvar_resid}.}
@@ -57,9 +53,6 @@ x <- bvar(data, lags = 2)
 # Get fitted values and adjust confidence bands to 10\%, 50\% and 90\%
 fitted(x, conf_bands = 0.10)
 
-# Only get the median and increase the iterations
-fitted(x, conf_bands = 0.5, n_thin = 10L)
-
 # Get residuals
 residuals(x)
 }
diff --git a/man/plot.bvar.Rd b/man/plot.bvar.Rd
index 8a7370d..21d842d 100644
--- a/man/plot.bvar.Rd
+++ b/man/plot.bvar.Rd
@@ -1,15 +1,18 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/80_plot.R
+% Please edit documentation in R/80_plot.R, R/deprecated.R
 \name{plot.bvar}
 \alias{plot.bvar}
+\alias{plot_bvar}
 \alias{bv_plot}
 \title{Plotting method for Bayesian VARs}
 \usage{
 \method{plot}{bvar}(x, type = c("full", "trace", "density"),
   vars = NULL, chains = list(), mar = c(2, 2, 2, 0.5), ...)
 
-bv_plot(x, type = c("full", "trace", "density"), vars = NULL,
+plot_bvar(x, type = c("full", "trace", "density"), vars = NULL,
   chains = list(), mar = c(2, 2, 2, 0.5), ...)
+
+bv_plot(x, mar = c(2, 2, 2, 0.5), ...)
 }
 \arguments{
 \item{x}{A \code{bvar} object, obtained from \code{\link{bvar}}.}
@@ -46,6 +49,9 @@ plot(x)
 
 # Compare with second chain
 plot(x, chains = y)
+
+# Only plot the marginal likelihood's density
+plot(x, "dens", "ml")
 }
 }
 \seealso{
diff --git a/man/plot.bvar_fcast.Rd b/man/plot.bvar_fcast.Rd
index 6feb8e8..b6745d9 100644
--- a/man/plot.bvar_fcast.Rd
+++ b/man/plot.bvar_fcast.Rd
@@ -1,7 +1,8 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/58_fcast_plot.R
+% Please edit documentation in R/58_fcast_plot.R, R/deprecated.R
 \name{plot.bvar_fcast}
 \alias{plot.bvar_fcast}
+\alias{plot_fcast}
 \alias{bv_plot_fcast}
 \title{Plotting method for Bayesian VAR forecasts}
 \usage{
@@ -9,7 +10,11 @@
   orientation = c("vertical", "horizontal"), mar = c(2, 2, 2, 0.5),
   ...)
 
-bv_plot_fcast(x, conf_bands, vars = NULL, variables = NULL,
+plot_fcast(x, conf_bands, vars = NULL, variables = NULL,
+  orientation = c("vertical", "horizontal"), mar = c(2, 2, 2, 0.5),
+  ...)
+
+bv_plot_fcast(x, conf_bands = 0.16, variables = NULL, vars = NULL,
   orientation = c("vertical", "horizontal"), mar = c(2, 2, 2, 0.5),
   ...)
 }
diff --git a/man/plot.bvar_irf.Rd b/man/plot.bvar_irf.Rd
index b5ac932..203cc24 100644
--- a/man/plot.bvar_irf.Rd
+++ b/man/plot.bvar_irf.Rd
@@ -1,15 +1,20 @@
 % Generated by roxygen2: do not edit by hand
-% Please edit documentation in R/68_irf_plot.R
+% Please edit documentation in R/68_irf_plot.R, R/deprecated.R
 \name{plot.bvar_irf}
 \alias{plot.bvar_irf}
+\alias{plot_irf}
 \alias{bv_plot_irf}
 \title{Plotting method for Bayesian VAR impulse responses}
 \usage{
 \method{plot}{bvar_irf}(x, conf_bands, vars_impulse = NULL,
   vars_response = NULL, variables = NULL, mar = c(2, 2, 2, 0.5), ...)
 
-bv_plot_irf(x, conf_bands, vars_impulse = NULL, vars_response = NULL,
+plot_irf(x, conf_bands, vars_impulse = NULL, vars_response = NULL,
   variables = NULL, mar = c(2, 2, 2, 0.5), ...)
+
+bv_plot_irf(x, conf_bands = 0.16, variables = NULL,
+  vars_impulse = NULL, vars_response = NULL, mar = c(2, 2, 2, 0.5),
+  ...)
 }
 \arguments{
 \item{x}{A \code{bvar} / \code{bvar_irf} object, obtained from