diff --git a/DESCRIPTION b/DESCRIPTION
index ec365a752a2c5a6d307e044c7305be753b280e7d..235db095aa05e28511ca813b881ee8fb3a074b7e 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,7 +1,7 @@
 Package: onlineforecast
 Type: Package
 Title: Forecast Modelling for Online Applications
-Version: 0.9.4
+Version: 0.10.0
 Description: A framework for fitting adaptive forecasting models. Provides a way to use forecasts as input to models, e.g. weather forecasts for energy related forecasting. The models can be fitted recursively and can easily be setup for updating parameters when new data arrives. See the included vignettes, the website <https://onlineforecasting.org> and the paper "Short-term heat load forecasting for single family houses" <doi:10.1016/j.enbuild.2013.04.022>.
 License: GPL-3
 Encoding: UTF-8
diff --git a/R/as.data.list.R b/R/as.data.list.R
index adeca1c8b6dc47cea32f3c3be674008da55f58f0..202d711927db6d727b6b89860f59f2e89432e8d6 100644
--- a/R/as.data.list.R
+++ b/R/as.data.list.R
@@ -21,7 +21,7 @@
 #' @title Convert to data.list class
 #' @param object The object to be converted into a data.list
 #' @return a value of class data.list
-#' @seealso \code{For specific detailed info see the children, e.g. \link{onlinefocast:::as.data.list.data.frame} }
+#' @seealso \code{For specific detailed info see the children, e.g. \link{as.data.list.data.frame} }
 #' 
 #' @rdname as.data.list
 #' @export
diff --git a/R/in_range.R b/R/in_range.R
index 4b9e72ccb144fe2b043fc8835ba91f1c17aa9781..b09d0910383ed0086c991374093506adb83a55e8 100644
--- a/R/in_range.R
+++ b/R/in_range.R
@@ -19,7 +19,6 @@
 #' @param tstart The start of the period.
 #' @param time The timestamps as POSIX.
 #' @param tend The end of the period. If not given then the period will have no end. 
-#' @param timezone The timezone of the timestamps, time.
 #' @return A logical vector indicating the selected period with TRUE
 #' @name in_range
 #' @examples
diff --git a/R/lagdf.R b/R/lagdf.R
index e3d326965318edde63620ee0af5c9ca3df15b613..a72012cdf0a3ef4911a706a1609bd2a023afd565 100644
--- a/R/lagdf.R
+++ b/R/lagdf.R
@@ -22,9 +22,9 @@
 #' lagdf(1:10, 3)
 #' # Back in time
 #' lagdf(1:10, -3)
-#' # Works but returns a numric
+#' # Works but returns a numeric column
 #' lagdf(as.factor(1:10), 3)
-#' # Works and returns a character
+#' # Works and returns a character column
 #' lagdf(as.character(1:10), 3)
 #' # Giving several lag values
 #' lagdf(1:10, c(1:3))
diff --git a/R/lapply.R b/R/lapply.R
index 9bffb1ca301fbf39cd1c73071be9f5e0d48e9cf5..6fb34ac169ffdbfd63e97e8fe85dcefc5cc2e45b 100644
--- a/R/lapply.R
+++ b/R/lapply.R
@@ -7,6 +7,7 @@
 #' Helper which does lapply and then cbind
 #' @param X object to apply on
 #' @param FUN function to apply
+#' @param ... passed on to lapply
 #' @export
 lapply_cbind <- function(X, FUN, ...){
   val <- lapply(X, FUN, ...)
@@ -16,6 +17,7 @@ lapply_cbind <- function(X, FUN, ...){
 #' Helper which does lapply and then rbind
 #' @param X object to apply on
 #' @param FUN function to apply
+#' @param ... passed on to lapply
 #' @export
 lapply_rbind <- function(X, FUN, ...){
   val <- lapply(X, FUN, ...)
@@ -25,6 +27,7 @@ lapply_rbind <- function(X, FUN, ...){
 #' Helper which does lapply, cbind and then as.data.frame
 #' @param X object to apply on
 #' @param FUN function to apply
+#' @param ... passed on to lapply
 #' @export
 lapply_cbind_df <- function(X, FUN, ...){
   val <- lapply(X, FUN, ...)
@@ -34,6 +37,7 @@ lapply_cbind_df <- function(X, FUN, ...){
 #' Helper which does lapply, rbind and then as.data.frame
 #' @param X object to apply on
 #' @param FUN function to apply
+#' @param ... passed on to lapply
 #' @export
 lapply_rbind_df <- function(X, FUN, ...){
   val <- lapply(X, FUN, ...)
diff --git a/R/step_optim.R b/R/step_optim.R
index a0b178831153bcfdee4b56df1249bc798b67bb15..971e2654f9cbbf5ddd3a9dfd1010708b3e4bf265 100644
--- a/R/step_optim.R
+++ b/R/step_optim.R
@@ -55,11 +55,13 @@
 #' the score will be calculated using only the complete cases across horizons
 #' and models in each step, see the last examples.
 #'
+#' Note, that either kseq or kseqopt must be set on the modelfull object. If kseqopt
+#' is set, then it is used no matter the value of kseq.
+#'
 #' @title Forward and backward model selection
 #' @param modelfull The full forecastmodel containing all inputs which will be
 #'     can be included in the selection.
 #' @param data The data.list which holds the data on which the model is fitted.
-#' @param kseq The horizons to fit for (if not set, then model$kseq is used)
 #' @param prm A list of integer parameters to be stepped. Given using the same
 #'     syntax as parameters for optimization, e.g. `list(I__degree = c(min=3,
 #'     max=7))` will step the "degree" for input "I".
@@ -78,6 +80,7 @@
 #'     in rls_optim()). Furthermore, information on complete cases are printed
 #'     and returned.
 #' @param scorefun The score function used.
+#' @param printout Logical. Passed on to fitting functions.
 #' @param mc.cores The mc.cores argument of mclapply. If debugging it can be
 #'     nessecary to set it to 1 to stop execution.
 #' @param ... Additional arguments which will be passed on to optimfun. For
@@ -117,7 +120,8 @@
 #' model$add_prmbounds(lambda = c(0.9, 0.99, 0.9999))
 #' 
 #' # Select a model, in the optimization just run it for a single horizon
-#' model$kseqopt <- 5
+#' # Note that kseqopt could also be set
+#' model$kseq <- 5
 #' # 
 #' prm <- list(mu_tday__nharmonics = c(min=3, max=7))
 #' 
@@ -125,7 +129,8 @@
 #' # Iterations in the prm optimization (MUST be increased in real applications)
 #' control <- list(maxit=1)
 #'
-#' # Run the default selection scheme, which is "both" (same as "backwardboth" if no start model is given)
+#' # Run the default selection scheme, which is "both"
+#' # (same as "backwardboth" if no start model is given)
 #' L <- step_optim(model, D, prm, control=control)
 #'
 #' # The optim value from each step is returned
@@ -142,7 +147,7 @@
 #' Lforwardboth <- step_optim(model, D, prm, "forwardboth", control=control, mc.cores=1)
 #'
 #' # It's possible avoid removing specified inputs
-#' L <- step_optim(model, D, prm, keepinputs = c("mu","mu_tday"), control=control)
+#' L <- step_optim(model, D, prm, keepinputs=c("mu","mu_tday"), control=control)
 #' 
 #' # Give a starting model
 #' modelstart <- model$clone_deep()
@@ -172,7 +177,7 @@
 #'
 #' @export
 
-step_optim <- function(modelfull, data, prm=list(NA), kseq = NA, direction = c("both","backward","forward","backwardboth","forwardboth"), modelstart=NA, keepinputs = FALSE, optimfun = rls_optim, fitfun = NA, scorefun = rmse, printout = FALSE, mc.cores = getOption("mc.cores", 2L), ...){
+step_optim <- function(modelfull, data, prm=list(NA), direction = c("both","backward","forward","backwardboth","forwardboth"), modelstart=NA, keepinputs = FALSE, optimfun = rls_optim, fitfun = NA, scorefun = rmse, printout = FALSE, mc.cores = getOption("mc.cores", 2L), ...){
     # Do:
     # - checking of input, model, ...
     # - Maybe have "cloneit" argument in optimfun, then don't clone inside optim.
@@ -191,7 +196,7 @@ step_optim <- function(modelfull, data, prm=list(NA), kseq = NA, direction = c("
     istep <- 1
     # Different start up, if a start model is given
     if( class(modelstart)[1] == "forecastmodel" ){
-        # The full model will not be changed from here, so don't need to clone it
+        # The full model will not be changed from here, so no need to clone it
         mfull <- modelfull
         m <- modelstart$clone()
     }else{
@@ -227,6 +232,10 @@ step_optim <- function(modelfull, data, prm=list(NA), kseq = NA, direction = c("
             scoreCurrent <- Inf
         }
     }
+    # If kseqopt is set, then make sure that it is used in all runs (also when only running fitfun)
+    if(!is.na(m$kseqopt)){
+        m$kseq <- m$kseqopt
+    }
     # Find the inputs to keep, if any
     if(class(keepinputs) == "logical"){
         if(keepinputs){
@@ -250,14 +259,13 @@ step_optim <- function(modelfull, data, prm=list(NA), kseq = NA, direction = c("
         # If the init model is not yet optimized
         if(istep == 1 & length(L) == 0){
             # Optimize
-            res <- optimfun(m, data, kseq, printout=printout, scorefun=scorefun, ...)
+            res <- optimfun(m, data, printout=printout, scorefun=scorefun, ...)
             # Should we forecast only on the complete cases?
             if(class(fitfun) == "function"){
                 # Forecast to get the complete cases
                 mtmp <- m$clone_deep()
-                mtmp$kseq <- kseq
                 Yhat <- fitfun(res$par, mtmp, data, printout=printout)$Yhat
-                scoreCurrent <- sum(score(residuals(Yhat,data[[model$output]]),data$scoreperiod))
+                scoreCurrent <- sum(score(residuals(Yhat,data[[m$output]]),data$scoreperiod))
                 casesCurrent <- complete_cases(Yhat)
             }else{
                 scoreCurrent <- res$value
@@ -360,7 +368,7 @@ step_optim <- function(modelfull, data, prm=list(NA), kseq = NA, direction = c("
 
             # Run the optimization
             Lstep <- mclapply(1:length(mStep), function(i, ...){
-                optimfun(mStep[[i]], data, kseq, printout=printout, ...)
+                optimfun(mStep[[i]], data, printout=printout, scorefun=scorefun, ...)
             }, mc.cores=mc.cores, ...)
             names(Lstep) <- names(mStep)
 
@@ -368,11 +376,10 @@ step_optim <- function(modelfull, data, prm=list(NA), kseq = NA, direction = c("
             if(class(fitfun) == "function"){
                 LYhat <- mclapply(1:length(mStep), function(i){
                     mtmp <- mStep[[i]]$clone_deep()
-                    mtmp$kseq <- kseq
                     fitfun(Lstep[[i]]$par, mtmp, data, printout=printout)$Yhat
                 }, mc.cores=mc.cores)
                 # Use complete cases across models and horizons per default
-                scoreStep <- apply(score(residuals(LYhat,data[[model$output]]), data$scoreperiod), 2, sum)
+                scoreStep <- apply(score(residuals(LYhat,data[[m$output]]), data$scoreperiod), 2, sum)
                 casesStep <- sapply(LYhat, complete_cases)
             }else{
                 # Use the scores from optimfun
@@ -393,7 +400,7 @@ step_optim <- function(modelfull, data, prm=list(NA), kseq = NA, direction = c("
                 tmp <- cbind(tmp, apply(casesStep != casesCurrent, 2, sum))
                 nams(tmp)[2] <- "CasesDiff"
             }
-            onlineforecast:::print_to_message(tmp)
+            print_to_message(tmp)
 
             # Compare scores: Is one the step models score smaller than the current ref?
             imin <- which.min(scoreStep)
diff --git a/cran-comments.md b/cran-comments.md
index a3024177f22c1ea47ce7f56edaff52697a199645..437dde459efb458e3fc75cff34c6d713f5e318bb 100644
--- a/cran-comments.md
+++ b/cran-comments.md
@@ -1,3 +1,15 @@
+#----------------------------------------------------------------
+# v0.10.0
+
+We have added features and done some small changes to functions. This version
+should be fully backward compatible.
+
+
+
+
+
+
+
 #----------------------------------------------------------------
 # v0.9.3
 # Response to review of v0.9.2 by Uwe Ligges
diff --git a/make.R b/make.R
index 70cc8b3c56bf3aaf8fdf5846758e0f2018511a50..28a07b7e4bc77291ce683ab7317def2def49b5e1 100644
--- a/make.R
+++ b/make.R
@@ -50,8 +50,8 @@ library(roxygen2)
 #test()
 
 # # Run the examples
-#load_all(as.package("../onlineforecast"))
-#run_examples()
+load_all(as.package("../onlineforecast"))
+run_examples()
 
 # # Run tests in a single file
 # test_file("tests/testthat/test-rls-heat-load.R")
@@ -60,10 +60,10 @@ library(roxygen2)
 # ----------------------------------------------------------------
 # Build the package
 document()
-build(".", vignettes=FALSE)
+build(".", vignettes=TRUE)
 
 # Install it
-install.packages("../onlineforecast_0.9.4.tar.gz")
+install.packages("../onlineforecast_0.10.0.tar.gz")
 
 library(onlineforecast)
 # ----------------------------------------------------------------
@@ -80,11 +80,11 @@ library(onlineforecast)
 # Test before release
 devtools::check()
 
-devtools::check_built("../onlineforecast_0.9.4.tar.gz")
+devtools::check_built("../onlineforecast_0.10.0.tar.gz")
 
 # Does give different results than check() above
 #system("R CMD check --as-cran ../onlineforecast_0.9.4.tar.gz")
-system("R CMD check ../onlineforecast_0.9.4.tar.gz")
+system("R CMD check ../onlineforecast_0.10.0.tar.gz")
 unlink("onlineforecast.Rcheck/", recursive=TRUE)
 
 # Use for more checking:
diff --git a/man/as.data.list.Rd b/man/as.data.list.Rd
index 14c29280e47b5f0b6ffadc4fe7b83652ee5262cf..9567fef28f93a2771f8bb7d8b6d642e7edbbb5c0 100644
--- a/man/as.data.list.Rd
+++ b/man/as.data.list.Rd
@@ -51,7 +51,7 @@ as.data.frame(as.data.list(X))
 
 }
 \seealso{
-\code{For specific detailed info see the children, e.g. \link{onlinefocast:::as.data.list.data.frame} }
+\code{For specific detailed info see the children, e.g. \link{as.data.list.data.frame} }
 
 as.data.list
 }
diff --git a/man/in_range.Rd b/man/in_range.Rd
index 658c200cea1dc1841fdb171bc2fbe16f12a3cc7e..c3412162da3a75c61bfd95ad9454014cbdfcf2ad 100644
--- a/man/in_range.Rd
+++ b/man/in_range.Rd
@@ -12,8 +12,6 @@ in_range(tstart, time, tend = NA)
 \item{time}{The timestamps as POSIX.}
 
 \item{tend}{The end of the period. If not given then the period will have no end.}
-
-\item{timezone}{The timezone of the timestamps, time.}
 }
 \value{
 A logical vector indicating the selected period with TRUE
diff --git a/man/lagdf.Rd b/man/lagdf.Rd
index a072d02c737abdebae65a12f312e0c48ec0a82f0..c13d5f56e47f5ae3443c01c6765bc1be5ca58494 100644
--- a/man/lagdf.Rd
+++ b/man/lagdf.Rd
@@ -36,9 +36,9 @@ This function lags the columns with the integer values specified with the argume
 lagdf(1:10, 3)
 # Back in time
 lagdf(1:10, -3)
-# Works but returns a numric
+# Works but returns a numeric column
 lagdf(as.factor(1:10), 3)
-# Works and returns a character
+# Works and returns a character column
 lagdf(as.character(1:10), 3)
 # Giving several lag values
 lagdf(1:10, c(1:3))
diff --git a/man/lapply_cbind.Rd b/man/lapply_cbind.Rd
index fd57b1d897ad439d3c7b501aaeadeba40b618207..f075fcc2e0fe7cbeb82d85b57bf7c0d788129ff2 100644
--- a/man/lapply_cbind.Rd
+++ b/man/lapply_cbind.Rd
@@ -10,6 +10,8 @@ lapply_cbind(X, FUN, ...)
 \item{X}{object to apply on}
 
 \item{FUN}{function to apply}
+
+\item{...}{passed on to lapply}
 }
 \description{
 Helper which does lapply and then cbind
diff --git a/man/lapply_cbind_df.Rd b/man/lapply_cbind_df.Rd
index 5b8e5a65d5037304f935647d8b50f3aa37b9a48a..4801b513be67a328f52a51bf39d871a8e200bb0a 100644
--- a/man/lapply_cbind_df.Rd
+++ b/man/lapply_cbind_df.Rd
@@ -10,6 +10,8 @@ lapply_cbind_df(X, FUN, ...)
 \item{X}{object to apply on}
 
 \item{FUN}{function to apply}
+
+\item{...}{passed on to lapply}
 }
 \description{
 Helper which does lapply, cbind and then as.data.frame
diff --git a/man/lapply_rbind.Rd b/man/lapply_rbind.Rd
index a4e3939e812ea59b57321d55633cfd4804d87b89..55d51afe182edc69c1459c6015eb33fdb81aa3e6 100644
--- a/man/lapply_rbind.Rd
+++ b/man/lapply_rbind.Rd
@@ -10,6 +10,8 @@ lapply_rbind(X, FUN, ...)
 \item{X}{object to apply on}
 
 \item{FUN}{function to apply}
+
+\item{...}{passed on to lapply}
 }
 \description{
 Helper which does lapply and then rbind
diff --git a/man/lapply_rbind_df.Rd b/man/lapply_rbind_df.Rd
index 83a114a5e8e30b56514ecd0980bdd1f7b1612227..3bd25abacc42a80e2f4e3097a8efca19c307689c 100644
--- a/man/lapply_rbind_df.Rd
+++ b/man/lapply_rbind_df.Rd
@@ -10,6 +10,8 @@ lapply_rbind_df(X, FUN, ...)
 \item{X}{object to apply on}
 
 \item{FUN}{function to apply}
+
+\item{...}{passed on to lapply}
 }
 \description{
 Helper which does lapply, rbind and then as.data.frame
diff --git a/man/step_optim.Rd b/man/step_optim.Rd
index a72de690f6fa0d7a8e853a36c131929a7b8c24c6..58e258c0f2977ba7edfe0b6b115e1ca28e536fc6 100644
--- a/man/step_optim.Rd
+++ b/man/step_optim.Rd
@@ -8,7 +8,6 @@ step_optim(
   modelfull,
   data,
   prm = list(NA),
-  kseq = NA,
   direction = c("both", "backward", "forward", "backwardboth", "forwardboth"),
   modelstart = NA,
   keepinputs = FALSE,
@@ -30,8 +29,6 @@ can be included in the selection.}
 syntax as parameters for optimization, e.g. `list(I__degree = c(min=3,
 max=7))` will step the "degree" for input "I".}
 
-\item{kseq}{The horizons to fit for (if not set, then model$kseq is used)}
-
 \item{direction}{The direction to be used in the selection process.}
 
 \item{modelstart}{A forecastmodel. If it's set then it will be used as the
@@ -53,6 +50,8 @@ and returned.}
 
 \item{scorefun}{The score function used.}
 
+\item{printout}{Logical. Passed on to fitting functions.}
+
 \item{mc.cores}{The mc.cores argument of mclapply. If debugging it can be
 nessecary to set it to 1 to stop execution.}
 
@@ -118,6 +117,9 @@ models, it can be very important to make sure that only complete cases are
 included when calculating the score. By providing the `fitfun` argument then
 the score will be calculated using only the complete cases across horizons
 and models in each step, see the last examples.
+
+Note, that either kseq or kseqopt must be set on the modelfull object. If kseqopt
+is set, then it is used no matter the value of kseq.
 }
 \examples{
 
@@ -147,7 +149,8 @@ model$add_regprm("rls_prm(lambda=0.9)")
 model$add_prmbounds(lambda = c(0.9, 0.99, 0.9999))
 
 # Select a model, in the optimization just run it for a single horizon
-model$kseqopt <- 5
+# Note that kseqopt could also be set
+model$kseq <- 5
 # 
 prm <- list(mu_tday__nharmonics = c(min=3, max=7))
 
@@ -155,7 +158,8 @@ prm <- list(mu_tday__nharmonics = c(min=3, max=7))
 # Iterations in the prm optimization (MUST be increased in real applications)
 control <- list(maxit=1)
 
-# Run the default selection scheme, which is "both" (same as "backwardboth" if no start model is given)
+# Run the default selection scheme, which is "both"
+# (same as "backwardboth" if no start model is given)
 L <- step_optim(model, D, prm, control=control)
 
 # The optim value from each step is returned
@@ -172,7 +176,7 @@ Lbackwardboth <- step_optim(model, D, prm, "backwardboth", control=control)
 Lforwardboth <- step_optim(model, D, prm, "forwardboth", control=control, mc.cores=1)
 
 # It's possible avoid removing specified inputs
-L <- step_optim(model, D, prm, keepinputs = c("mu","mu_tday"), control=control)
+L <- step_optim(model, D, prm, keepinputs=c("mu","mu_tday"), control=control)
 
 # Give a starting model
 modelstart <- model$clone_deep()
diff --git a/vignettes/forecast-evaluation.Rmd b/vignettes/forecast-evaluation.Rmd
index 33f3f810eed224dc4d3fbbeea7266060c99764fa..1b44e287aec1d5c3985f523cd6c2c1b0249d9d62 100644
--- a/vignettes/forecast-evaluation.Rmd
+++ b/vignettes/forecast-evaluation.Rmd
@@ -7,7 +7,7 @@ output:
     toc: true
     toc_debth: 3
 vignette: >
-  %\VignetteIndexEntry{Online updating of onlineforecast models}
+  %\VignetteIndexEntry{Forecast evaluation}
   %\VignetteEngine{knitr::rmarkdown}
   %\VignetteEncoding{UTF-8}
 ---
diff --git a/vignettes/model-selection.Rmd b/vignettes/model-selection.Rmd
index 2d7e0cbdf029a5d801177e68069fd6bf440780b2..2fc3e8e0ffef498b39e64fa1b9f11335575fc822 100644
--- a/vignettes/model-selection.Rmd
+++ b/vignettes/model-selection.Rmd
@@ -150,8 +150,8 @@ model$add_prmbounds(lambda = c(0.9, 0.99, 0.9999))
 ```
 Finally, set the horizons to run (just keep a vector for later use):
 ```{r}
-# Select a model, just run it for a single horizon
-kseq <- 5
+# Select a model, just run optimization and score for a single horizon
+model$kseq <- 5
 ```
 
 Now we can use the `step_optim()` function for the selection.  In each step new models are generated, with either one removed input or one added input, and     then all the generated models are optimized and their scores compared. If     any new model have an improved score compared to the currently selected     model, then the new is selected and the process is repeated until no new     improvement is achieved.
@@ -184,7 +184,7 @@ The default procedure is backward selection with stepping in both directions:
 # Run the default selection, which is "both" and equivalent to "backwadboth"
 # Note the control argument, which is passed to optim, it's now set to few
 # iterations in the prm optimization
-Lboth <- step_optim(model, D, kseq, prm, direction="both", control=list(maxit=1))
+Lboth <- step_optim(model, D, prm, direction="both", control=list(maxit=1))
 ```
 We now have the models selected in each step in and we see that the final model
 is decreased:
@@ -194,14 +194,14 @@ getse(Lboth, "model")
 
 Forward selection:
 ```{r, message=FALSE, results="hide"}
-Lforward <- step_optim(model, D, kseq, prm, "forward", control=list(maxit=1))
+Lforward <- step_optim(model, D, prm, "forward", control=list(maxit=1))
 ```
 ```{r}
 getse(Lforward, "model")
 ```
 Same model is selected, which is also the case in backwards selection:
 ```{r, message=FALSE, results="hide"}
-Lbackward <- step_optim(model, D, kseq, prm, "backward", control=list(maxit=1))
+Lbackward <- step_optim(model, D, prm, "backward", control=list(maxit=1))
 ```
 ```{r}
 getse(Lbackward, "model")
@@ -215,7 +215,7 @@ modelstart <- model$clone_deep()
 # Remove two inputs
 modelstart$inputs[2:3] <- NULL
 # Run the selection
-L <- step_optim(model, D, kseq, prm, modelstart=modelstart, control=control)
+L <- step_optim(model, D, prm, modelstart=modelstart, control=list(maxit=1))
 ```
 ```{r}
 getse(L, "model")