diff --git a/DESCRIPTION b/DESCRIPTION
index b102dfd..0c77cb5 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,6 +1,6 @@
 Package: PaleoSpec
 Title: Spectral tools for the ECUS group
-Version: 0.32
+Version: 0.33
 Authors@R: c(
   person("Thomas", "Laepple", email = "tlaepple@awi.de", role = c("aut", "cre")),
   person("Thomas", "Muench", email = "tmuench@awi.de", role = c("aut")),
diff --git a/NEWS.md b/NEWS.md
index f4c4208..8d43f92 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -1,3 +1,14 @@
+# PaleoSpec 0.33
+
+* AddConfInterval calculation changed calculated confidence intervals give the 
+  interval within which the true spectrum should lie with frequency 1-p. 
+  Previous versions returned the interval which, when applied to the 'true' 
+  spectrum, would contain new spectral estimates with frequency 1-p.
+  
+* Improvements to calculation of DOF in FilterSpec and FilterSpecLog. DOF 
+  calculation now accounts for non independence of ordinates from existing 
+  filter or tapering.
+
 # PaleoSpec 0.32
 
 * SpecACF can now use slepian tapers
diff --git a/R/AddConfInterval.R b/R/AddConfInterval.R
index 9f17c7d..70b29e0 100644
--- a/R/AddConfInterval.R
+++ b/R/AddConfInterval.R
@@ -14,33 +14,62 @@
 #' @return the input object including the new list elements \code{lim.1} and
 #'   \code{lim.2} giving the upper and lower bound of the confidence interval,
 #'   respectively.
-#' @author Thomas Laepple
+#' @description
+#' Calculated confidence intervals give the interval within which the true
+#' spectrum should lie with frequency 1-p. This behaviour changed with version
+#' 0.33.
+#'
+#' Previous versions returned the interval which, when applied to the 'true' spectrum, would contain new spectral estimates
+#' with frequency 1-p.
+#'
+#' @author Thomas Laepple, Andrew Dolman
 #' @examples
 #'
-#' N.R <- 1000
-#' N.T <- 100
-#' save.spec <- matrix(NA, N.T / 2, N.R)
-#' for (i.R in 1 : N.R) {
-#'   save.spec[, i.R] <- SpecMTM(ts(SimPowerlaw(1, N.T)))$spec
-#' }
+#' alpha <- 0.1
+#' beta <- 1
+#'
+#' spec_sim <- SpecMTM(ts(SimPLS(N = 1e03, beta = beta, alpha = alpha)))
+#'
+#' # Using a large nominal p value of 0.25 to reduce variation between random
+#' # timeseries
+#'
+#' # The true spec should be p/2 times above the lower CI, and p/2 below the upper CI
+#'
+#' nominal_p <- 0.25
+#' spec_sim <- AddConfInterval(spec_sim, pval = nominal_p)
+#' true_spec <- alpha * spec_sim$freq^-beta
+#'
+#' LPlot(spec_sim)
+#' abline(a = log10(alpha), b = -beta, lty = 2, col = "red")
 #'
-#' q.empirical <- apply(save.spec, 1, quantile, c(0.025, 0.975))
-#' testspec <- SpecMTM(ts(SimPowerlaw(1, N.T)))
+#' SpecCoverage <- function(spec, true_spec){
 #'
-#' LPlot(AddConfInterval(testspec), ylim = c(0.05, 10))
-#' lines(testspec$freq, q.empirical[1, ], col = "red")
-#' lines(testspec$freq, q.empirical[2, ],col = "red")
-#' legend("bottomleft", lwd = 2, col = c("black", "red"),
-#'        legend = c("one realization with chisq conf intervals",
-#'                   "MC confidence intervals"))
+#'   stopifnot("lim.1" %in% names(spec))
+#'
+#'   n <- length(spec$freq)
+#'   below <- sum(spec$lim.1 < true_spec)
+#'   above <- sum(spec$lim.2 > true_spec)
+#'
+#'   total <- below + above
+#'
+#'   count <- list(n = n, above = above, below = below, total = total,
+#'                 p_above = above / n, p_below = below / n,
+#'                 p_total = total / n,
+#'                 nominal_pval = spec$pval)
+#'
+#'   as.data.frame(count)
+#' }
+#'
+#' SpecCoverage(spec_sim, true_spec)
 #'
 #' @export
 AddConfInterval <- function(spec, pval = 0.05, MINVALUE = 1e-10) {
 
   is.spectrum(spec)
 
-  spec$lim.1 <- spec$spec * qchisq(c(1 - pval / 2), spec$dof) / (spec$dof)
-  spec$lim.2 <- spec$spec * qchisq(c(pval / 2), spec$dof) / (spec$dof)
+  spec$lim.1 <- spec$spec * 1 / (qchisq(c(pval / 2), spec$dof) / (spec$dof))
+  spec$lim.2 <- spec$spec * 1 / (qchisq(c(1 - pval / 2), spec$dof) / (spec$dof))
+
   spec$lim.1[spec$lim.1 < MINVALUE] <- MINVALUE
   spec$lim.2[spec$lim.2 < MINVALUE] <- MINVALUE
 
diff --git a/R/FilterSpec.R b/R/FilterSpec.R
index a01914d..0212f8b 100644
--- a/R/FilterSpec.R
+++ b/R/FilterSpec.R
@@ -1,11 +1,14 @@
 #' Filter a Power Spectrum Object
 #'
 #' @param spec A spec object
+#' @param method single integer for choosing an endpoint constraint method.
+#' Available choices are integers 0-4, see details of \code{\link{ApplyFilter}}
 #' @param keep_low_f Keep filtered (smoothed) low frequencies or replace with unfiltered
 #' @inheritParams stats::spec.pgram
 #' @inheritParams ApplyFilter
 #' @return A spec object (list)
 #' @family functions to filter / smooth spectra
+#' @author Andrew Dolman
 #' @export
 #'
 #' @examples
@@ -57,34 +60,47 @@ FilterSpec <- function(spec, spans, method = 3, keep_low_f = TRUE) {
   kernel <- stats::kernel("modified.daniell", spans %/% 2)
   filter <- kernel[-kernel$m:kernel$m]
 
-  spec_filt <- ApplyFilter(spec$spec, filter = filter, method = method)
+  spec_filt <- ApplyFilter(spec$spec, filter = filter,
+                           method = method)
 
   if (keep_low_f == FALSE) {
     # replace filtered spec with original in area where freqs have been reflected
     i <- 1:ceiling(length(filter) / 2)
     spec_filt[i] <- spec$spec[i]
-
     iend <- length(spec$freq) - (i-1)
-
     spec_filt[iend] <- spec$spec[iend]
-
   }
 
   spec$spec <- as.numeric(spec_filt)
 
-  # degrees of freedom of the kernel
-  df.kern <- stats::df.kernel(kernel)
 
-  spec$dof <- df.kern * spec$dof / 2
+  ## New degrees of freedom of the filtered spec
+  ## This can be estimated using the convolution of the new filter and an
+  ## approximate equivalent filter to the existing DOF
 
-  if (keep_low_f == FALSE) {
+  # Approximate an equivalent filter to the existing DOF
+  fl <- round(mean(spec$dof)/2)
+  if ((fl%%2) == 0) fl <- fl-1 #ensure that the filter length is odd
+
+  prior_filter <- rep(1/fl, fl)
+
+  combined_filter <- convolve(filter, prior_filter, type = "open")
+
+  m3 <- floor(length(combined_filter)/2)
+  kernel3 <- combined_filter[1:ceiling(length(combined_filter)/2)]
+  kernel3 <- kernel(coef = sort(kernel3, decreasing = TRUE), m = m3)
+
+  new.dof <- stats::df.kernel(kernel3)
+
+  spec$dof <- rep(new.dof, length(spec$freq))
 
+
+  if (keep_low_f == FALSE) {
     i <- 1:ceiling(length(filter) / 2)
     spec$dof[i] <- dof0[i]
 
     iend <- length(spec$freq) - (i-1)
     spec$dof[iend] <- dof0[iend]
-
   }
 
   # Adjust DOF in reflected filter region
@@ -116,7 +132,6 @@ FilterSpec <- function(spec, spans, method = 3, keep_low_f = TRUE) {
 
   spec <- AddConfInterval(spec)
 
-
   return(spec)
 }
 
@@ -126,11 +141,12 @@ FilterSpec <- function(spec, spans, method = 3, keep_low_f = TRUE) {
 #'
 #' @param spec A spec object
 #' @inheritParams LogSmooth
-#' @inheritParams ApplyFilter
-#'
+#' @param method single integer for choosing an endpoint constraint method.
+#' Available choices are integers 0-4, see details of \code{\link{ApplyFilter}}
 #' @return A spec object (list)
 #' @family functions to filter / smooth spectra
 #' @export
+#' @author Andrew Dolman
 #' @examples
 #' library(PaleoSpec)
 #'
@@ -219,7 +235,7 @@ FilterSpecLog <- function(spec,
   dof0 <- spec$dof
 
   # Gets the difference in delta_f for the log and standard freq axis
-  NpF <- function(freq, fw, df){
+  NpF2 <- function(freq, fw, df){
 
     posdiff <- (exp(log(freq) + df) - freq)
     negdiff <- (freq - exp(log(freq) - df))
@@ -228,9 +244,12 @@ FilterSpecLog <- function(spec,
 
     2 * fw * (fdiff/df) * 1/(2*max(freq))
   }
-  df.logkern <- NpF(spec$freq, length(filter), df = diff(freq_logspace[1:2]))
 
-  spec$dof <- spec$dof + df.mod * df.logkern * spec$dof/2
+  deltaf.logkern <- NpF2(spec$freq, length(filter), df = diff(freq_logspace[1:2]))
+
+  # this would be more accurate as the dof of the convolution of the filter not sum of the old and new dofs
+  spec$dof <- df.mod * deltaf.logkern + spec$dof
+
   spec$shape <- spec$dof/2
   spec$spans <- paste(spans, collapse = ",")
 
@@ -244,7 +263,5 @@ FilterSpecLog <- function(spec,
 
   spec <- AddConfInterval(spec)
 
-
-
   return(spec)
 }
diff --git a/man/AddConfInterval.Rd b/man/AddConfInterval.Rd
index 09efe07..6ec4df9 100644
--- a/man/AddConfInterval.Rd
+++ b/man/AddConfInterval.Rd
@@ -23,30 +23,58 @@ the input object including the new list elements \code{lim.1} and
   respectively.
 }
 \description{
+Calculated confidence intervals give the interval within which the true
+spectrum should lie with frequency 1-p. This behaviour changed with version
+0.33.
+
+Previous versions returned the interval which, when applied to the 'true' spectrum, would contain new spectral estimates
+with frequency 1-p.
+}
+\details{
 Add a confidence interval to a spectral estimate assuming the uncertainty of
 the spectral estimates follows a chi-squared distribution with the degrees of
 the freedom from the spectral estimates.
 }
 \examples{
 
-N.R <- 1000
-N.T <- 100
-save.spec <- matrix(NA, N.T / 2, N.R)
-for (i.R in 1 : N.R) {
-  save.spec[, i.R] <- SpecMTM(ts(SimPowerlaw(1, N.T)))$spec
-}
+alpha <- 0.1
+beta <- 1
+
+spec_sim <- SpecMTM(ts(SimPLS(N = 1e03, beta = beta, alpha = alpha)))
+
+# Using a large nominal p value of 0.25 to reduce variation between random
+# timeseries
+
+# The true spec should be p/2 times above the lower CI, and p/2 below the upper CI
+
+nominal_p <- 0.25
+spec_sim <- AddConfInterval(spec_sim, pval = nominal_p)
+true_spec <- alpha * spec_sim$freq^-beta
 
-q.empirical <- apply(save.spec, 1, quantile, c(0.025, 0.975))
-testspec <- SpecMTM(ts(SimPowerlaw(1, N.T)))
+LPlot(spec_sim)
+abline(a = log10(alpha), b = -beta, lty = 2, col = "red")
+
+SpecCoverage <- function(spec, true_spec){
+
+  stopifnot("lim.1" \%in\% names(spec))
+
+  n <- length(spec$freq)
+  below <- sum(spec$lim.1 < true_spec)
+  above <- sum(spec$lim.2 > true_spec)
+
+  total <- below + above
+
+  count <- list(n = n, above = above, below = below, total = total,
+                p_above = above / n, p_below = below / n,
+                p_total = total / n,
+                nominal_pval = spec$pval)
+
+  as.data.frame(count)
+}
 
-LPlot(AddConfInterval(testspec), ylim = c(0.05, 10))
-lines(testspec$freq, q.empirical[1, ], col = "red")
-lines(testspec$freq, q.empirical[2, ],col = "red")
-legend("bottomleft", lwd = 2, col = c("black", "red"),
-       legend = c("one realization with chisq conf intervals",
-                  "MC confidence intervals"))
+SpecCoverage(spec_sim, true_spec)
 
 }
 \author{
-Thomas Laepple
+Thomas Laepple, Andrew Dolman
 }
diff --git a/man/FilterSpec.Rd b/man/FilterSpec.Rd
index e7c168c..61f93ba 100644
--- a/man/FilterSpec.Rd
+++ b/man/FilterSpec.Rd
@@ -12,8 +12,8 @@ FilterSpec(spec, spans, method = 3, keep_low_f = TRUE)
 \item{spans}{vector of odd integers giving the widths of modified
     Daniell smoothers to be used to smooth the periodogram.}
 
-\item{method}{single integer for choosing an endpoint constraint method;
-available choices are integers 0-4, see details.}
+\item{method}{single integer for choosing an endpoint constraint method.
+Available choices are integers 0-4, see details of \code{\link{ApplyFilter}}}
 
 \item{keep_low_f}{Keep filtered (smoothed) low frequencies or replace with unfiltered}
 }
@@ -64,4 +64,7 @@ Other functions to filter / smooth spectra:
 \code{\link{FilterSpecLog}()},
 \code{\link{LogSmooth}()}
 }
+\author{
+Andrew Dolman
+}
 \concept{functions to filter / smooth spectra}
diff --git a/man/FilterSpecLog.Rd b/man/FilterSpecLog.Rd
index a82e2eb..82f33a6 100644
--- a/man/FilterSpecLog.Rd
+++ b/man/FilterSpecLog.Rd
@@ -11,8 +11,8 @@ FilterSpecLog(spec, df.log = 0.05, spans = NULL, method = 3, f.res = 10)
 
 \item{df.log}{width of the smoother in log units}
 
-\item{method}{single integer for choosing an endpoint constraint method;
-available choices are integers 0-4, see details.}
+\item{method}{single integer for choosing an endpoint constraint method.
+Available choices are integers 0-4, see details of \code{\link{ApplyFilter}}}
 }
 \value{
 A spec object (list)
@@ -63,4 +63,7 @@ Other functions to filter / smooth spectra:
 \code{\link{FilterSpec}()},
 \code{\link{LogSmooth}()}
 }
+\author{
+Andrew Dolman
+}
 \concept{functions to filter / smooth spectra}
diff --git a/man/figures/README-unnamed-chunk-10-1.png b/man/figures/README-unnamed-chunk-10-1.png
index 4561da1..30f544a 100644
Binary files a/man/figures/README-unnamed-chunk-10-1.png and b/man/figures/README-unnamed-chunk-10-1.png differ
diff --git a/man/figures/README-unnamed-chunk-5-1.png b/man/figures/README-unnamed-chunk-5-1.png
index bd6a23e..3209ba7 100644
Binary files a/man/figures/README-unnamed-chunk-5-1.png and b/man/figures/README-unnamed-chunk-5-1.png differ
diff --git a/man/figures/README-unnamed-chunk-6-1.png b/man/figures/README-unnamed-chunk-6-1.png
index 79f9046..650aaf0 100644
Binary files a/man/figures/README-unnamed-chunk-6-1.png and b/man/figures/README-unnamed-chunk-6-1.png differ
diff --git a/man/figures/README-unnamed-chunk-9-1.png b/man/figures/README-unnamed-chunk-9-1.png
index 70409b2..adde0d8 100644
Binary files a/man/figures/README-unnamed-chunk-9-1.png and b/man/figures/README-unnamed-chunk-9-1.png differ
diff --git a/tests/testthat/test-AddConfInterval.R b/tests/testthat/test-AddConfInterval.R
index bb8fdb9..e1451c4 100644
--- a/tests/testthat/test-AddConfInterval.R
+++ b/tests/testthat/test-AddConfInterval.R
@@ -10,8 +10,10 @@ test_that("adding the confidence interval is correct.", {
 
   pval <- 0.05
 
-  lim1 <- spec * qchisq(1 - pval / 2, dof) / dof
-  lim2 <- spec * qchisq(pval / 2, dof) / dof
+  lim1 <- spec * 1 / (qchisq(c(pval / 2), dof) / (dof))
+  lim2 <- spec * 1 / (qchisq(c(1 - pval / 2), dof) / (dof))
+
+
   lim1[lim1 < MINVALUE] <- MINVALUE
   lim2[lim2 < MINVALUE] <- MINVALUE
 
@@ -23,8 +25,9 @@ test_that("adding the confidence interval is correct.", {
   MINVALUE <- 1.
   pval <- 0.1
 
-  lim1 <- spec * qchisq(1 - pval / 2, dof) / dof
-  lim2 <- spec * qchisq(pval / 2, dof) / dof
+  lim1 <- spec * 1 / (qchisq(c(pval / 2), dof) / (dof))
+  lim2 <- spec * 1 / (qchisq(c(1 - pval / 2), dof) / (dof))
+
   lim1[lim1 < MINVALUE] <- MINVALUE
   lim2[lim2 < MINVALUE] <- MINVALUE
 
diff --git a/tests/testthat/test-gg_spec.R b/tests/testthat/test-gg_spec.R
index efe38d8..899deaf 100644
--- a/tests/testthat/test-gg_spec.R
+++ b/tests/testthat/test-gg_spec.R
@@ -11,7 +11,7 @@ sp1 <- SpecMTM(ts1)
 sp2 <- SpecMTM(ts2)
 
 # 1 as spec_df
-sp1_df <- as.data.frame(sp1)
+sp1_df <- Spec2DF(sp1)
 
 # list of spectra, 1 is a spec_df
 sp_lst <- list(sp2 = sp2, sp1_df = sp1_df)
@@ -64,3 +64,4 @@ test_that("removeFirst and removeLast operate on a per spec_id basis", {
   expect_that(nrow(g0$data) - nrow(gl3$data), equals(6))
 
 })
+