diff --git a/CRAN-SUBMISSION b/CRAN-SUBMISSION
deleted file mode 100644
index b61d8e6..0000000
--- a/CRAN-SUBMISSION
+++ /dev/null
@@ -1,3 +0,0 @@
-Version: 0.4.7
-Date: 2025-01-14 02:18:38 UTC
-SHA: 17a1f29170c6dddf83f8daa57344c5cfc106673a
diff --git a/NEWS.md b/NEWS.md
index e193c58..9e6981d 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -24,8 +24,6 @@
 
 -   Add R-hub workflow to enable checks for platforms not covered already.
 
--   Vignettes now force `Sys.setenv("OMP_THREAD_LIMIT" = 1)` and `options(Ncpus = 1)` in a bid to get around CRAN CMD changes.
-
 # joineRML 0.4.6
 
 ## Housekeeping
diff --git a/cran-comments.md b/cran-comments.md
index 9e44c78..cbf80c6 100644
--- a/cran-comments.md
+++ b/cran-comments.md
@@ -13,7 +13,7 @@
 
 ## R CMD check results
 
-0 errors | 0 warnings | 2 notes
+0 errors | 0 warnings | 1 notes
 
 Win-Builder NOTE: "checking CRAN incoming feasibility ... NOTE"
 
diff --git a/vignettes/joineRML-tidy.Rmd b/vignettes/joineRML-tidy.Rmd
index 7681c47..ae1f3b2 100644
--- a/vignettes/joineRML-tidy.Rmd
+++ b/vignettes/joineRML-tidy.Rmd
@@ -17,6 +17,7 @@ library(knitr)
 Sys.setenv("OMP_THREAD_LIMIT" = 1)
 Sys.setenv("OMP_NUM_THREADS" = 1)
 options(Ncpus = 1)
+options(cores = 2)
 ```
 
 # Introduction
@@ -41,7 +42,7 @@ These methods are specifically useful when plotting results of a joint model or
 
 We use the sample example from the introductory vignette to `joineRML` using the heart valve data.
 
-```{r vignette, eval=FALSE}
+```{r vignette, eval=FALSE, purl=FALSE}
 vignette("joineRML", package = "joineRML")
 help("heart.valve", package = "joineRML")
 ```
@@ -99,7 +100,7 @@ tidy(fit, ci = TRUE, conf.level = 0.99)
 
 The standard errors reported by default are based on the empirical information matrix, as in `mjoint`. It is of course possible to use bootstrapped standard errors as follows:
 
-```{r tidy-boot, eval=FALSE}
+```{r tidy-boot, eval=FALSE, purl=FALSE}
 bSE <- bootSE(fit, nboot = 100, safe.boot = TRUE, progress = FALSE)
 tidy(fit, boot_se = bSE, conf.int = TRUE)
 ```
@@ -177,6 +178,6 @@ glance(fit2)
 
 Several examples of how to use `broom` including more details are available on its introductory vignette:
 
-```{r vignette-broom, eval = FALSE}
+```{r vignette-broom, eval=FALSE, purl=FALSE}
 vignette(topic = "broom", package = "broom")
 ```
diff --git a/vignettes/joineRML.Rmd b/vignettes/joineRML.Rmd
index 3779ccf..6c74a2d 100644
--- a/vignettes/joineRML.Rmd
+++ b/vignettes/joineRML.Rmd
@@ -25,13 +25,14 @@ if (requireNamespace('joineR', quietly = TRUE)) {
 Sys.setenv("OMP_THREAD_LIMIT" = 1)
 Sys.setenv("OMP_NUM_THREADS" = 1)
 options(Ncpus = 1)
+options(cores = 2)
 ```
 
 # Introduction
 
 The `joineRML` package implements methods for analyzing data from *multiple* longitudinal studies in which the responses from each subject consists of time-sequences of repeated measurements and a possibly censored time-to-event outcome. The modelling framework for the repeated measurements is the multivariate linear mixed effects model. The model for the time-to-event outcome is a Cox proportional hazards model with log-Gaussian frailty. Stochastic dependence is captured by allowing the Gaussian random effects of the linear model to be correlated with the frailty term of the Cox proportional hazards model. For full details of the model, please consult the technical vignette by running
 
-```{r vignette, eval=FALSE}
+```{r vignette, eval=FALSE, purl=FALSE}
 vignette("technical", package = "joineRML")
 ```
 
@@ -41,7 +42,7 @@ vignette("technical", package = "joineRML")
 
 The simplest way to explain the concepts of the package is through an example. `joineRML` comes with the data set `heart.valve`. Details of this data can be found in the help file by running the command
 
-```{r heart.valve_help, eval=FALSE}
+```{r heart.valve_help, eval=FALSE, purl=FALSE}
 help("heart.valve", package = "joineRML")
 ```
 
@@ -133,7 +134,7 @@ plot(fit, params = "beta")
 
 Once an `mjoint` model has converged, and assuming the `pfs` argument is `TRUE` (default), then approximated standard errors are calculated based on the empirical information matrix of the profile likelihood at the maximizer. Theoretically, these standard errors will be underestimated (see the technical vignette). In principle, residual Monte Carlo error will oppose this through an increase in uncertainty.
 
-```{r hvd_model_boot, eval=FALSE}
+```{r hvd_model_boot, eval=FALSE, purl=FALSE}
 fit.se <- bootSE(fit, nboot = 100)
 ```
 
@@ -141,13 +142,13 @@ Bootstrapping is a computationally intensive method, possibly taking many hours
 
 We can call the `bootSE` object to interrogate it
 
-```{r hvd_model_boot_print, eval=FALSE}
+```{r hvd_model_boot_print, eval=FALSE, purl=FALSE}
 fit.se
 ```
 
 or alternatively re-run the `summary` command, passing the additional argument of `bootSE = fit.se`
 
-```{r hvd_model_boot_summary, eval=FALSE}
+```{r hvd_model_boot_summary, eval=FALSE, purl=FALSE}
 summary(fit, bootSE = fit.se)
 ```
 
@@ -171,7 +172,7 @@ To fit a univariate model in `joineR` we run the following code for the `hvd` da
 joineR_available <- require(joineR)
 ```
 
-```{r joineR, cache=TRUE, eval=joineR_available}
+```{r joineR, cache=TRUE, eval=joineR_available, purl=joineR_available}
 library(joineR, quietly = TRUE)
 
 hvd.surv <- UniqueVariables(hvd, var.col = c("fuyrs", "status"), id.col = "num")
@@ -207,7 +208,7 @@ summary(fit.joinerml)
 
 In addition to just comparing model parameter estimates, we can also extract the predicted (or posterior) random effects from each model and plot them.
 
-```{r re_comp_plot, fig.width=7.25, fig.height=4, eval=joineR_available}
+```{r re_comp_plot, fig.width=7.25, fig.height=4, eval=joineR_available, purl=joineR_available}
 id <- as.numeric(row.names(fit.joiner$coefficients$random))
 id.ord <- order(id) # joineR rearranges patient ordering during EM fit
 par(mfrow = c(1, 2))