From 36657396bfc37cf5853c7729d7576ef717fb4cc1 Mon Sep 17 00:00:00 2001 From: Andrew Johnson Date: Tue, 14 Mar 2023 12:26:17 +0200 Subject: [PATCH 01/13] Update structure for rstantools --- .Rbuildignore | 2 + .github/.gitignore | 1 + .github/workflows/R-CMD-check.yaml | 50 +++++++++++++++++++ R/stanmodels.R | 37 -------------- R/zzz.R | 4 -- configure | 3 +- configure.win | 3 +- inst/include/meta_header.hpp | 1 - {src/stan_files => inst/stan}/color.stan | 0 .../chunks => inst/stan/include}/license.stan | 0 {src/stan_files => inst/stan}/linear.stan | 0 .../stan}/reaction_time.stan | 0 .../stan}/success_rate.stan | 0 {src/stan_files => inst/stan}/ttest.stan | 0 src/Makevars | 26 ---------- src/Makevars.win | 26 ---------- src/init.cpp | 20 -------- tools/make_cc.R | 48 ------------------ 18 files changed, 55 insertions(+), 166 deletions(-) create mode 100644 .github/.gitignore create mode 100644 .github/workflows/R-CMD-check.yaml delete mode 100644 R/stanmodels.R delete mode 100644 R/zzz.R delete mode 100644 inst/include/meta_header.hpp rename {src/stan_files => inst/stan}/color.stan (100%) rename {src/stan_files/chunks => inst/stan/include}/license.stan (100%) rename {src/stan_files => inst/stan}/linear.stan (100%) rename {src/stan_files => inst/stan}/reaction_time.stan (100%) rename {src/stan_files => inst/stan}/success_rate.stan (100%) rename {src/stan_files => inst/stan}/ttest.stan (100%) delete mode 100644 src/Makevars delete mode 100644 src/Makevars.win delete mode 100644 src/init.cpp delete mode 100644 tools/make_cc.R diff --git a/.Rbuildignore b/.Rbuildignore index 6f55499..f2fec2f 100644 --- a/.Rbuildignore +++ b/.Rbuildignore @@ -6,3 +6,5 @@ cran-comments.md ^doc$ ^Meta$ ^CRAN-RELEASE$ +^src$ +^\.github$ diff --git a/.github/.gitignore b/.github/.gitignore new file mode 100644 index 0000000..2d19fc7 --- /dev/null +++ b/.github/.gitignore @@ -0,0 +1 @@ +*.html diff --git a/.github/workflows/R-CMD-check.yaml b/.github/workflows/R-CMD-check.yaml new file mode 100644 index 0000000..85c5068 --- /dev/null +++ b/.github/workflows/R-CMD-check.yaml @@ -0,0 +1,50 @@ +# Workflow derived from https://github.com/r-lib/actions/tree/v2/examples +# Need help debugging build failures? Start at https://github.com/r-lib/actions#where-to-find-help +on: + push: + branches: [main, master] + pull_request: + branches: [main, master] + +name: R-CMD-check + +jobs: + R-CMD-check: + runs-on: ${{ matrix.config.os }} + + name: ${{ matrix.config.os }} (${{ matrix.config.r }}) + + strategy: + fail-fast: false + matrix: + config: + - {os: macos-latest, r: 'release'} + - {os: windows-latest, r: 'release'} + - {os: ubuntu-latest, r: 'devel', http-user-agent: 'release'} + - {os: ubuntu-latest, r: 'release'} + - {os: ubuntu-latest, r: 'oldrel-1'} + + env: + GITHUB_PAT: ${{ secrets.GITHUB_TOKEN }} + R_KEEP_PKG_SOURCE: yes + + steps: + - uses: actions/checkout@v3 + + - uses: r-lib/actions/setup-pandoc@v2 + + - uses: r-lib/actions/setup-r@v2 + with: + r-version: ${{ matrix.config.r }} + http-user-agent: ${{ matrix.config.http-user-agent }} + use-public-rspm: true + + - uses: r-lib/actions/setup-r-dependencies@v2 + with: + extra-packages: any::rcmdcheck + needs: check + + - uses: r-lib/actions/check-r-package@v2 + with: + upload-snapshots: true + error-on: "error" diff --git a/R/stanmodels.R b/R/stanmodels.R deleted file mode 100644 index a92a6ec..0000000 --- a/R/stanmodels.R +++ /dev/null @@ -1,37 +0,0 @@ -# Part of the bayes4psy package for estimating model parameters -# Copyright (C) 2015, 2016, 2017 Trustees of Columbia University -# -# This program is free software; you can redistribute it and/or -# modify it under the terms of the GNU General Public License -# as published by the Free Software Foundation; either version 3 -# of the License, or (at your option) any later version. -# -# This program is distributed in the hope that it will be useful, -# but WITHOUT ANY WARRANTY; without even the implied warranty of -# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the -# GNU General Public License for more details. -# -# You should have received a copy of the GNU General Public License -# along with this program; if not, write to the Free Software -# Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. - -# This file is only intended to be used during the installation process -# nocov start -MODELS_HOME <- "src" -if (!file.exists(MODELS_HOME)) MODELS_HOME <- sub("R$", "src", getwd()) - -stan_files <- dir(file.path(MODELS_HOME, "stan_files"), - pattern="stan$", full.names=TRUE) -stanmodels <- lapply(stan_files, function(f) { - model_cppname <- sub("\\.stan$", "", basename(f)) - stanfit <- rstan::stanc(f, allow_undefined=TRUE, - obfuscate_model_name=FALSE) - stanfit$model_cpp <- list(model_cppname=stanfit$model_name, - model_cppcode=stanfit$cppcode) - return(do.call(methods::new, args=c(stanfit[-(1:3)], Class="stanmodel", - mk_cppmodule=function(x) get(paste0("model_", model_cppname))))) - } -) -names(stanmodels) <- sub("\\.stan$", "", basename(stan_files)) -rm(MODELS_HOME) -# nocov end diff --git a/R/zzz.R b/R/zzz.R deleted file mode 100644 index eed32cf..0000000 --- a/R/zzz.R +++ /dev/null @@ -1,4 +0,0 @@ -.onLoad <- function(libname, pkgname) { - modules <- paste0("stan_fit4", names(stanmodels), "_mod") - for (m in modules) loadModule(m, what=TRUE) -} diff --git a/configure b/configure index 0304fc5..1c04798 100644 --- a/configure +++ b/configure @@ -1,5 +1,4 @@ -#! /bin/sh - # Generated by rstantools. Do not edit by hand. +#! /bin/sh "${R_HOME}/bin/Rscript" -e "rstantools::rstan_config()" diff --git a/configure.win b/configure.win index 5e2dceb..94d77bd 100644 --- a/configure.win +++ b/configure.win @@ -1,5 +1,4 @@ -#! /bin/sh - # Generated by rstantools. Do not edit by hand. +#! /bin/sh "${R_HOME}/bin${R_ARCH_BIN}/Rscript.exe" -e "rstantools::rstan_config()" diff --git a/inst/include/meta_header.hpp b/inst/include/meta_header.hpp deleted file mode 100644 index 3b914da..0000000 --- a/inst/include/meta_header.hpp +++ /dev/null @@ -1 +0,0 @@ -// Insert all #include statements here diff --git a/src/stan_files/color.stan b/inst/stan/color.stan similarity index 100% rename from src/stan_files/color.stan rename to inst/stan/color.stan diff --git a/src/stan_files/chunks/license.stan b/inst/stan/include/license.stan similarity index 100% rename from src/stan_files/chunks/license.stan rename to inst/stan/include/license.stan diff --git a/src/stan_files/linear.stan b/inst/stan/linear.stan similarity index 100% rename from src/stan_files/linear.stan rename to inst/stan/linear.stan diff --git a/src/stan_files/reaction_time.stan b/inst/stan/reaction_time.stan similarity index 100% rename from src/stan_files/reaction_time.stan rename to inst/stan/reaction_time.stan diff --git a/src/stan_files/success_rate.stan b/inst/stan/success_rate.stan similarity index 100% rename from src/stan_files/success_rate.stan rename to inst/stan/success_rate.stan diff --git a/src/stan_files/ttest.stan b/inst/stan/ttest.stan similarity index 100% rename from src/stan_files/ttest.stan rename to inst/stan/ttest.stan diff --git a/src/Makevars b/src/Makevars deleted file mode 100644 index 810342b..0000000 --- a/src/Makevars +++ /dev/null @@ -1,26 +0,0 @@ -STANHEADERS_SRC = $(shell "$(R_HOME)/bin$(R_ARCH_BIN)/Rscript" -e "message()" -e "cat(system.file('include', 'src', package = 'StanHeaders', mustWork = TRUE))" -e "message()" | grep "StanHeaders") - -STANC_FLAGS = $(shell "$(R_HOME)/bin$(R_ARCH_BIN)/Rscript" -e "cat(ifelse(utils::packageVersion('rstan') >= 2.26, '-DUSE_STANC3',''))") -PKG_CPPFLAGS = -I"../inst/include" -I"$(STANHEADERS_SRC)" -DBOOST_DISABLE_ASSERTS -DEIGEN_NO_DEBUG -DBOOST_MATH_OVERFLOW_ERROR_POLICY=errno_on_error $(STANC_FLAGS) -PKG_CXXFLAGS = $(shell "$(R_HOME)/bin$(R_ARCH_BIN)/Rscript" -e "RcppParallel::CxxFlags()") $(shell "$(R_HOME)/bin$(R_ARCH_BIN)/Rscript" -e "StanHeaders:::CxxFlags()") -PKG_LIBS = $(shell "$(R_HOME)/bin$(R_ARCH_BIN)/Rscript" -e "RcppParallel::RcppParallelLibs()") $(shell "$(R_HOME)/bin$(R_ARCH_BIN)/Rscript" -e "StanHeaders:::LdFlags()") - - -CXX_STD = CXX14 - -SOURCES = $(wildcard stan_files/*.stan) -OBJECTS = $(SOURCES:.stan=.o) init.o - -all: $(SHLIB) - @if test -e "/usr/bin/install_name_tool" && test -e "/usr/local/clang4/lib/libc++.1.dylib" && test -e "/usr/lib/libc++.1.dylib"; then /usr/bin/install_name_tool -change /usr/local/clang4/lib/libc++.1.dylib /usr/lib/libc++.1.dylib $(SHLIB); fi - -clean: - rm -rf stan_files/*.o - rm -rf *.so *.o - rm -rf stan_files/*.cc - rm -rf stan_files/*.hpp - -%.cc: %.stan - "$(R_HOME)/bin$(R_ARCH_BIN)/Rscript" -e "source(file.path('..', 'tools', 'make_cc.R')); make_cc(commandArgs(TRUE))" $< - -.phony: all clean diff --git a/src/Makevars.win b/src/Makevars.win deleted file mode 100644 index 895a819..0000000 --- a/src/Makevars.win +++ /dev/null @@ -1,26 +0,0 @@ -STANHEADERS_SRC = $(shell "$(R_HOME)/bin$(R_ARCH_BIN)/Rscript" -e "message()" -e "cat(system.file('include', 'src', package = 'StanHeaders', mustWork = TRUE))" -e "message()" | grep "StanHeaders") - -STANC_FLAGS = $(shell "$(R_HOME)/bin$(R_ARCH_BIN)/Rscript" -e "cat(ifelse(utils::packageVersion('rstan') >= 2.26, '-DUSE_STANC3',''))") -PKG_CPPFLAGS = -I"../inst/include" -I"$(STANHEADERS_SRC)" -DBOOST_DISABLE_ASSERTS -DEIGEN_NO_DEBUG -DRCPP_PARALLEL_USE_TBB=1 $(STANC_FLAGS) -PKG_CXXFLAGS = $(shell "$(R_HOME)/bin$(R_ARCH_BIN)/Rscript" -e "RcppParallel::CxxFlags()") $(shell "$(R_HOME)/bin$(R_ARCH_BIN)/Rscript" -e "StanHeaders:::CxxFlags()") -PKG_LIBS = $(shell "$(R_HOME)/bin$(R_ARCH_BIN)/Rscript" -e "RcppParallel::RcppParallelLibs()") $(shell "$(R_HOME)/bin$(R_ARCH_BIN)/Rscript" -e "StanHeaders:::LdFlags()") - -CXX_STD = CXX14 - -SOURCES = $(wildcard stan_files/*.stan) -OBJECTS = $(SOURCES:.stan=.o) init.o - -all: $(SHLIB) - -clean: - RM -rf stan_files/*.o - RM -rf *.so *.o - RM -rf stan_files/*.cc - RM -rf stan_files/*.hpp - -%.cc: %.stan - "$(R_HOME)/bin$(R_ARCH_BIN)/Rscript" -e "source(file.path('..', 'tools', 'make_cc.R')); make_cc(commandArgs(TRUE))" $< - - - -.phony: clean diff --git a/src/init.cpp b/src/init.cpp deleted file mode 100644 index d3ad850..0000000 --- a/src/init.cpp +++ /dev/null @@ -1,20 +0,0 @@ -// Generated by the rstantools package - - -#include -#include -#include -#include -#include - - -static const R_CallMethodDef CallEntries[] = { - {NULL, NULL, 0} -}; - - -void attribute_visible R_init_bayes4psy(DllInfo *dll) { - // next line is necessary to avoid a NOTE from R CMD check - R_registerRoutines(dll, NULL, CallEntries, NULL, NULL); - R_useDynamicSymbols(dll, TRUE); // necessary for .onLoad() to work -} diff --git a/tools/make_cc.R b/tools/make_cc.R deleted file mode 100644 index 2c0e4c1..0000000 --- a/tools/make_cc.R +++ /dev/null @@ -1,48 +0,0 @@ -# Part of the rstanarm package for estimating model parameters -# Copyright (C) 2015, 2016, 2017 Trustees of Columbia University -# -# This program is free software; you can redistribute it and/or -# modify it under the terms of the GNU General Public License -# as published by the Free Software Foundation; either version 3 -# of the License, or (at your option) any later version. -# -# This program is distributed in the hope that it will be useful, -# but WITHOUT ANY WARRANTY; without even the implied warranty of -# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the -# GNU General Public License for more details. -# -# You should have received a copy of the GNU General Public License -# along with this program; if not, write to the Free Software -# Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. -options(warn = 3L) -options("useFancyQuotes" = FALSE) - -make_cc <- function(file) { - file <- sub("\\.cc$", ".stan", file) - cppcode <- rstan::stanc(file, allow_undefined = TRUE, - obfuscate_model_name = FALSE)$cppcode - cppcode <- sub("(class[[:space:]]+[A-Za-z_][A-Za-z0-9_]*[[:space:]]*: public prob_grad \\{)", - paste("#include \n", "\\1"), cppcode) - - cat(readLines(dir("stan_files", pattern = "license.stan", recursive = TRUE, full.names = TRUE)), - "#ifndef MODELS_HPP", "#define MODELS_HPP", "#define STAN__SERVICES__COMMAND_HPP", - "#include ", - cppcode, "#endif", file = sub("\\.stan$", ".hpp", file), - sep = "\n", append = FALSE) - - f <- sub("\\.stan$", "", basename(file)) - Rcpp::exposeClass(class = paste0("model_", f), - constructors = list(c("SEXP", "SEXP", "SEXP")), fields = character(), - methods = c("call_sampler", - "param_names", "param_names_oi", "param_fnames_oi", - "param_dims", "param_dims_oi", "update_param_oi", "param_oi_tidx", - "grad_log_prob", "log_prob", - "unconstrain_pars", "constrain_pars", "num_pars_unconstrained", - "unconstrained_param_names", "constrained_param_names"), - file = file.path("stan_files", paste0(f, ".cc")), - header = paste0('#include "', f, '.hpp"'), - module = paste0("stan_fit4", f, "_mod"), - CppClass = "rstan::stan_fit ", - Rfile = FALSE) - return(invisible(NULL)) -} From 0d85bbc0affd447065488ced122b8186e721058b Mon Sep 17 00:00:00 2001 From: Andrew Johnson Date: Tue, 14 Mar 2023 12:30:33 +0200 Subject: [PATCH 02/13] Update arg --- .github/workflows/R-CMD-check.yaml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.github/workflows/R-CMD-check.yaml b/.github/workflows/R-CMD-check.yaml index 85c5068..2021bd7 100644 --- a/.github/workflows/R-CMD-check.yaml +++ b/.github/workflows/R-CMD-check.yaml @@ -47,4 +47,4 @@ jobs: - uses: r-lib/actions/check-r-package@v2 with: upload-snapshots: true - error-on: "error" + error-on: "'error'" From 2a7eba783a50be7f82014599d38fd80876dea5c4 Mon Sep 17 00:00:00 2001 From: Andrew Johnson Date: Tue, 14 Mar 2023 12:39:40 +0200 Subject: [PATCH 03/13] Executable --- configure | 0 configure.win | 0 2 files changed, 0 insertions(+), 0 deletions(-) mode change 100644 => 100755 configure mode change 100644 => 100755 configure.win diff --git a/configure b/configure old mode 100644 new mode 100755 diff --git a/configure.win b/configure.win old mode 100644 new mode 100755 From ccea5b58d4580c8484732052d35b3e401ee7d969 Mon Sep 17 00:00:00 2001 From: Andrew Johnson Date: Tue, 14 Mar 2023 12:53:28 +0200 Subject: [PATCH 04/13] Testing? --- .Rbuildignore | 1 - 1 file changed, 1 deletion(-) diff --git a/.Rbuildignore b/.Rbuildignore index f2fec2f..e206297 100644 --- a/.Rbuildignore +++ b/.Rbuildignore @@ -6,5 +6,4 @@ cran-comments.md ^doc$ ^Meta$ ^CRAN-RELEASE$ -^src$ ^\.github$ From d9a532c2f1f57ff5e76968b74690b9e2cae6cb87 Mon Sep 17 00:00:00 2001 From: Andrew Johnson Date: Tue, 14 Mar 2023 12:58:33 +0200 Subject: [PATCH 05/13] Cleanup workflow --- .Rbuildignore | 1 - .github/.gitignore | 1 - .github/workflows/R-CMD-check.yaml | 50 ------------------------------ 3 files changed, 52 deletions(-) delete mode 100644 .github/.gitignore delete mode 100644 .github/workflows/R-CMD-check.yaml diff --git a/.Rbuildignore b/.Rbuildignore index e206297..6f55499 100644 --- a/.Rbuildignore +++ b/.Rbuildignore @@ -6,4 +6,3 @@ cran-comments.md ^doc$ ^Meta$ ^CRAN-RELEASE$ -^\.github$ diff --git a/.github/.gitignore b/.github/.gitignore deleted file mode 100644 index 2d19fc7..0000000 --- a/.github/.gitignore +++ /dev/null @@ -1 +0,0 @@ -*.html diff --git a/.github/workflows/R-CMD-check.yaml b/.github/workflows/R-CMD-check.yaml deleted file mode 100644 index 2021bd7..0000000 --- a/.github/workflows/R-CMD-check.yaml +++ /dev/null @@ -1,50 +0,0 @@ -# Workflow derived from https://github.com/r-lib/actions/tree/v2/examples -# Need help debugging build failures? Start at https://github.com/r-lib/actions#where-to-find-help -on: - push: - branches: [main, master] - pull_request: - branches: [main, master] - -name: R-CMD-check - -jobs: - R-CMD-check: - runs-on: ${{ matrix.config.os }} - - name: ${{ matrix.config.os }} (${{ matrix.config.r }}) - - strategy: - fail-fast: false - matrix: - config: - - {os: macos-latest, r: 'release'} - - {os: windows-latest, r: 'release'} - - {os: ubuntu-latest, r: 'devel', http-user-agent: 'release'} - - {os: ubuntu-latest, r: 'release'} - - {os: ubuntu-latest, r: 'oldrel-1'} - - env: - GITHUB_PAT: ${{ secrets.GITHUB_TOKEN }} - R_KEEP_PKG_SOURCE: yes - - steps: - - uses: actions/checkout@v3 - - - uses: r-lib/actions/setup-pandoc@v2 - - - uses: r-lib/actions/setup-r@v2 - with: - r-version: ${{ matrix.config.r }} - http-user-agent: ${{ matrix.config.http-user-agent }} - use-public-rspm: true - - - uses: r-lib/actions/setup-r-dependencies@v2 - with: - extra-packages: any::rcmdcheck - needs: check - - - uses: r-lib/actions/check-r-package@v2 - with: - upload-snapshots: true - error-on: "'error'" From 1eae54d158ad7739e16c80d89c27e8d3fa47774d Mon Sep 17 00:00:00 2001 From: Jure Demsar Date: Thu, 16 Mar 2023 12:23:23 +0100 Subject: [PATCH 06/13] CRAN prep. --- DESCRIPTION | 2 +- NEWS.md | 4 ++++ cran-comments.md | 4 ++++ 3 files changed, 9 insertions(+), 1 deletion(-) diff --git a/DESCRIPTION b/DESCRIPTION index a82cf1c..583fcca 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,5 +1,5 @@ Package: bayes4psy -Version: 1.2.10 +Version: 1.2.11 Title: User Friendly Bayesian Data Analysis for Psychology Description: Contains several Bayesian models for data analysis of psychological tests. A user friendly interface for these models should enable students and researchers to perform professional level Bayesian data analysis without advanced knowledge in programming and Bayesian statistics. This package is based on the Stan platform (Carpenter et el. 2017 ). Authors@R: diff --git a/NEWS.md b/NEWS.md index cb26a44..83f587f 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,5 +1,9 @@ # bayes4psy change log +## bayes4psy 1.2.11 + +Optimized building of RStan models. + ## bayes4psy 1.2.10 Tidied up the README file. diff --git a/cran-comments.md b/cran-comments.md index 5127ef9..242bdba 100644 --- a/cran-comments.md +++ b/cran-comments.md @@ -1,5 +1,9 @@ # Revisions +## CRAN submission 16. 3. 2023 + +Optimized building of RStan models. + ## CRAN submission 27. 9. 2021 Additional fixes required for newer versions of the RStan package. From 3a33e1492d813daa7469ada3a2c2b7ba61c8bfc4 Mon Sep 17 00:00:00 2001 From: Jure Demsar Date: Thu, 16 Mar 2023 13:33:52 +0100 Subject: [PATCH 07/13] Updated some dependencies. --- DESCRIPTION | 2 +- man/bayes4psy-datasets.Rd | 6 ++++-- 2 files changed, 5 insertions(+), 3 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index 583fcca..1e7072d 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -44,6 +44,6 @@ LinkingTo: SystemRequirements: GNU make VignetteBuilder: knitr NeedsCompilation: yes -RoxygenNote: 7.1.1 +RoxygenNote: 7.2.1 URL: https://github.com/bstatcomp/bayes4psy BugReports: https://github.com/bstatcomp/bayes4psy/issues diff --git a/man/bayes4psy-datasets.Rd b/man/bayes4psy-datasets.Rd index e112620..6055046 100644 --- a/man/bayes4psy-datasets.Rd +++ b/man/bayes4psy-datasets.Rd @@ -15,7 +15,8 @@ Example datasets for use in \pkg{rstanarm} examples and vignettes. The datasets were extracted from the internal MBLab \url{http://www.mblab.si} repository. MBLab is a research lab at the Faculty of Arts, Department of Psychology, University of Ljubljana, Slovenia.} -\format{\describe{ +\format{ +\describe{ \item{\code{adaptation_level_small}}{ Small dataset on subjects picking up weights and determining their weights from 1..10. @@ -154,7 +155,8 @@ Source: Internal MBLab repository. \item \code{naming_incongruent} average response time for naming incongruent stimuli. } } -}} +} +} \description{ Datasets for bayes4psy examples Example datasets for use in \pkg{rstanarm} examples and vignettes. From 08222b82475c3b0e1fa7d74c2a00632257df628a Mon Sep 17 00:00:00 2001 From: Andrew Johnson Date: Wed, 22 Mar 2023 10:38:28 +0200 Subject: [PATCH 08/13] Update gitignore and DESCRIPTION, add testing workflow --- .Rbuildignore | 1 + .github/.gitignore | 1 + .github/workflows/R-CMD-check.yaml | 50 ++++++++++++++++++++++++++++++ .gitignore | 1 - DESCRIPTION | 11 ++++--- NAMESPACE | 3 +- R/bayes4psy-package.R | 3 +- README.md | 3 ++ inst/include/stan_meta_header.hpp | 1 + 9 files changed, 67 insertions(+), 7 deletions(-) create mode 100644 .github/.gitignore create mode 100644 .github/workflows/R-CMD-check.yaml create mode 100644 inst/include/stan_meta_header.hpp diff --git a/.Rbuildignore b/.Rbuildignore index 6f55499..e206297 100644 --- a/.Rbuildignore +++ b/.Rbuildignore @@ -6,3 +6,4 @@ cran-comments.md ^doc$ ^Meta$ ^CRAN-RELEASE$ +^\.github$ diff --git a/.github/.gitignore b/.github/.gitignore new file mode 100644 index 0000000..2d19fc7 --- /dev/null +++ b/.github/.gitignore @@ -0,0 +1 @@ +*.html diff --git a/.github/workflows/R-CMD-check.yaml b/.github/workflows/R-CMD-check.yaml new file mode 100644 index 0000000..600f7ac --- /dev/null +++ b/.github/workflows/R-CMD-check.yaml @@ -0,0 +1,50 @@ +# Workflow derived from https://github.com/r-lib/actions/tree/v2/examples +# Need help debugging build failures? Start at https://github.com/r-lib/actions#where-to-find-help +on: + push: + branches: [main, master, develop] + pull_request: + branches: [main, master, develop] + +name: R-CMD-check + +jobs: + R-CMD-check: + runs-on: ${{ matrix.config.os }} + + name: ${{ matrix.config.os }} (${{ matrix.config.r }}) + + strategy: + fail-fast: false + matrix: + config: + - {os: macos-latest, r: 'release'} + - {os: windows-latest, r: 'release'} + - {os: ubuntu-latest, r: 'devel', http-user-agent: 'release'} + - {os: ubuntu-latest, r: 'release'} + - {os: ubuntu-latest, r: 'oldrel-1'} + + env: + GITHUB_PAT: ${{ secrets.GITHUB_TOKEN }} + R_KEEP_PKG_SOURCE: yes + + steps: + - uses: actions/checkout@v3 + + - uses: r-lib/actions/setup-pandoc@v2 + + - uses: r-lib/actions/setup-r@v2 + with: + r-version: ${{ matrix.config.r }} + http-user-agent: ${{ matrix.config.http-user-agent }} + use-public-rspm: true + + - uses: r-lib/actions/setup-r-dependencies@v2 + with: + extra-packages: any::rcmdcheck + needs: check + + - uses: r-lib/actions/check-r-package@v2 + with: + upload-snapshots: true + error-on: '"error"' diff --git a/.gitignore b/.gitignore index 06bacdc..cf4f1c4 100644 --- a/.gitignore +++ b/.gitignore @@ -41,6 +41,5 @@ vignettes/*.pdf *.o *.so *.cc -*.hpp doc Meta diff --git a/DESCRIPTION b/DESCRIPTION index 1e7072d..d4506f7 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -30,8 +30,9 @@ Imports: mcmcse (>= 1.4.1), scales (>= 1.1.1), stats (>= 4.0.0), - utils (>= 4.0.0) -Suggests: + utils (>= 4.0.0), + RcppParallel +Suggests: testthat (>= 3.0.0), rmarkdown (>= 2.5.0), knitr (>= 1.30.0) @@ -40,10 +41,12 @@ LinkingTo: rstan (>= 2.21.0), BH (>= 1.72.0.3), Rcpp (>= 1.0.5), - RcppEigen (>= 0.3.3.7.0) + RcppEigen (>= 0.3.3.7.0), + RcppParallel SystemRequirements: GNU make VignetteBuilder: knitr NeedsCompilation: yes -RoxygenNote: 7.2.1 +UseLTO: true +RoxygenNote: 7.2.3 URL: https://github.com/bstatcomp/bayes4psy BugReports: https://github.com/bstatcomp/bayes4psy/issues diff --git a/NAMESPACE b/NAMESPACE index 1d848f0..cf2e3c6 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -32,6 +32,7 @@ import(dplyr) import(ggplot2) import(methods) import(rstan) -import(rstantools) +importFrom(RcppParallel,RcppParallelLibs) importFrom(rstan,sampling) +importFrom(rstantools,rstan_config) useDynLib(bayes4psy, .registration = TRUE) diff --git a/R/bayes4psy-package.R b/R/bayes4psy-package.R index 504881e..3074171 100644 --- a/R/bayes4psy-package.R +++ b/R/bayes4psy-package.R @@ -8,7 +8,8 @@ #' @useDynLib bayes4psy, .registration = TRUE #' @import methods #' @import Rcpp -#' @import rstantools +#' @importFrom rstantools rstan_config +#' @importFrom RcppParallel RcppParallelLibs #' @importFrom rstan sampling #' #' @references diff --git a/README.md b/README.md index 1ac550f..d91309f 100644 --- a/README.md +++ b/README.md @@ -1,6 +1,9 @@ # bayes4psy—User Friendly Bayesian Data Analysis for Psychology + +[![R-CMD-check](https://github.com/bstatcomp/bayes4psy/actions/workflows/R-CMD-check.yaml/badge.svg)](https://github.com/bstatcomp/bayes4psy/actions/workflows/R-CMD-check.yaml) [![CRAN\_Status\_Badge](http://www.r-pkg.org/badges/version/bayes4psy?color=blue)](https://CRAN.R-project.org/package=bayes4psy)[![Downloads](http://cranlogs.r-pkg.org/badges/bayes4psy?color=blue)](https://CRAN.R-project.org/package=bayes4psy) + ## Authors diff --git a/inst/include/stan_meta_header.hpp b/inst/include/stan_meta_header.hpp new file mode 100644 index 0000000..3b914da --- /dev/null +++ b/inst/include/stan_meta_header.hpp @@ -0,0 +1 @@ +// Insert all #include statements here From 9a178a8cf7170b102787deb436b9ed5b373f482c Mon Sep 17 00:00:00 2001 From: Jure Demsar Date: Wed, 22 Mar 2023 15:36:04 +0100 Subject: [PATCH 09/13] Removed unnecessary notes in compilation. --- DESCRIPTION | 2 - R/reaction_time_class.R | 12 +- R/success_rate_class.R | 269 ++++++++++++---------- R/ttest_class.R | 395 +++++++++++++++++--------------- man/success_rate_class-class.Rd | 56 ++--- man/ttest_class-class.Rd | 76 +++--- 6 files changed, 424 insertions(+), 386 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index d4506f7..0f004b5 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -28,9 +28,7 @@ Imports: rstan (>= 2.21.2), rstantools (>= 2.1.1), mcmcse (>= 1.4.1), - scales (>= 1.1.1), stats (>= 4.0.0), - utils (>= 4.0.0), RcppParallel Suggests: testthat (>= 3.0.0), diff --git a/R/reaction_time_class.R b/R/reaction_time_class.R index da20352..381ca55 100644 --- a/R/reaction_time_class.R +++ b/R/reaction_time_class.R @@ -485,7 +485,7 @@ setMethod(f="compare_means", signature(object="reaction_time_class"), definition # provided a list of fits i <- 2 for (fit in arguments$fits) { - if (class(fit) != "reaction_time_class") { + if (!("reaction_time_class" %in% class(fit))) { stop("One of the fits in the fits list is not a valid reaction_time_class object.") } @@ -608,7 +608,7 @@ setMethod(f="plot_means_difference", signature(object="reaction_time_class"), de } else if (!is.null(arguments$fits)) { i <- 2 for (fit in arguments$fits) { - if (class(fit) != "reaction_time_class") { + if (!("reaction_time_class" %in% class(fit))) { stop("One of the fits in the fits list is not a valid reaction_time_class object.") } @@ -737,7 +737,7 @@ setMethod(f="plot_means", signature(object="reaction_time_class"), definition=fu } else if (!is.null(arguments$fits)) { i <- 2 for (fit in arguments$fits) { - if (class(fit) != "reaction_time_class") { + if (!("reaction_time_class" %in% class(fit))) { stop("One of the fits in the fits list is not a valid reaction_time_class object.") } @@ -836,7 +836,7 @@ setMethod(f="compare_distributions", signature(object="reaction_time_class"), de # provided a list of fits i <- 2 for (fit in arguments$fits) { - if (class(fit) != "reaction_time_class") { + if (!("reaction_time_class" %in% class(fit))) { stop("One of the fits in the fits list is not a valid reaction_time_class object.") } @@ -927,7 +927,7 @@ setMethod(f="plot_distributions", signature(object="reaction_time_class"), defin } else if (!is.null(arguments$fits)) { i <- 2 for (fit in arguments$fits) { - if (class(fit) != "reaction_time_class") { + if (!("reaction_time_class" %in% class(fit))) { stop("One of the fits in the fits list is not a valid reaction_time_class object.") } @@ -1040,7 +1040,7 @@ setMethod(f="plot_distributions_difference", signature(object="reaction_time_cla } else if (!is.null(arguments$fits)) { i <- 2 for (fit in arguments$fits) { - if (class(fit) != "reaction_time_class") { + if (!("reaction_time_class" %in% class(fit))) { stop("One of the fits in the fits list is not a valid reaction_time_class object.") } diff --git a/R/success_rate_class.R b/R/success_rate_class.R index 8292bc8..8b27cb8 100644 --- a/R/success_rate_class.R +++ b/R/success_rate_class.R @@ -60,29 +60,31 @@ #' #' @examples #' \donttest{ -#' #priors -#' p_prior <- b_prior(family="beta", pars=c(1, 1)) -#' tau_prior <- b_prior(family="uniform", pars=c(0, 500)) +#' # priors +#' p_prior <- b_prior(family = "beta", pars = c(1, 1)) +#' tau_prior <- b_prior(family = "uniform", pars = c(0, 500)) #' #' # attach priors to relevant parameters -#' priors <- list(c("p", p_prior), -#' c("tau", tau_prior)) +#' priors <- list( +#' c("p", p_prior), +#' c("tau", tau_prior) +#' ) #' #' # subjects #' s <- rep(1:5, 20) #' #' # generate data and fit -#' data1 <- rbinom(100, size=1, prob=0.6) -#' fit1 <- b_success_rate(r=data1, s=s, priors=priors, chains=1) +#' data1 <- rbinom(100, size = 1, prob = 0.6) +#' fit1 <- b_success_rate(r = data1, s = s, priors = priors, chains = 1) #' -#' data2 <- rbinom(100, size=1, prob=0.1) -#' fit2 <- b_success_rate(r=data2, s=s, priors=priors, chains=1) +#' data2 <- rbinom(100, size = 1, prob = 0.1) +#' fit2 <- b_success_rate(r = data2, s = s, priors = priors, chains = 1) #' -#' data3 <- rbinom(100, size=1, prob=0.5) -#' fit3 <- b_success_rate(r=data3, s=s, priors=priors, chains=1) +#' data3 <- rbinom(100, size = 1, prob = 0.5) +#' fit3 <- b_success_rate(r = data3, s = s, priors = priors, chains = 1) #' -#' data4 <- rbinom(100, size=1, prob=0.9) -#' fit4 <- b_success_rate(r=data4, s=s, priors=priors, chains=1) +#' data4 <- rbinom(100, size = 1, prob = 0.9) +#' fit4 <- b_success_rate(r = data4, s = s, priors = priors, chains = 1) #' #' # fit list #' fit_list <- list(fit2, fit3, fit4) @@ -100,8 +102,8 @@ #' #' # plot the fitted distribution against the data, #' # plot on the top (group) level -#' plot(fit1, subjects=FALSE) -#' plot_fit(fit1, subjects=FALSE) +#' plot(fit1, subjects = FALSE) +#' plot_fit(fit1, subjects = FALSE) #' #' # traceplot of the fitted parameters #' plot_trace(fit1) @@ -113,49 +115,49 @@ #' subject_parameters <- get_subject_parameters(fit1) #' #' # compare means between two fits, use a rope interval -#' compare_means(fit1, fit2=fit2, rope=0.05) +#' compare_means(fit1, fit2 = fit2, rope = 0.05) #' #' # compare means between multiple fits -#' compare_means(fit1, fits=fit_list) +#' compare_means(fit1, fits = fit_list) #' #' # visualize difference in means between two fits, #' # specify number of histogram bins and rope interval -#' plot_means_difference(fit1, fit2=fit2, bins=40, rope=0.05) +#' plot_means_difference(fit1, fit2 = fit2, bins = 40, rope = 0.05) #' #' # visualize difference in means between multiple fits -#' plot_means_difference(fit1, fits=fit_list) +#' plot_means_difference(fit1, fits = fit_list) #' #' # visualize means of a single fit #' plot_means(fit1) #' #' # visualize means of two fits -#' plot_means(fit1, fit2=fit2) +#' plot_means(fit1, fit2 = fit2) #' #' # visualize means of multiple fits -#' plot_means(fit1, fits=fit_list) +#' plot_means(fit1, fits = fit_list) #' #' # draw samples from distributions underlying two fits and compare them, #' # use a rope interval -#' compare_distributions(fit1, fit2=fit2, rope=0.05) +#' compare_distributions(fit1, fit2 = fit2, rope = 0.05) #' #' # draw samples from distributions underlying multiple fits and compare them -#' compare_distributions(fit1, fits=fit_list) +#' compare_distributions(fit1, fits = fit_list) #' #' # visualize the distribution underlying a fit #' plot_distributions(fit1) #' #' # visualize distributions underlying two fits -#' plot_distributions(fit1, fit2=fit2) +#' plot_distributions(fit1, fit2 = fit2) #' #' # visualize distributions underlying multiple fits -#' plot_distributions(fit1, fits=fit_list) +#' plot_distributions(fit1, fits = fit_list) #' #' # visualize difference between distributions underlying two fits, #' # use a rope interval -#' plot_distributions_difference(fit1, fit2=fit2, rope=0.05) +#' plot_distributions_difference(fit1, fit2 = fit2, rope = 0.05) #' #' # visualize difference between distributions underlying multiple fits -#' plot_distributions_difference(fit1, fits=fit_list) +#' plot_distributions_difference(fit1, fits = fit_list) #' } #' success_rate_class <- setClass( @@ -180,7 +182,7 @@ success_rate_class <- setClass( #' # along with an example of how to use this function #' ?success_rate_class #' -setMethod(f="summary", signature(object="success_rate_class"), definition=function(object) { +setMethod(f = "summary", signature(object = "success_rate_class"), definition = function(object) { # get means p <- mean(object@extract$p0) tau <- mean(object@extract$tau) @@ -190,10 +192,14 @@ setMethod(f="summary", signature(object="success_rate_class"), definition=functi tau_hdi <- mcmc_hdi(object@extract$tau) # print - cat(sprintf("Success rate:\t%.2f +/- %.5f\t95%% HDI: [%.2f, %.2f]\n", - p, mcmcse::mcse(object@extract$p0)$se, p_hdi[1], p_hdi[2])) - cat(sprintf("Tau:\t\t%.2f +/- %.5f\t95%% HDI: [%.2f, %.2f]\n", - tau, mcmcse::mcse(object@extract$tau)$se, tau_hdi[1], tau_hdi[2])) + cat(sprintf( + "Success rate:\t%.2f +/- %.5f\t95%% HDI: [%.2f, %.2f]\n", + p, mcmcse::mcse(object@extract$p0)$se, p_hdi[1], p_hdi[2] + )) + cat(sprintf( + "Tau:\t\t%.2f +/- %.5f\t95%% HDI: [%.2f, %.2f]\n", + tau, mcmcse::mcse(object@extract$tau)$se, tau_hdi[1], tau_hdi[2] + )) }) @@ -208,7 +214,7 @@ setMethod(f="summary", signature(object="success_rate_class"), definition=functi #' # along with an example of how to use this function #' ?success_rate_class #' -setMethod(f="show", signature(object="success_rate_class"), definition=function(object) { +setMethod(f = "show", signature(object = "success_rate_class"), definition = function(object) { # print show(object@fit) }) @@ -227,8 +233,8 @@ setMethod(f="show", signature(object="success_rate_class"), definition=function( #' # along with an example of how to use this function #' ?success_rate_class #' -setMethod(f="plot", signature(x="success_rate_class", y="missing"), definition=function(x, ...) { - return(plot_fit(object=x, ...)) +setMethod(f = "plot", signature(x = "success_rate_class", y = "missing"), definition = function(x, ...) { + return(plot_fit(object = x, ...)) }) @@ -246,7 +252,7 @@ setMethod(f="plot", signature(x="success_rate_class", y="missing"), definition=f #' # along with an example of how to use this function #' ?success_rate_class #' -setMethod(f="plot_fit", signature(object="success_rate_class"), definition=function(object, ...) { +setMethod(f = "plot_fit", signature(object = "success_rate_class"), definition = function(object, ...) { # init local varibales for CRAN check variable <- value <- NULL @@ -258,34 +264,36 @@ setMethod(f="plot_fit", signature(object="success_rate_class"), definition=funct subjects <- arguments$subjects } - df_data <- data.frame(value=object@data$r, variable=object@data$s) + df_data <- data.frame(value = object@data$r, variable = object@data$s) if (!subjects) { mean_p <- mean(df_data$value) - df_fit <- data.frame(value=object@extract$p0) + df_fit <- data.frame(value = object@extract$p0) graph <- ggplot() + - geom_vline(xintercept=mean_p, color="#ff4e3f") + - geom_density(data=df_fit, aes(x=value), fill="#3182bd", alpha=0.4, color=NA) + + geom_vline(xintercept = mean_p, color = "#ff4e3f") + + geom_density(data = df_fit, aes(x = value), fill = "#3182bd", alpha = 0.4, color = NA) + xlim(0, 1) + xlab("success rate") } else { - df_data <- df_data %>% group_by(variable) %>% summarize(value=mean(value)) + df_data <- df_data %>% + group_by(variable) %>% + summarize(value = mean(value)) df_fit <- object@extract$p colnames(df_fit) <- seq(1:ncol(df_fit)) - df_fit <- reshape::melt(as.data.frame(df_fit), id=NULL) + df_fit <- reshape::melt(as.data.frame(df_fit), id = NULL) # ncol n_col <- ceiling(sqrt(nrow(df_data))) # density per subject graph <- ggplot() + - geom_vline(data=df_data, aes(xintercept=value), color="#ff4e3f") + - geom_density(data=df_fit, aes(x=value), fill="#3182bd", alpha=0.4, color=NA) + + geom_vline(data = df_data, aes(xintercept = value), color = "#ff4e3f") + + geom_density(data = df_fit, aes(x = value), fill = "#3182bd", alpha = 0.4, color = NA) + xlim(0, 1) + - facet_wrap(~ variable, ncol=n_col) + + facet_wrap(~variable, ncol = n_col) + xlab("success rate") } @@ -306,8 +314,8 @@ setMethod(f="plot_fit", signature(object="success_rate_class"), definition=funct #' # along with an example of how to use this function #' ?success_rate_class #' -setMethod(f="plot_trace", signature(object="success_rate_class"), definition=function(object) { - rstan::traceplot(object@fit, pars=c("p"), inc_warmup=TRUE) +setMethod(f = "plot_trace", signature(object = "success_rate_class"), definition = function(object) { + rstan::traceplot(object@fit, pars = c("p"), inc_warmup = TRUE) }) @@ -324,9 +332,11 @@ setMethod(f="plot_trace", signature(object="success_rate_class"), definition=fun #' # along with an example of how to use this function #' ?success_rate_class #' -setMethod(f="get_parameters", signature(object="success_rate_class"), definition=function(object) { - df <- data.frame(p=object@extract$p0, - tau=object@extract$tau) +setMethod(f = "get_parameters", signature(object = "success_rate_class"), definition = function(object) { + df <- data.frame( + p = object@extract$p0, + tau = object@extract$tau + ) return(df) }) @@ -345,14 +355,16 @@ setMethod(f="get_parameters", signature(object="success_rate_class"), definition #' # along with an example of how to use this function #' ?success_rate_class #' -setMethod(f="get_subject_parameters", signature(object="success_rate_class"), definition=function(object) { - df <- data.frame(p=numeric(), subject=numeric()) +setMethod(f = "get_subject_parameters", signature(object = "success_rate_class"), definition = function(object) { + df <- data.frame(p = numeric(), subject = numeric()) n <- length(unique(object@data$s)) for (i in 1:n) { - df_subject <- data.frame(p = object@extract$p[,i], - subject = i) + df_subject <- data.frame( + p = object@extract$p[, i], + subject = i + ) df <- rbind(df, df_subject) } @@ -375,7 +387,7 @@ setMethod(f="get_subject_parameters", signature(object="success_rate_class"), de #' # along with an example of how to use this function #' ?success_rate_class #' -setMethod(f="compare_means", signature(object="success_rate_class"), definition=function(object, ...) { +setMethod(f = "compare_means", signature(object = "success_rate_class"), definition = function(object, ...) { arguments <- list(...) wrong_arguments <- "The provided arguments for the compare_means function are invalid, compare_means(success_rate_class, fit2=success_rate_class) or compare_means(success_rate_class, fits=list) is required! You can also provide the rope parameter, e.g. compare_means(success_rate_class, fit2=success_rate_class, rope=numeric)." @@ -410,7 +422,7 @@ setMethod(f="compare_means", signature(object="success_rate_class"), definition= # provided a list of fits i <- 2 for (fit in arguments$fits) { - if (class(fit) != "success_rate_class") { + if (!("success_rate_class" %in% class(fit))) { stop("One of the fits in the fits list is not a valid success_rate_class object.") } y[[i]] <- fit@extract$p0 @@ -422,12 +434,12 @@ setMethod(f="compare_means", signature(object="success_rate_class"), definition= n <- length(y) comparison_matrix <- matrix(nrow = n, ncol = n) - for (i in 1:(n-1)) { - for (j in (i+1):n) { + for (i in 1:(n - 1)) { + for (j in (i + 1):n) { cat(sprintf("\n---------- Group %d vs Group %d ----------\n", i, j)) - result <- difference(y1=y[[i]], y2=y[[j]], rope=rope, group1=i, group2=j) - comparison_matrix[j,i] <- result[1] - comparison_matrix[i,j] <- result[2] + result <- difference(y1 = y[[i]], y2 = y[[j]], rope = rope, group1 = i, group2 = j) + comparison_matrix[j, i] <- result[1] + comparison_matrix[i, j] <- result[2] cat("\n") } } @@ -436,10 +448,10 @@ setMethod(f="compare_means", signature(object="success_rate_class"), definition= if (n > 2) { cat("-----------------------------------------") cat("\nProbabilities that a certain group is\nsmallest/largest or equal to all others:\n\n") - smallest_largest <- is_smallest_or_largest(data=y, rope=rope) + smallest_largest <- is_smallest_or_largest(data = y, rope = rope) print(smallest_largest) cat("\n\n") - return(list(comparison_matrix=comparison_matrix, smallest_largest=smallest_largest)) + return(list(comparison_matrix = comparison_matrix, smallest_largest = smallest_largest)) } else { return(comparison_matrix) } @@ -460,7 +472,7 @@ setMethod(f="compare_means", signature(object="success_rate_class"), definition= #' # along with an example of how to use this function #' ?success_rate_class #' -setMethod(f="plot_means_difference", signature(object="success_rate_class"), definition=function(object, ...) { +setMethod(f = "plot_means_difference", signature(object = "success_rate_class"), definition = function(object, ...) { # init local varibales for CRAN check value <- NULL @@ -495,7 +507,7 @@ setMethod(f="plot_means_difference", signature(object="success_rate_class"), def } else if (!is.null(arguments$fits)) { i <- 2 for (fit in arguments$fits) { - if (class(fit) != "success_rate_class") { + if (!("success_rate_class" %in% class(fit))) { stop("One of the fits in the fits list is not a valid success_rate_class object.") } y[[i]] <- fit@extract$p0 @@ -515,7 +527,7 @@ setMethod(f="plot_means_difference", signature(object="success_rate_class"), def # if no list is provided if (is.null(arguments$fits)) { # call plot difference shared function - graph <- plot_difference(y1=y[[1]], y2=y[[2]], rope=rope, bins=bins) + graph <- plot_difference(y1 = y[[1]], y2 = y[[2]], rope = rope, bins = bins) return(graph) } else { graphs <- list() @@ -524,24 +536,24 @@ setMethod(f="plot_means_difference", signature(object="success_rate_class"), def for (j in i:n) { # if both are equal plot means, else plot difference if (i == j) { - df <- data.frame(value=y[[i]]) - index <- (i-1)*n + i + df <- data.frame(value = y[[i]]) + index <- (i - 1) * n + i graphs[[index]] <- ggplot() + - geom_density(data=df, aes(x=value), fill="#3182bd", color=NA, alpha=0.4) + + geom_density(data = df, aes(x = value), fill = "#3182bd", color = NA, alpha = 0.4) + xlab("probability") + xlim(0, 1) } else { - index1 <- (i-1)*n + j - graphs[[index1]] <- plot_difference(y1=y[[i]], y2=y[[j]], rope=rope, bins=bins, nrow=n) + index1 <- (i - 1) * n + j + graphs[[index1]] <- plot_difference(y1 = y[[i]], y2 = y[[j]], rope = rope, bins = bins, nrow = n) - index2 <- (j-1)*n + i - graphs[[index2]] <- plot_difference(y1=y[[j]], y2=y[[i]], rope=rope, bins=bins, nrow=n) + index2 <- (j - 1) * n + i + graphs[[index2]] <- plot_difference(y1 = y[[j]], y2 = y[[i]], rope = rope, bins = bins, nrow = n) } } } # cowplot - graph <- suppressWarnings(cowplot::plot_grid(plotlist=graphs, nrow=n, ncol=n, scale=0.9)) + graph <- suppressWarnings(cowplot::plot_grid(plotlist = graphs, nrow = n, ncol = n, scale = 0.9)) return(graph) } }) @@ -561,12 +573,12 @@ setMethod(f="plot_means_difference", signature(object="success_rate_class"), def #' # along with an example of how to use this function #' ?success_rate_class #' -setMethod(f="plot_means", signature(object="success_rate_class"), definition=function(object, ...) { +setMethod(f = "plot_means", signature(object = "success_rate_class"), definition = function(object, ...) { # init local varibales for CRAN check group <- value <- NULL # first group data - df <- data.frame(value=object@extract$p0, group="1") + df <- data.frame(value = object@extract$p0, group = "1") # second group data arguments <- list(...) @@ -579,14 +591,14 @@ setMethod(f="plot_means", signature(object="success_rate_class"), definition=fun fit2 <- arguments[[1]] } - df <- rbind(df, data.frame(value=fit2@extract$p0, group="2")) + df <- rbind(df, data.frame(value = fit2@extract$p0, group = "2")) } else if (!is.null(arguments$fits)) { i <- 2 for (fit in arguments$fits) { - if (class(fit) != "success_rate_class") { + if (!("success_rate_class" %in% class(fit))) { stop("One of the fits in the fits list is not a valid success_rate_class object.") } - df <- rbind(df, data.frame(value=fit@extract$p0, group=as.factor(i))) + df <- rbind(df, data.frame(value = fit@extract$p0, group = as.factor(i))) i <- i + 1 } } @@ -594,21 +606,21 @@ setMethod(f="plot_means", signature(object="success_rate_class"), definition=fun # plot graph <- ggplot() + - geom_density(data=df, aes(x=value, fill=group), color=NA, alpha=0.4) + + geom_density(data = df, aes(x = value, fill = group), color = NA, alpha = 0.4) + xlab("probability") + xlim(0, 1) n_groups <- max(as.numeric(df$group)) if (n_groups == 2) { graph <- graph + - scale_fill_manual(values=c("#3182bd", "#ff4e3f")) + scale_fill_manual(values = c("#3182bd", "#ff4e3f")) } else if (n_groups > 2) { graph <- graph + scale_fill_hue() } else { graph <- graph + - scale_fill_manual(values=c("#3182bd")) + - theme(legend.position="none") + scale_fill_manual(values = c("#3182bd")) + + theme(legend.position = "none") } return(suppressWarnings(graph)) @@ -629,7 +641,7 @@ setMethod(f="plot_means", signature(object="success_rate_class"), definition=fun #' # along with an example of how to use this function #' ?success_rate_class #' -setMethod(f="compare_distributions", signature(object="success_rate_class"), definition=function(object, ...) { +setMethod(f = "compare_distributions", signature(object = "success_rate_class"), definition = function(object, ...) { arguments <- list(...) wrong_arguments <- "The provided arguments for the compare_distributions function are invalid, compare_distributions(success_rate_class, fit2=success_rate_class) or compare_distributions(success_rate_class, fits=list) is required!. You can also provide the rope parameter, e.g. compare_distributions(success_rate_class, fit2=success_rate_class, rope=numeric)." @@ -650,7 +662,7 @@ setMethod(f="compare_distributions", signature(object="success_rate_class"), def n <- 100000 p0 <- mean(object@extract$p0) tau <- mean(object@extract$tau) - y[[1]] <- stats::rbeta(n, p0*tau, (1 - p0)*tau) + y[[1]] <- stats::rbeta(n, p0 * tau, (1 - p0) * tau) # second group data if (!is.null(arguments$fit2) || class(arguments[[1]])[1] == "success_rate_class") { @@ -662,17 +674,17 @@ setMethod(f="compare_distributions", signature(object="success_rate_class"), def } p0 <- mean(fit2@extract$p0) tau <- mean(fit2@extract$tau) - y[[2]] <- stats::rbeta(n, p0*tau, (1 - p0)*tau) + y[[2]] <- stats::rbeta(n, p0 * tau, (1 - p0) * tau) } else if (!is.null(arguments$fits)) { # provided a list of fits i <- 2 for (fit in arguments$fits) { - if (class(fit) != "success_rate_class") { + if (!("success_rate_class" %in% class(fit))) { stop("One of the fits in the fits list is not a valid success_rate_class object.") } p0 <- mean(fit@extract$p0) tau <- mean(fit@extract$tau) - y[[i]] <- stats::rbeta(n, p0*tau, (1 - p0)*tau) + y[[i]] <- stats::rbeta(n, p0 * tau, (1 - p0) * tau) i <- i + 1 } } else { @@ -681,12 +693,12 @@ setMethod(f="compare_distributions", signature(object="success_rate_class"), def n <- length(y) comparison_matrix <- matrix(nrow = n, ncol = n) - for (i in 1:(n-1)) { - for (j in (i+1):n) { + for (i in 1:(n - 1)) { + for (j in (i + 1):n) { cat(sprintf("\n---------- Group %d vs Group %d ----------\n", i, j)) - result <- difference(y1=y[[i]], y2=y[[j]], rope=rope, group1=i, group2=j) - comparison_matrix[j,i] <- result[1] - comparison_matrix[i,j] <- result[2] + result <- difference(y1 = y[[i]], y2 = y[[j]], rope = rope, group1 = i, group2 = j) + comparison_matrix[j, i] <- result[1] + comparison_matrix[i, j] <- result[2] cat("\n") } } @@ -695,10 +707,10 @@ setMethod(f="compare_distributions", signature(object="success_rate_class"), def if (n > 2) { cat("-----------------------------------------") cat("\nProbabilities that a certain group is\nsmallest/largest or equal to all others:\n\n") - smallest_largest <- is_smallest_or_largest(data=y, rope=rope) + smallest_largest <- is_smallest_or_largest(data = y, rope = rope) print(smallest_largest) cat("\n\n") - return(list(comparison_matrix=comparison_matrix, smallest_largest=smallest_largest)) + return(list(comparison_matrix = comparison_matrix, smallest_largest = smallest_largest)) } else { return(comparison_matrix) } @@ -719,7 +731,7 @@ setMethod(f="compare_distributions", signature(object="success_rate_class"), def #' # along with an example of how to use this function #' ?success_rate_class #' -setMethod(f="plot_distributions", signature(object="success_rate_class"), definition=function(object, ...) { +setMethod(f = "plot_distributions", signature(object = "success_rate_class"), definition = function(object, ...) { # init local varibales for CRAN check group <- x <- y <- NULL @@ -728,8 +740,8 @@ setMethod(f="plot_distributions", signature(object="success_rate_class"), defini betas <- vector() p0 <- mean(object@extract$p0) tau <- mean(object@extract$tau) - alphas[[1]] <- p0*tau - betas[[1]] <- (1 - p0)*tau + alphas[[1]] <- p0 * tau + betas[[1]] <- (1 - p0) * tau # second group data arguments <- list(...) @@ -743,18 +755,18 @@ setMethod(f="plot_distributions", signature(object="success_rate_class"), defini } p0 <- mean(fit2@extract$p0) tau <- mean(fit2@extract$tau) - alphas[[2]] <- p0*tau - betas[[2]] <- (1 - p0)*tau + alphas[[2]] <- p0 * tau + betas[[2]] <- (1 - p0) * tau } else if (!is.null(arguments$fits)) { i <- 2 for (fit in arguments$fits) { - if (class(fit) != "success_rate_class") { + if (!("success_rate_class" %in% class(fit))) { stop("One of the fits in the fits list is not a valid success_rate_class object.") } p0 <- mean(fit@extract$p0) tau <- mean(fit@extract$tau) - alphas[[i]] <- p0*tau - betas[[i]] <- (1 - p0)*tau + alphas[[i]] <- p0 * tau + betas[[i]] <- (1 - p0) * tau i <- i + 1 } } @@ -762,34 +774,37 @@ setMethod(f="plot_distributions", signature(object="success_rate_class"), defini # calculate data points step <- 1 / 1000 - df <- data.frame(x=numeric(), y=numeric(), group=factor()) + df <- data.frame(x = numeric(), y = numeric(), group = factor()) n_groups <- length(alphas) for (i in 1:n_groups) { - df_group <- data.frame(x = seq(0, 1, step), - y = stats::dbeta(seq(0, 1, step), - shape1 = alphas[i], - shape2 = betas[i]), - group=as.factor(i)) + df_group <- data.frame( + x = seq(0, 1, step), + y = stats::dbeta(seq(0, 1, step), + shape1 = alphas[i], + shape2 = betas[i] + ), + group = as.factor(i) + ) df <- rbind(df, df_group) } # plot graph <- ggplot() + - geom_area(data=df, aes(x=x, y=y, fill=group), alpha=0.4, position="identity") + + geom_area(data = df, aes(x = x, y = y, fill = group), alpha = 0.4, position = "identity") + xlab("probability") + ylab("density") if (n_groups == 2) { graph <- graph + - scale_fill_manual(values=c("#3182bd", "#ff4e3f")) + scale_fill_manual(values = c("#3182bd", "#ff4e3f")) } else if (n_groups > 2) { graph <- graph + scale_fill_hue() } else { graph <- graph + - scale_fill_manual(values=c("#3182bd")) + - theme(legend.position="none") + scale_fill_manual(values = c("#3182bd")) + + theme(legend.position = "none") } return(suppressWarnings(graph)) @@ -810,7 +825,7 @@ setMethod(f="plot_distributions", signature(object="success_rate_class"), defini #' # along with an example of how to use this function #' ?success_rate_class #' -setMethod(f="plot_distributions_difference", signature(object="success_rate_class"), definition=function(object, ...) { +setMethod(f = "plot_distributions_difference", signature(object = "success_rate_class"), definition = function(object, ...) { # init local varibales for CRAN check value <- NULL @@ -834,7 +849,7 @@ setMethod(f="plot_distributions_difference", signature(object="success_rate_clas n <- 100000 p0 <- mean(object@extract$p0) tau <- mean(object@extract$tau) - y[[1]] <- stats::rbeta(n, p0*tau, (1 - p0)*tau) + y[[1]] <- stats::rbeta(n, p0 * tau, (1 - p0) * tau) # second group data if (!is.null(arguments$fit2) || class(arguments[[1]])[1] == "success_rate_class") { @@ -846,16 +861,16 @@ setMethod(f="plot_distributions_difference", signature(object="success_rate_clas } p0 <- mean(fit2@extract$p0) tau <- mean(fit2@extract$tau) - y[[2]] <- stats::rbeta(n, p0*tau, (1 - p0)*tau) + y[[2]] <- stats::rbeta(n, p0 * tau, (1 - p0) * tau) } else if (!is.null(arguments$fits)) { i <- 2 for (fit in arguments$fits) { - if (class(fit) != "success_rate_class") { + if (!("success_rate_class" %in% class(fit))) { stop("One of the fits in the fits list is not a valid success_rate_class object.") } p0 <- mean(fit@extract$p0) tau <- mean(fit@extract$tau) - y[[i]] <- stats::rbeta(n, p0*tau, (1 - p0)*tau) + y[[i]] <- stats::rbeta(n, p0 * tau, (1 - p0) * tau) i <- i + 1 } } else { @@ -871,7 +886,7 @@ setMethod(f="plot_distributions_difference", signature(object="success_rate_clas # if no list is provided if (is.null(arguments$fits)) { # call plot difference shared function - graph <- plot_difference(y1=y[[1]], y2=y[[2]], rope=rope, bins=bins) + graph <- plot_difference(y1 = y[[1]], y2 = y[[2]], rope = rope, bins = bins) return(graph) } else { graphs <- list() @@ -880,24 +895,24 @@ setMethod(f="plot_distributions_difference", signature(object="success_rate_clas for (j in i:n) { # if both are equal plot samples, else plot difference if (i == j) { - df <- data.frame(value=y[[i]]) - index <- (i-1)*n + i + df <- data.frame(value = y[[i]]) + index <- (i - 1) * n + i graphs[[index]] <- ggplot() + - geom_density(data=df, aes(x=value), fill="#3182bd", color=NA, alpha=0.4) + + geom_density(data = df, aes(x = value), fill = "#3182bd", color = NA, alpha = 0.4) + xlab("probability") + xlim(0, 1) } else { - index1 <- (i-1)*n + j - graphs[[index1]] <- plot_difference(y1=y[[i]], y2=y[[j]], rope=rope, bins=bins, nrow=n) + index1 <- (i - 1) * n + j + graphs[[index1]] <- plot_difference(y1 = y[[i]], y2 = y[[j]], rope = rope, bins = bins, nrow = n) - index2 <- (j-1)*n + i - graphs[[index2]] <- plot_difference(y1=y[[j]], y2=y[[i]], rope=rope, bins=bins, nrow=n) + index2 <- (j - 1) * n + i + graphs[[index2]] <- plot_difference(y1 = y[[j]], y2 = y[[i]], rope = rope, bins = bins, nrow = n) } } } # cowplot - graph <- suppressWarnings(cowplot::plot_grid(plotlist=graphs, nrow=n, ncol=n, scale=0.9)) + graph <- suppressWarnings(cowplot::plot_grid(plotlist = graphs, nrow = n, ncol = n, scale = 0.9)) return(graph) } }) diff --git a/R/ttest_class.R b/R/ttest_class.R index 7c45ac9..c59c2c5 100644 --- a/R/ttest_class.R +++ b/R/ttest_class.R @@ -73,27 +73,29 @@ #' @examples #' \donttest{ #' # priors -#' mu_prior <- b_prior(family="normal", pars=c(0, 1000)) -#' sigma_prior <- b_prior(family="uniform", pars=c(0, 500)) -#' nu_prior <- b_prior(family="normal", pars=c(2000, 1000)) +#' mu_prior <- b_prior(family = "normal", pars = c(0, 1000)) +#' sigma_prior <- b_prior(family = "uniform", pars = c(0, 500)) +#' nu_prior <- b_prior(family = "normal", pars = c(2000, 1000)) #' #' # attach priors to relevant parameters -#' priors <- list(c("mu", mu_prior), -#' c("sigma", sigma_prior), -#' c("nu", nu_prior)) +#' priors <- list( +#' c("mu", mu_prior), +#' c("sigma", sigma_prior), +#' c("nu", nu_prior) +#' ) #' #' # generate data and fit -#' data1 <- rnorm(20, mean=150, sd=20) -#' fit1 <- b_ttest(data=data1, priors=priors, chains=1) +#' data1 <- rnorm(20, mean = 150, sd = 20) +#' fit1 <- b_ttest(data = data1, priors = priors, chains = 1) #' -#' data2 <- rnorm(20, mean=200, sd=20) -#' fit2 <- b_ttest(data=data2, priors=priors, chains=1) +#' data2 <- rnorm(20, mean = 200, sd = 20) +#' fit2 <- b_ttest(data = data2, priors = priors, chains = 1) #' -#' data3 <- rnorm(20, mean=150, sd=40) -#' fit3 <- b_ttest(data=data3, priors=priors, chains=1) +#' data3 <- rnorm(20, mean = 150, sd = 40) +#' fit3 <- b_ttest(data = data3, priors = priors, chains = 1) #' -#' data4 <- rnorm(20, mean=50, sd=10) -#' fit4 <- b_ttest(data=data4, priors=priors, chains=1) +#' data4 <- rnorm(20, mean = 50, sd = 10) +#' fit4 <- b_ttest(data = data4, priors = priors, chains = 1) #' #' # fit list #' fit_list <- list(fit2, fit3, fit4) @@ -116,92 +118,92 @@ #' parameters <- get_parameters(fit1) #' #' # compare means between two fits -#' compare_means(fit1, fit2=fit2) +#' compare_means(fit1, fit2 = fit2) #' #' # compare means between two fits, use a rope interval -#' compare_means(fit1, fit2=fit2, rope=2) +#' compare_means(fit1, fit2 = fit2, rope = 2) #' #' # compare means between a fit and a constant value -#' compare_means(fit1, mu=150) +#' compare_means(fit1, mu = 150) #' #' # compare means between a fit and a distribution, #' # sigma is used for calculating Cohen's d -#' compare_means(fit1, mu=150, sigma=20) +#' compare_means(fit1, mu = 150, sigma = 20) #' #' # compare means between multiple fits -#' compare_means(fit1, fits=fit_list) +#' compare_means(fit1, fits = fit_list) #' #' # visualize difference in means between two fits, #' # specify number of histogram bins -#' plot_means_difference(fit1, fit2=fit2, bins=20) +#' plot_means_difference(fit1, fit2 = fit2, bins = 20) #' #' # visualize difference in means between a fit and a constant value -#' plot_means_difference(fit1, mu=150) +#' plot_means_difference(fit1, mu = 150) #' #' # visualize difference in means between multiple fits, use a rope interval -#' plot_means_difference(fit1, fits=fit_list, rope=2) +#' plot_means_difference(fit1, fits = fit_list, rope = 2) #' #' # visualize means of a single fit #' plot_means(fit1) #' #' # visualize means of two fits -#' plot_means(fit1, fit2=fit2) +#' plot_means(fit1, fit2 = fit2) #' #' # visualize means of a fit and a constant value -#' plot_means(fit1, mu=150) +#' plot_means(fit1, mu = 150) #' #' # visualize means of multiple fits -#' plot_means(fit1, fits=fit_list) +#' plot_means(fit1, fits = fit_list) #' #' # draw samples from distributions underlying two fits and compare them -#' compare_distributions(fit1, fit2=fit2) +#' compare_distributions(fit1, fit2 = fit2) #' #' # draw samples from a distribution underlying the fit #' # and compare them with a constant, use a rope interval -#' compare_distributions(fit1, mu=150, rope=2) +#' compare_distributions(fit1, mu = 150, rope = 2) #' #' # draw samples from a distribution underlying the fit and #' # compare them with a user defined distribution -#' compare_distributions(fit1, mu=150, sigma=20) +#' compare_distributions(fit1, mu = 150, sigma = 20) #' #' # draw samples from distributions underlying multiple fits and compare them -#' compare_distributions(fit1, fits=fit_list) +#' compare_distributions(fit1, fits = fit_list) #' #' # visualize the distribution underlying a fit #' plot_distributions(fit1) #' #' # visualize distributions underlying two fits -#' plot_distributions(fit1, fit2=fit2) +#' plot_distributions(fit1, fit2 = fit2) #' #' # visualize the distribution underlying a fit and a constant value -#' plot_distributions(fit1, mu=150) +#' plot_distributions(fit1, mu = 150) #' #' # visualize the distribution underlying a fit and a user defined distribution -#' plot_distributions(fit1, mu=150, sigma=20) +#' plot_distributions(fit1, mu = 150, sigma = 20) #' #' # visualize distributions underlying multiple fits -#' plot_distributions(fit1, fits=fit_list) +#' plot_distributions(fit1, fits = fit_list) #' #' # visualize difference between distributions underlying two fits, #' # use a rope interval -#' plot_distributions_difference(fit1, fit2=fit2, rope=2) +#' plot_distributions_difference(fit1, fit2 = fit2, rope = 2) #' #' # visualize difference between a distribution underlying the fit #' # and a constant value -#' plot_distributions_difference(fit1, mu=150) +#' plot_distributions_difference(fit1, mu = 150) #' #' # visualize difference between a distribution underlying the fits #' # and a user defined distribution -#' plot_distributions_difference(fit1, mu=150, sigma=20) +#' plot_distributions_difference(fit1, mu = 150, sigma = 20) #' #' # visualize difference between distributions underlying multiple fits -#' plot_distributions_difference(fit1, fits=fit_list) +#' plot_distributions_difference(fit1, fits = fit_list) #' } #' ttest_class <- setClass( "ttest_class", slots = c( - extract = "list", + extract = "list", fit = "stanfit", data = "numeric" ), @@ -219,7 +221,7 @@ ttest_class <- setClass( #' # along with an example of how to use this function #' ?ttest_class #' -setMethod(f="summary", signature(object="ttest_class"), definition=function(object) { +setMethod(f = "summary", signature(object = "ttest_class"), definition = function(object) { # get means mu <- mean(object@extract$mu) sigma <- mean(object@extract$sigma) @@ -231,12 +233,18 @@ setMethod(f="summary", signature(object="ttest_class"), definition=function(obje nu_hdi <- mcmc_hdi(object@extract$nu) # print) - cat(sprintf("mu:\t\t%.2f +/- %.5f\t95%% HDI: [%.2f, %.2f]\n", - mu, mcmcse::mcse(object@extract$mu)$se, mu_hdi[1], mu_hdi[2])) - cat(sprintf("sigma:\t\t%.2f +/- %.5f\t95%% HDI: [%.2f, %.2f]\n", - sigma, mcmcse::mcse(object@extract$sigma)$se, sigma_hdi[1], sigma_hdi[2])) - cat(sprintf("nu:\t\t%.2f +/- %.5f\t95%% HDI: [%.2f, %.2f]\n", nu, - mcmcse::mcse(object@extract$nu)$se, nu_hdi[1], nu_hdi[2])) + cat(sprintf( + "mu:\t\t%.2f +/- %.5f\t95%% HDI: [%.2f, %.2f]\n", + mu, mcmcse::mcse(object@extract$mu)$se, mu_hdi[1], mu_hdi[2] + )) + cat(sprintf( + "sigma:\t\t%.2f +/- %.5f\t95%% HDI: [%.2f, %.2f]\n", + sigma, mcmcse::mcse(object@extract$sigma)$se, sigma_hdi[1], sigma_hdi[2] + )) + cat(sprintf( + "nu:\t\t%.2f +/- %.5f\t95%% HDI: [%.2f, %.2f]\n", nu, + mcmcse::mcse(object@extract$nu)$se, nu_hdi[1], nu_hdi[2] + )) }) @@ -251,7 +259,7 @@ setMethod(f="summary", signature(object="ttest_class"), definition=function(obje #' # along with an example of how to use this function #' ?ttest_class #' -setMethod(f="show", signature(object="ttest_class"), definition=function(object) { +setMethod(f = "show", signature(object = "ttest_class"), definition = function(object) { # print show(object@fit) }) @@ -268,8 +276,8 @@ setMethod(f="show", signature(object="ttest_class"), definition=function(object) #' # along with an example of how to use this function #' ?ttest_class #' -setMethod(f="plot", signature(x="ttest_class", y="missing"), definition=function(x) { - return(plot_fit(object=x)) +setMethod(f = "plot", signature(x = "ttest_class", y = "missing"), definition = function(x) { + return(plot_fit(object = x)) }) @@ -286,29 +294,29 @@ setMethod(f="plot", signature(x="ttest_class", y="missing"), definition=function #' # along with an example of how to use this function #' ?ttest_class #' -setMethod(f="plot_fit", signature(object="ttest_class"), definition=function(object) { +setMethod(f = "plot_fit", signature(object = "ttest_class"), definition = function(object) { # init local varibales for CRAN check value <- NULL - df_data <- data.frame(value=object@data) + df_data <- data.frame(value = object@data) nu <- mean(object@extract$nu) mu <- mean(object@extract$mu) sigma <- mean(object@extract$sigma) # get x range - x_min <- min(mu - 4*sigma, df_data$value) - x_max <- max(mu + 4*sigma, df_data$value) + x_min <- min(mu - 4 * sigma, df_data$value) + x_max <- max(mu + 4 * sigma, df_data$value) diff <- x_max - x_min - x_min <- x_min - 0.1*diff - x_max <- x_max + 0.1*diff + x_min <- x_min - 0.1 * diff + x_max <- x_max + 0.1 * diff - df_x <- data.frame(x=c(x_min, x_max)) + df_x <- data.frame(x = c(x_min, x_max)) - graph <- ggplot(data=df_x) + - geom_density(data=df_data, aes(x=value), fill="#3182bd", alpha=0.4, color=NA) + - stat_function(fun=metRology::dt.scaled, n=10000, args=list(df=nu, mean=mu, sd=sigma), colour="#3182bd", size=1) + + graph <- ggplot(data = df_x) + + geom_density(data = df_data, aes(x = value), fill = "#3182bd", alpha = 0.4, color = NA) + + stat_function(fun = metRology::dt.scaled, n = 10000, args = list(df = nu, mean = mu, sd = sigma), colour = "#3182bd", size = 1) + xlab("value") + xlim(x_min, x_max) @@ -329,8 +337,8 @@ setMethod(f="plot_fit", signature(object="ttest_class"), definition=function(obj #' # along with an example of how to use this function #' ?ttest_class #' -setMethod(f="plot_trace", signature(object="ttest_class"), definition=function(object) { - traceplot(object@fit, pars=c("mu", "sigma", "nu"), inc_warmup=TRUE) +setMethod(f = "plot_trace", signature(object = "ttest_class"), definition = function(object) { + traceplot(object@fit, pars = c("mu", "sigma", "nu"), inc_warmup = TRUE) }) @@ -347,10 +355,12 @@ setMethod(f="plot_trace", signature(object="ttest_class"), definition=function(o #' # along with an example of how to use this function #' ?ttest_class #' -setMethod(f="get_parameters", signature(object="ttest_class"), definition=function(object) { - df <- data.frame(mu=object@extract$mu, - sigma=object@extract$sigma, - nu=object@extract$nu) +setMethod(f = "get_parameters", signature(object = "ttest_class"), definition = function(object) { + df <- data.frame( + mu = object@extract$mu, + sigma = object@extract$sigma, + nu = object@extract$nu + ) return(df) }) @@ -370,7 +380,7 @@ setMethod(f="get_parameters", signature(object="ttest_class"), definition=functi #' # along with an example of how to use this function #' ?ttest_class #' -setMethod(f="compare_means", signature(object="ttest_class"), definition=function(object, ...) { +setMethod(f = "compare_means", signature(object = "ttest_class"), definition = function(object, ...) { arguments <- list(...) wrong_arguments <- "The provided arguments for the compare_means function are invalid, compare_means(ttest_class, fit2=ttest_class), compare_means(ttest_class, mu=numeric), compare_means(ttest_class, mu=numeric, sigma=numeric) or compare_means(ttest_class, fits=list) is required! You provide the rope parameter, e.g. compare_means(ttest_class, fit2=ttest_class, rope=numeric) or execute the comparison through the sigma or nu parameter, e.g. compare_means(ttest_class, fit2=ttest_class, par=\"sigma\")." @@ -437,7 +447,7 @@ setMethod(f="compare_means", signature(object="ttest_class"), definition=functio sigma2 <- mean(fit2@extract$sigma) } else if (!is.null(arguments$mu) && par == "mu") { # provided mu - y[[2]] <- arguments$mu; + y[[2]] <- arguments$mu # provided also sigma? if (!is.null(arguments$sigma)) { @@ -445,13 +455,13 @@ setMethod(f="compare_means", signature(object="ttest_class"), definition=functio } } else if (!is.null(arguments$sigma) && par == "sigma") { # provided sigma - y[[2]] <- arguments$sigma; + y[[2]] <- arguments$sigma sigma2 <- arguments$sigma } else if (!is.null(arguments$fits)) { # provided a list of fits i <- 2 for (fit in arguments$fits) { - if (class(fit) != "ttest_class") { + if (!("ttest_class" %in% class(fit))) { stop("One of the fits in the fits list is not a valid ttest_class object.") } @@ -471,12 +481,12 @@ setMethod(f="compare_means", signature(object="ttest_class"), definition=functio n <- length(y) comparison_matrix <- matrix(nrow = n, ncol = n) - for (i in 1:(n-1)) { - for (j in (i+1):n) { + for (i in 1:(n - 1)) { + for (j in (i + 1):n) { cat(sprintf("\n---------- Group %d vs Group %d ----------\n", i, j)) - result <- difference(y1=y[[i]], y2=y[[j]], rope=rope, group1=i, group2=j) - comparison_matrix[j,i] <- result[1] - comparison_matrix[i,j] <- result[2] + result <- difference(y1 = y[[i]], y2 = y[[j]], rope = rope, group1 = i, group2 = j) + comparison_matrix[j, i] <- result[1] + comparison_matrix[i, j] <- result[2] cat("\n") } } @@ -484,7 +494,7 @@ setMethod(f="compare_means", signature(object="ttest_class"), definition=functio # add Cohen's d if 2 groups are provided if (!is.null(sigma2)) { diff <- mean(y[[1]]) - mean(y[[2]]) - cohens_d <- diff / sqrt((n*sigma1^2 + n*sigma2^2) / (n + n - 2)); + cohens_d <- diff / sqrt((n * sigma1^2 + n * sigma2^2) / (n + n - 2)) cat(sprintf("\nCohen's d: %.2f\n", cohens_d)) } @@ -492,10 +502,10 @@ setMethod(f="compare_means", signature(object="ttest_class"), definition=functio if (n > 2) { cat("-----------------------------------------") cat("\nProbabilities that a certain group is\nsmallest/largest or equal to all others:\n\n") - smallest_largest <- is_smallest_or_largest(data=y, rope=rope) + smallest_largest <- is_smallest_or_largest(data = y, rope = rope) print(smallest_largest) cat("\n\n") - return(list(comparison_matrix=comparison_matrix, smallest_largest=smallest_largest)) + return(list(comparison_matrix = comparison_matrix, smallest_largest = smallest_largest)) } else { return(comparison_matrix) } @@ -516,7 +526,7 @@ setMethod(f="compare_means", signature(object="ttest_class"), definition=functio #' # along with an example of how to use this function #' ?ttest_class #' -setMethod(f="plot_means_difference", signature(object="ttest_class"), definition=function(object, ...) { +setMethod(f = "plot_means_difference", signature(object = "ttest_class"), definition = function(object, ...) { # init local varibales for CRAN check value <- NULL @@ -582,13 +592,13 @@ setMethod(f="plot_means_difference", signature(object="ttest_class"), definition y[[2]] <- fit2@extract$nu } } else if (!is.null(arguments$mu) && par == "mu") { - y[[2]] <- arguments$mu; + y[[2]] <- arguments$mu } else if (!is.null(arguments$sigma) && par == "sigma") { - y[[2]] <- arguments$sigma; + y[[2]] <- arguments$sigma } else if (!is.null(arguments$fits)) { i <- 2 for (fit in arguments$fits) { - if (class(fit) != "ttest_class") { + if (!("ttest_class" %in% class(fit))) { stop("One of the fits in the fits list is not a valid ttest_class object.") } @@ -619,12 +629,12 @@ setMethod(f="plot_means_difference", signature(object="ttest_class"), definition # if no list is provided if (is.null(arguments$fits)) { # call plot difference shared function - graph <- plot_difference(y1=y[[1]], y2=y[[2]], rope=rope, bins=bins) + graph <- plot_difference(y1 = y[[1]], y2 = y[[2]], rope = rope, bins = bins) return(graph) } else { diff <- x_max - x_min - x_min <- x_min - 0.1*diff - x_max <- x_max + 0.1*diff + x_min <- x_min - 0.1 * diff + x_max <- x_max + 0.1 * diff graphs <- list() n <- length(y) @@ -632,24 +642,24 @@ setMethod(f="plot_means_difference", signature(object="ttest_class"), definition for (j in i:n) { # if both are equal plot means, else plot difference if (i == j) { - df <- data.frame(value=y[[i]]) - index <- (i-1)*n + i + df <- data.frame(value = y[[i]]) + index <- (i - 1) * n + i graphs[[index]] <- ggplot() + - geom_density(data=df, aes(x=value), fill="#3182bd", color=NA, alpha=0.4) + + geom_density(data = df, aes(x = value), fill = "#3182bd", color = NA, alpha = 0.4) + xlab("value") + xlim(x_min, x_max) } else { - index1 <- (i-1)*n + j - graphs[[index1]] <- plot_difference(y1=y[[i]], y2=y[[j]], rope=rope, bins=bins, nrow=n) + index1 <- (i - 1) * n + j + graphs[[index1]] <- plot_difference(y1 = y[[i]], y2 = y[[j]], rope = rope, bins = bins, nrow = n) - index2 <- (j-1)*n + i - graphs[[index2]] <- plot_difference(y1=y[[j]], y2=y[[i]], rope=rope, bins=bins, nrow=n) + index2 <- (j - 1) * n + i + graphs[[index2]] <- plot_difference(y1 = y[[j]], y2 = y[[i]], rope = rope, bins = bins, nrow = n) } } } # cowplot - graph <- suppressWarnings(cowplot::plot_grid(plotlist=graphs, nrow=n, ncol=n, scale=0.9)) + graph <- suppressWarnings(cowplot::plot_grid(plotlist = graphs, nrow = n, ncol = n, scale = 0.9)) return(graph) } }) @@ -669,7 +679,7 @@ setMethod(f="plot_means_difference", signature(object="ttest_class"), definition #' # along with an example of how to use this function #' ?ttest_class #' -setMethod(f="plot_means", signature(object="ttest_class"), definition=function(object, ...) { +setMethod(f = "plot_means", signature(object = "ttest_class"), definition = function(object, ...) { # init local varibales for CRAN check group <- value <- NULL @@ -691,11 +701,11 @@ setMethod(f="plot_means", signature(object="ttest_class"), definition=function(o # first group data df <- NULL if (par == "mu") { - df <- data.frame(value=object@extract$mu, group="1") + df <- data.frame(value = object@extract$mu, group = "1") } else if (par == "sigma") { - df <- data.frame(value=object@extract$sigma, group="1") + df <- data.frame(value = object@extract$sigma, group = "1") } else if (par == "nu") { - df <- data.frame(value=object@extract$nu, group="1") + df <- data.frame(value = object@extract$nu, group = "1") } # second group data @@ -711,32 +721,31 @@ setMethod(f="plot_means", signature(object="ttest_class"), definition=function(o } if (par == "mu") { - df <- rbind(df, data.frame(value=fit2@extract$mu, group="2")) + df <- rbind(df, data.frame(value = fit2@extract$mu, group = "2")) } else if (par == "sigma") { - df <- rbind(df, data.frame(value=fit2@extract$sigma, group="2")) + df <- rbind(df, data.frame(value = fit2@extract$sigma, group = "2")) } else if (par == "nu") { - df <- rbind(df, data.frame(value=fit2@extract$nu, group="2")) + df <- rbind(df, data.frame(value = fit2@extract$nu, group = "2")) } - } else if (!is.null(arguments$mu) && par == "mu") { # provided mu and sigma - par2 <- arguments$mu; + par2 <- arguments$mu } else if (!is.null(arguments$sigma) && par == "sigma") { # provided mu and sigma - par2 <- arguments$sigma; + par2 <- arguments$sigma } else if (!is.null(arguments$fits)) { i <- 2 for (fit in arguments$fits) { - if (class(fit) != "ttest_class") { + if (!("ttest_class" %in% class(fit))) { stop("One of the fits in the fits list is not a valid ttest_class object.") } if (par == "mu") { - df <- rbind(df, data.frame(value=fit@extract$mu, group=as.factor(i))) + df <- rbind(df, data.frame(value = fit@extract$mu, group = as.factor(i))) } else if (par == "sigma") { - df <- rbind(df, data.frame(value=fit@extract$sigma, group=as.factor(i))) + df <- rbind(df, data.frame(value = fit@extract$sigma, group = as.factor(i))) } else if (par == "nu") { - df <- rbind(df, data.frame(value=fit@extract$nu, group=as.factor(i))) + df <- rbind(df, data.frame(value = fit@extract$nu, group = as.factor(i))) } i <- i + 1 @@ -746,13 +755,13 @@ setMethod(f="plot_means", signature(object="ttest_class"), definition=function(o # plot graph <- ggplot() + - geom_density(data=df, aes(x=value, fill=group), color=NA, alpha=0.4) + + geom_density(data = df, aes(x = value, fill = group), color = NA, alpha = 0.4) + xlab("value") n_groups <- max(as.numeric(df$group)) if (n_groups == 2) { graph <- graph + - scale_fill_manual(values=c("#3182bd", "#ff4e3f")) + scale_fill_manual(values = c("#3182bd", "#ff4e3f")) } else if (n_groups > 2) { graph <- graph + scale_fill_hue() @@ -760,14 +769,14 @@ setMethod(f="plot_means", signature(object="ttest_class"), definition=function(o y_max <- ggplot_build(graph)$layout$panel_scales_y[[1]]$range$range graph <- graph + - geom_segment(aes(x=par2, xend=par2, y=0, yend=y_max[2]*1.05), size=1.5, color="#ff4e3f", alpha=0.4) + - geom_text(aes(label=sprintf("%.2f", par2), x=par2, y=y_max[2]*1.08), size=4) + - scale_fill_manual(values=c("#3182bd")) + - theme(legend.position="none") + geom_segment(aes(x = par2, xend = par2, y = 0, yend = y_max[2] * 1.05), size = 1.5, color = "#ff4e3f", alpha = 0.4) + + geom_text(aes(label = sprintf("%.2f", par2), x = par2, y = y_max[2] * 1.08), size = 4) + + scale_fill_manual(values = c("#3182bd")) + + theme(legend.position = "none") } else { graph <- graph + - scale_fill_manual(values=c("#3182bd")) + - theme(legend.position="none") + scale_fill_manual(values = c("#3182bd")) + + theme(legend.position = "none") } # limits @@ -782,8 +791,8 @@ setMethod(f="plot_means", signature(object="ttest_class"), definition=function(o } diff <- x_max - x_min - x_min <- x_min - 0.1*diff - x_max <- x_max + 0.1*diff + x_min <- x_min - 0.1 * diff + x_max <- x_max + 0.1 * diff graph <- graph + xlim(x_min, x_max) @@ -805,7 +814,7 @@ setMethod(f="plot_means", signature(object="ttest_class"), definition=function(o #' # along with an example of how to use this function #' ?ttest_class #' -setMethod(f="compare_distributions", signature(object="ttest_class"), definition=function(object, ...) { +setMethod(f = "compare_distributions", signature(object = "ttest_class"), definition = function(object, ...) { arguments <- list(...) wrong_arguments <- "The provided arguments for the compare_distributions function are invalid, compare_distributions(ttest_class, fit2=ttest_class), compare_distributions(ttest_class, fits=list), compare_distributions(ttest_class, mu=numeric) or compare_distributions(ttest_class, mu=numeric, sigma=numeric) is required! You can also provide the rope parameter, e.g. compare_distributions(ttest_class, fit2=ttest_class, rope=numeric)." @@ -826,9 +835,10 @@ setMethod(f="compare_distributions", signature(object="ttest_class"), definition y <- NULL sigma1 <- mean(object@extract$sigma) y[[1]] <- metRology::rt.scaled(n, - df=mean(object@extract$nu), - mean=mean(object@extract$mu), - sd=sigma1) + df = mean(object@extract$nu), + mean = mean(object@extract$mu), + sd = sigma1 + ) # second group data sigma2 <- 0 @@ -842,13 +852,14 @@ setMethod(f="compare_distributions", signature(object="ttest_class"), definition sigma2 <- mean(fit2@extract$sigma) y[[2]] <- metRology::rt.scaled(n, - df=mean(fit2@extract$nu), - mean=mean(fit2@extract$mu), - sd=sigma2) + df = mean(fit2@extract$nu), + mean = mean(fit2@extract$mu), + sd = sigma2 + ) } else if (!is.null(arguments$mu)) { # provided mu and sigma if (!is.null(arguments$sigma)) { - sigma2 <- arguments$sigma; + sigma2 <- arguments$sigma y[[2]] <- stats::rnorm(n, arguments$mu, sigma2) } else { y[[2]] <- stats::rnorm(n, arguments$mu, 0) @@ -856,14 +867,15 @@ setMethod(f="compare_distributions", signature(object="ttest_class"), definition } else if (!is.null(arguments$fits)) { i <- 2 for (fit in arguments$fits) { - if (class(fit) != "ttest_class") { + if (!("ttest_class" %in% class(fit))) { stop("One of the fits in the fits list is not a valid ttest_class object.") } y[[i]] <- metRology::rt.scaled(n, - df=mean(fit@extract$nu), - mean=mean(fit@extract$mu), - sd=mean(fit@extract$sigma)) + df = mean(fit@extract$nu), + mean = mean(fit@extract$mu), + sd = mean(fit@extract$sigma) + ) i <- i + 1 } @@ -873,12 +885,12 @@ setMethod(f="compare_distributions", signature(object="ttest_class"), definition n <- length(y) comparison_matrix <- matrix(nrow = n, ncol = n) - for (i in 1:(n-1)) { - for (j in (i+1):n) { + for (i in 1:(n - 1)) { + for (j in (i + 1):n) { cat(sprintf("\n---------- Group %d vs Group %d ----------\n", i, j)) - result <- difference(y1=y[[i]], y2=y[[j]], rope=rope, group1=i, group2=j) - comparison_matrix[j,i] <- result[1] - comparison_matrix[i,j] <- result[2] + result <- difference(y1 = y[[i]], y2 = y[[j]], rope = rope, group1 = i, group2 = j) + comparison_matrix[j, i] <- result[1] + comparison_matrix[i, j] <- result[2] cat("\n") } } @@ -886,7 +898,7 @@ setMethod(f="compare_distributions", signature(object="ttest_class"), definition # add Cohen's d if 2 groups are provided if (!is.null(sigma2)) { diff <- mean(y[[1]]) - mean(y[[2]]) - cohens_d <- diff / sqrt((n*sigma1^2 + n*sigma2^2) / (n + n - 2)); + cohens_d <- diff / sqrt((n * sigma1^2 + n * sigma2^2) / (n + n - 2)) cat(sprintf("\nCohen's d: %.2f\n", cohens_d)) } @@ -894,10 +906,10 @@ setMethod(f="compare_distributions", signature(object="ttest_class"), definition if (n > 2) { cat("-----------------------------------------") cat("\nProbabilities that a certain group is\nsmallest/largest or equal to all others:\n\n") - smallest_largest <- is_smallest_or_largest(data=y, rope=rope) + smallest_largest <- is_smallest_or_largest(data = y, rope = rope) print(smallest_largest) cat("\n\n") - return(list(comparison_matrix=comparison_matrix, smallest_largest=smallest_largest)) + return(list(comparison_matrix = comparison_matrix, smallest_largest = smallest_largest)) } else { return(comparison_matrix) } @@ -918,7 +930,7 @@ setMethod(f="compare_distributions", signature(object="ttest_class"), definition #' # along with an example of how to use this function #' ?ttest_class #' -setMethod(f="plot_distributions", signature(object="ttest_class"), definition=function(object, ...) { +setMethod(f = "plot_distributions", signature(object = "ttest_class"), definition = function(object, ...) { # init local varibales for CRAN check x <- y <- group <- NULL @@ -947,15 +959,15 @@ setMethod(f="plot_distributions", signature(object="ttest_class"), definition=fu sigmas[2] <- mean(fit2@extract$sigma) } else if (!is.null(arguments$mu)) { # provided mu and/or sigma - mu2 <- arguments$mu; + mu2 <- arguments$mu if (!is.null(arguments$sigma)) { - sigma2 <- arguments$sigma; + sigma2 <- arguments$sigma } } else if (!is.null(arguments$fits)) { i <- 2 for (fit in arguments$fits) { - if (class(fit) != "ttest_class") { + if (!("ttest_class" %in% class(fit))) { stop("One of the fits in the fits list is not a valid ttest_class object.") } @@ -968,43 +980,49 @@ setMethod(f="plot_distributions", signature(object="ttest_class"), definition=fu } # get boundaries - x_min <- min(mus - 4*sigmas) - x_max <- max(mus + 4*sigmas) + x_min <- min(mus - 4 * sigmas) + x_max <- max(mus + 4 * sigmas) if (!is.null(mu2) & !is.null(sigma2)) { - x_min <- min(x_min, mu2 - 4*sigma2) - x_max <- max(x_max, mu2 + 4*sigma2) + x_min <- min(x_min, mu2 - 4 * sigma2) + x_max <- max(x_max, mu2 + 4 * sigma2) } else if (!is.null(mu2)) { if (mu2 < x_min) { x_min <- mu2 - x_min <- x_min - (0.1*(x_max-x_min)) + x_min <- x_min - (0.1 * (x_max - x_min)) } else if (mu2 > x_max) { x_max <- mu2 - x_max <- x_max + (0.1*(x_max-x_min)) + x_max <- x_max + (0.1 * (x_max - x_min)) } } # calculate data points step <- (x_max - x_min) / 1000 - df <- data.frame(x=numeric(), y=numeric(), group=factor()) + df <- data.frame(x = numeric(), y = numeric(), group = factor()) n_groups <- length(mus) for (i in 1:n_groups) { - df_group <- data.frame(x = seq(x_min, x_max, step), - y = metRology::dt.scaled(seq(x_min, x_max, step), - df = nus[i], - mean = mus[i], - sd = sigmas[i]), - group=as.factor(i)) + df_group <- data.frame( + x = seq(x_min, x_max, step), + y = metRology::dt.scaled(seq(x_min, x_max, step), + df = nus[i], + mean = mus[i], + sd = sigmas[i] + ), + group = as.factor(i) + ) df <- rbind(df, df_group) } if (!is.null(mu2) & !is.null(sigma2)) { - df_group <- data.frame(x = seq(x_min, x_max, step), - y = stats::dnorm(seq(x_min, x_max, step), - mean = mu2, - sd = sigma2), - group="2") + df_group <- data.frame( + x = seq(x_min, x_max, step), + y = stats::dnorm(seq(x_min, x_max, step), + mean = mu2, + sd = sigma2 + ), + group = "2" + ) df <- rbind(df, df_group) @@ -1013,13 +1031,13 @@ setMethod(f="plot_distributions", signature(object="ttest_class"), definition=fu # plot graph <- ggplot() + - geom_area(data=df, aes(x=x, y=y, fill=group), alpha=0.4, position="identity") + + geom_area(data = df, aes(x = x, y = y, fill = group), alpha = 0.4, position = "identity") + xlab("value") + ylab("density") if (n_groups == 2) { graph <- graph + - scale_fill_manual(values=c("#3182bd", "#ff4e3f")) + scale_fill_manual(values = c("#3182bd", "#ff4e3f")) } else if (n_groups > 2) { graph <- graph + scale_fill_hue() @@ -1027,14 +1045,14 @@ setMethod(f="plot_distributions", signature(object="ttest_class"), definition=fu y_max <- ggplot_build(graph)$layout$panel_scales_y[[1]]$range$range graph <- graph + - geom_segment(aes(x=mu2, xend=mu2, y=0, yend=y_max[2]*1.05), size=1.5, color="#ff4e3f", alpha=0.4) + - geom_text(aes(label=sprintf("%.2f", mu2), x=mu2, y=y_max[2]*1.08), size=4) + - scale_fill_manual(values=c("#3182bd")) + - theme(legend.position="none") + geom_segment(aes(x = mu2, xend = mu2, y = 0, yend = y_max[2] * 1.05), size = 1.5, color = "#ff4e3f", alpha = 0.4) + + geom_text(aes(label = sprintf("%.2f", mu2), x = mu2, y = y_max[2] * 1.08), size = 4) + + scale_fill_manual(values = c("#3182bd")) + + theme(legend.position = "none") } else { graph <- graph + - scale_fill_manual(values=c("#3182bd")) + - theme(legend.position="none") + scale_fill_manual(values = c("#3182bd")) + + theme(legend.position = "none") } return(suppressWarnings(graph)) @@ -1055,7 +1073,7 @@ setMethod(f="plot_distributions", signature(object="ttest_class"), definition=fu #' # along with an example of how to use this function #' ?ttest_class #' -setMethod(f="plot_distributions_difference", signature(object="ttest_class"), definition=function(object, ...) { +setMethod(f = "plot_distributions_difference", signature(object = "ttest_class"), definition = function(object, ...) { # init local varibales for CRAN check value <- NULL @@ -1078,9 +1096,10 @@ setMethod(f="plot_distributions_difference", signature(object="ttest_class"), de y <- list() n <- 100000 y[[1]] <- metRology::rt.scaled(n, - df=mean(object@extract$nu), - mean=mean(object@extract$mu), - sd=mean(object@extract$sigma)) + df = mean(object@extract$nu), + mean = mean(object@extract$mu), + sd = mean(object@extract$sigma) + ) # limits x_min <- min(y[[1]]) @@ -1095,9 +1114,10 @@ setMethod(f="plot_distributions_difference", signature(object="ttest_class"), de fit2 <- arguments[[1]] } y[[2]] <- metRology::rt.scaled(n, - df=mean(fit2@extract$nu), - mean=mean(fit2@extract$mu), - sd=mean(fit2@extract$sigma)) + df = mean(fit2@extract$nu), + mean = mean(fit2@extract$mu), + sd = mean(fit2@extract$sigma) + ) } else if (!is.null(arguments$mu)) { # provided mu and sigma sigma2 <- 0 @@ -1110,14 +1130,15 @@ setMethod(f="plot_distributions_difference", signature(object="ttest_class"), de } else if (!is.null(arguments$fits)) { i <- 2 for (fit in arguments$fits) { - if (class(fit) != "ttest_class") { + if (!("ttest_class" %in% class(fit))) { stop("One of the fits in the fits list is not a valid ttest_class object.") } y[[i]] <- metRology::rt.scaled(n, - df=mean(fit@extract$nu), - mean=mean(fit@extract$mu), - sd=mean(fit@extract$sigma)) + df = mean(fit@extract$nu), + mean = mean(fit@extract$mu), + sd = mean(fit@extract$sigma) + ) # limits x_min <- min(x_min, y[[i]]) @@ -1138,12 +1159,12 @@ setMethod(f="plot_distributions_difference", signature(object="ttest_class"), de # if no list is provided if (is.null(arguments$fits)) { # call plot difference shared function - graph <- plot_difference(y1=y[[1]], y2=y[[2]], rope=rope, bins=bins) + graph <- plot_difference(y1 = y[[1]], y2 = y[[2]], rope = rope, bins = bins) return(graph) } else { diff <- x_max - x_min - x_min <- x_min - 0.1*diff - x_max <- x_max + 0.1*diff + x_min <- x_min - 0.1 * diff + x_max <- x_max + 0.1 * diff graphs <- list() n <- length(y) @@ -1151,24 +1172,24 @@ setMethod(f="plot_distributions_difference", signature(object="ttest_class"), de for (j in i:n) { # if both are equal plot samples, else plot difference if (i == j) { - df <- data.frame(value=y[[i]]) - index <- (i-1)*n + i + df <- data.frame(value = y[[i]]) + index <- (i - 1) * n + i graphs[[index]] <- ggplot() + - geom_density(data=df, aes(x=value), fill="#3182bd", color=NA, alpha=0.4) + + geom_density(data = df, aes(x = value), fill = "#3182bd", color = NA, alpha = 0.4) + xlab("value") + xlim(x_min, x_max) } else { - index1 <- (i-1)*n + j - graphs[[index1]] <- plot_difference(y1=y[[i]], y2=y[[j]], rope=rope, bins=bins, nrow=n) + index1 <- (i - 1) * n + j + graphs[[index1]] <- plot_difference(y1 = y[[i]], y2 = y[[j]], rope = rope, bins = bins, nrow = n) - index2 <- (j-1)*n + i - graphs[[index2]] <- plot_difference(y1=y[[j]], y2=y[[i]], rope=rope, bins=bins, nrow=n) + index2 <- (j - 1) * n + i + graphs[[index2]] <- plot_difference(y1 = y[[j]], y2 = y[[i]], rope = rope, bins = bins, nrow = n) } } } # cowplot - graph <- suppressWarnings(cowplot::plot_grid(plotlist=graphs, nrow=n, ncol=n, scale=0.9)) + graph <- suppressWarnings(cowplot::plot_grid(plotlist = graphs, nrow = n, ncol = n, scale = 0.9)) return(graph) } }) diff --git a/man/success_rate_class-class.Rd b/man/success_rate_class-class.Rd index 4f3e50b..95ed14d 100644 --- a/man/success_rate_class-class.Rd +++ b/man/success_rate_class-class.Rd @@ -72,29 +72,31 @@ plot_fit(`success_rate_class`): plots fitted model against the data. Use this fu \examples{ \donttest{ -#priors -p_prior <- b_prior(family="beta", pars=c(1, 1)) -tau_prior <- b_prior(family="uniform", pars=c(0, 500)) +# priors +p_prior <- b_prior(family = "beta", pars = c(1, 1)) +tau_prior <- b_prior(family = "uniform", pars = c(0, 500)) # attach priors to relevant parameters -priors <- list(c("p", p_prior), - c("tau", tau_prior)) +priors <- list( + c("p", p_prior), + c("tau", tau_prior) +) # subjects s <- rep(1:5, 20) # generate data and fit -data1 <- rbinom(100, size=1, prob=0.6) -fit1 <- b_success_rate(r=data1, s=s, priors=priors, chains=1) +data1 <- rbinom(100, size = 1, prob = 0.6) +fit1 <- b_success_rate(r = data1, s = s, priors = priors, chains = 1) -data2 <- rbinom(100, size=1, prob=0.1) -fit2 <- b_success_rate(r=data2, s=s, priors=priors, chains=1) +data2 <- rbinom(100, size = 1, prob = 0.1) +fit2 <- b_success_rate(r = data2, s = s, priors = priors, chains = 1) -data3 <- rbinom(100, size=1, prob=0.5) -fit3 <- b_success_rate(r=data3, s=s, priors=priors, chains=1) +data3 <- rbinom(100, size = 1, prob = 0.5) +fit3 <- b_success_rate(r = data3, s = s, priors = priors, chains = 1) -data4 <- rbinom(100, size=1, prob=0.9) -fit4 <- b_success_rate(r=data4, s=s, priors=priors, chains=1) +data4 <- rbinom(100, size = 1, prob = 0.9) +fit4 <- b_success_rate(r = data4, s = s, priors = priors, chains = 1) # fit list fit_list <- list(fit2, fit3, fit4) @@ -112,8 +114,8 @@ plot_fit(fit1) # plot the fitted distribution against the data, # plot on the top (group) level -plot(fit1, subjects=FALSE) -plot_fit(fit1, subjects=FALSE) +plot(fit1, subjects = FALSE) +plot_fit(fit1, subjects = FALSE) # traceplot of the fitted parameters plot_trace(fit1) @@ -125,49 +127,49 @@ parameters <- get_parameters(fit1) subject_parameters <- get_subject_parameters(fit1) # compare means between two fits, use a rope interval -compare_means(fit1, fit2=fit2, rope=0.05) +compare_means(fit1, fit2 = fit2, rope = 0.05) # compare means between multiple fits -compare_means(fit1, fits=fit_list) +compare_means(fit1, fits = fit_list) # visualize difference in means between two fits, # specify number of histogram bins and rope interval -plot_means_difference(fit1, fit2=fit2, bins=40, rope=0.05) +plot_means_difference(fit1, fit2 = fit2, bins = 40, rope = 0.05) # visualize difference in means between multiple fits -plot_means_difference(fit1, fits=fit_list) +plot_means_difference(fit1, fits = fit_list) # visualize means of a single fit plot_means(fit1) # visualize means of two fits -plot_means(fit1, fit2=fit2) +plot_means(fit1, fit2 = fit2) # visualize means of multiple fits -plot_means(fit1, fits=fit_list) +plot_means(fit1, fits = fit_list) # draw samples from distributions underlying two fits and compare them, # use a rope interval -compare_distributions(fit1, fit2=fit2, rope=0.05) +compare_distributions(fit1, fit2 = fit2, rope = 0.05) # draw samples from distributions underlying multiple fits and compare them -compare_distributions(fit1, fits=fit_list) +compare_distributions(fit1, fits = fit_list) # visualize the distribution underlying a fit plot_distributions(fit1) # visualize distributions underlying two fits -plot_distributions(fit1, fit2=fit2) +plot_distributions(fit1, fit2 = fit2) # visualize distributions underlying multiple fits -plot_distributions(fit1, fits=fit_list) +plot_distributions(fit1, fits = fit_list) # visualize difference between distributions underlying two fits, # use a rope interval -plot_distributions_difference(fit1, fit2=fit2, rope=0.05) +plot_distributions_difference(fit1, fit2 = fit2, rope = 0.05) # visualize difference between distributions underlying multiple fits -plot_distributions_difference(fit1, fits=fit_list) +plot_distributions_difference(fit1, fits = fit_list) } } diff --git a/man/ttest_class-class.Rd b/man/ttest_class-class.Rd index 080eb99..10efa8d 100644 --- a/man/ttest_class-class.Rd +++ b/man/ttest_class-class.Rd @@ -85,27 +85,29 @@ plot_distributions_difference(`ttest_class`, fits=`list`): a visualization of th \examples{ \donttest{ # priors -mu_prior <- b_prior(family="normal", pars=c(0, 1000)) -sigma_prior <- b_prior(family="uniform", pars=c(0, 500)) -nu_prior <- b_prior(family="normal", pars=c(2000, 1000)) +mu_prior <- b_prior(family = "normal", pars = c(0, 1000)) +sigma_prior <- b_prior(family = "uniform", pars = c(0, 500)) +nu_prior <- b_prior(family = "normal", pars = c(2000, 1000)) # attach priors to relevant parameters -priors <- list(c("mu", mu_prior), - c("sigma", sigma_prior), - c("nu", nu_prior)) +priors <- list( + c("mu", mu_prior), + c("sigma", sigma_prior), + c("nu", nu_prior) +) # generate data and fit -data1 <- rnorm(20, mean=150, sd=20) -fit1 <- b_ttest(data=data1, priors=priors, chains=1) +data1 <- rnorm(20, mean = 150, sd = 20) +fit1 <- b_ttest(data = data1, priors = priors, chains = 1) -data2 <- rnorm(20, mean=200, sd=20) -fit2 <- b_ttest(data=data2, priors=priors, chains=1) +data2 <- rnorm(20, mean = 200, sd = 20) +fit2 <- b_ttest(data = data2, priors = priors, chains = 1) -data3 <- rnorm(20, mean=150, sd=40) -fit3 <- b_ttest(data=data3, priors=priors, chains=1) +data3 <- rnorm(20, mean = 150, sd = 40) +fit3 <- b_ttest(data = data3, priors = priors, chains = 1) -data4 <- rnorm(20, mean=50, sd=10) -fit4 <- b_ttest(data=data4, priors=priors, chains=1) +data4 <- rnorm(20, mean = 50, sd = 10) +fit4 <- b_ttest(data = data4, priors = priors, chains = 1) # fit list fit_list <- list(fit2, fit3, fit4) @@ -128,86 +130,86 @@ plot_trace(fit1) parameters <- get_parameters(fit1) # compare means between two fits -compare_means(fit1, fit2=fit2) +compare_means(fit1, fit2 = fit2) # compare means between two fits, use a rope interval -compare_means(fit1, fit2=fit2, rope=2) +compare_means(fit1, fit2 = fit2, rope = 2) # compare means between a fit and a constant value -compare_means(fit1, mu=150) +compare_means(fit1, mu = 150) # compare means between a fit and a distribution, # sigma is used for calculating Cohen's d -compare_means(fit1, mu=150, sigma=20) +compare_means(fit1, mu = 150, sigma = 20) # compare means between multiple fits -compare_means(fit1, fits=fit_list) +compare_means(fit1, fits = fit_list) # visualize difference in means between two fits, # specify number of histogram bins -plot_means_difference(fit1, fit2=fit2, bins=20) +plot_means_difference(fit1, fit2 = fit2, bins = 20) # visualize difference in means between a fit and a constant value -plot_means_difference(fit1, mu=150) +plot_means_difference(fit1, mu = 150) # visualize difference in means between multiple fits, use a rope interval -plot_means_difference(fit1, fits=fit_list, rope=2) +plot_means_difference(fit1, fits = fit_list, rope = 2) # visualize means of a single fit plot_means(fit1) # visualize means of two fits -plot_means(fit1, fit2=fit2) +plot_means(fit1, fit2 = fit2) # visualize means of a fit and a constant value -plot_means(fit1, mu=150) +plot_means(fit1, mu = 150) # visualize means of multiple fits -plot_means(fit1, fits=fit_list) +plot_means(fit1, fits = fit_list) # draw samples from distributions underlying two fits and compare them -compare_distributions(fit1, fit2=fit2) +compare_distributions(fit1, fit2 = fit2) # draw samples from a distribution underlying the fit # and compare them with a constant, use a rope interval -compare_distributions(fit1, mu=150, rope=2) +compare_distributions(fit1, mu = 150, rope = 2) # draw samples from a distribution underlying the fit and # compare them with a user defined distribution -compare_distributions(fit1, mu=150, sigma=20) +compare_distributions(fit1, mu = 150, sigma = 20) # draw samples from distributions underlying multiple fits and compare them -compare_distributions(fit1, fits=fit_list) +compare_distributions(fit1, fits = fit_list) # visualize the distribution underlying a fit plot_distributions(fit1) # visualize distributions underlying two fits -plot_distributions(fit1, fit2=fit2) +plot_distributions(fit1, fit2 = fit2) # visualize the distribution underlying a fit and a constant value -plot_distributions(fit1, mu=150) +plot_distributions(fit1, mu = 150) # visualize the distribution underlying a fit and a user defined distribution -plot_distributions(fit1, mu=150, sigma=20) +plot_distributions(fit1, mu = 150, sigma = 20) # visualize distributions underlying multiple fits -plot_distributions(fit1, fits=fit_list) +plot_distributions(fit1, fits = fit_list) # visualize difference between distributions underlying two fits, # use a rope interval -plot_distributions_difference(fit1, fit2=fit2, rope=2) +plot_distributions_difference(fit1, fit2 = fit2, rope = 2) # visualize difference between a distribution underlying the fit # and a constant value -plot_distributions_difference(fit1, mu=150) +plot_distributions_difference(fit1, mu = 150) # visualize difference between a distribution underlying the fits # and a user defined distribution -plot_distributions_difference(fit1, mu=150, sigma=20) +plot_distributions_difference(fit1, mu = 150, sigma = 20) # visualize difference between distributions underlying multiple fits -plot_distributions_difference(fit1, fits=fit_list) +plot_distributions_difference(fit1, fits = fit_list) } } From ef815e5c8ec8e6c2ea2d326184a17032f9771de8 Mon Sep 17 00:00:00 2001 From: Jure Demsar Date: Wed, 22 Mar 2023 17:31:30 +0100 Subject: [PATCH 10/13] Fixed a http => https issue. --- .Rbuildignore | 1 + CRAN-SUBMISSION | 3 +++ R/b_bootstrap.R | 20 ++++++++++---------- 3 files changed, 14 insertions(+), 10 deletions(-) create mode 100644 CRAN-SUBMISSION diff --git a/.Rbuildignore b/.Rbuildignore index e206297..966aa16 100644 --- a/.Rbuildignore +++ b/.Rbuildignore @@ -7,3 +7,4 @@ cran-comments.md ^Meta$ ^CRAN-RELEASE$ ^\.github$ +^CRAN-SUBMISSION$ diff --git a/CRAN-SUBMISSION b/CRAN-SUBMISSION new file mode 100644 index 0000000..eabbd28 --- /dev/null +++ b/CRAN-SUBMISSION @@ -0,0 +1,3 @@ +Version: 1.2.11 +Date: 2023-03-22 14:49:36 UTC +SHA: 9a178a8cf7170b102787deb436b9ed5b373f482c diff --git a/R/b_bootstrap.R b/R/b_bootstrap.R index cfa39c4..3d468b0 100644 --- a/R/b_bootstrap.R +++ b/R/b_bootstrap.R @@ -1,7 +1,7 @@ #' @title b_bootstrap #' @description Performs a Bayesian bootstrap and returns a sample of size n1 representing the posterior distribution of the statistic. Returns a vector if the statistic is one-dimensional (like for mean(...)) or a data.frame if the statistic is multi-dimensional (like for the coefficients of lm). #' @author Rasmus Baath -#' @references \url{http://www.sumsar.net/blog/2015/07/easy-bayesian-bootstrap-in-r/} +#' @references \url{https://www.sumsar.net/blog/2015/07/easy-bayesian-bootstrap-in-r/} #' @references Rubin, D. B. (1981). The Bayesian Bootstrap. The annals of statistics, 9(1), 130-134. #' @export #' @param data The data as either a vector, matrix or data.frame. @@ -24,35 +24,35 @@ #' data <- adaptation_level_small #' #' # bootstrap -#' data_bootstrap <- b_bootstrap(data, lm_statistic, n1=1000, n2=1000) +#' data_bootstrap <- b_bootstrap(data, lm_statistic, n1 = 1000, n2 = 1000) #' -b_bootstrap <- function(data, statistic, n1=1000, n2=1000, use_weights=FALSE, weight_arg=NULL, ...) { +b_bootstrap <- function(data, statistic, n1 = 1000, n2 = 1000, use_weights = FALSE, weight_arg = NULL, ...) { # Draw from a uniform Dirichlet dist. with alpha set to rep(1, n_dim). # Using the facts that you can transform gamma distributed draws into # Dirichlet draws and that rgamma(n, 1) <=> rexp(n, 1) - dirichlet_weights <- matrix(stats::rexp(NROW(data) * n1, 1) , ncol=NROW(data), byrow=TRUE) + dirichlet_weights <- matrix(stats::rexp(NROW(data) * n1, 1), ncol = NROW(data), byrow = TRUE) dirichlet_weights <- dirichlet_weights / rowSums(dirichlet_weights) - if(use_weights) { + if (use_weights) { stat_call <- quote(statistic(data, w, ...)) names(stat_call)[3] <- weight_arg boot_sample <- apply(dirichlet_weights, 1, function(w) { eval(stat_call) }) } else { - if(is.null(dim(data)) || length(dim(data)) < 2) { # data is a list type of object + if (is.null(dim(data)) || length(dim(data)) < 2) { # data is a list type of object boot_sample <- apply(dirichlet_weights, 1, function(w) { - data_sample <- sample(data, size=n2, replace=TRUE, prob=w) + data_sample <- sample(data, size = n2, replace = TRUE, prob = w) statistic(data_sample, ...) }) } else { # data is a table type of object boot_sample <- apply(dirichlet_weights, 1, function(w) { - index_sample <- sample(nrow(data), size=n2, replace=TRUE, prob=w) - statistic(data[index_sample, , drop=FALSE], ...) + index_sample <- sample(nrow(data), size = n2, replace = TRUE, prob = w) + statistic(data[index_sample, , drop = FALSE], ...) }) } } - if(is.null(dim(boot_sample)) || length(dim(boot_sample)) < 2) { + if (is.null(dim(boot_sample)) || length(dim(boot_sample)) < 2) { # If the bootstrap sample is just a simple vector return it. boot_sample } else { From 6080a4b5f55024ce2adbf776a463c45ccd7c5d7c Mon Sep 17 00:00:00 2001 From: Jure Demsar Date: Wed, 22 Mar 2023 18:00:39 +0100 Subject: [PATCH 11/13] CRAN submission --- CRAN-SUBMISSION | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/CRAN-SUBMISSION b/CRAN-SUBMISSION index eabbd28..639277b 100644 --- a/CRAN-SUBMISSION +++ b/CRAN-SUBMISSION @@ -1,3 +1,3 @@ Version: 1.2.11 -Date: 2023-03-22 14:49:36 UTC -SHA: 9a178a8cf7170b102787deb436b9ed5b373f482c +Date: 2023-03-22 17:00:29 UTC +SHA: ef815e5c8ec8e6c2ea2d326184a17032f9771de8 From ebd5258c187e9b08ce60238e50168d25c2b5e910 Mon Sep 17 00:00:00 2001 From: Jure Demsar Date: Wed, 22 Mar 2023 20:47:46 +0100 Subject: [PATCH 12/13] Documentation update. --- man/b_bootstrap.Rd | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/man/b_bootstrap.Rd b/man/b_bootstrap.Rd index 7f9785b..5992608 100644 --- a/man/b_bootstrap.Rd +++ b/man/b_bootstrap.Rd @@ -46,11 +46,11 @@ lm_statistic <- function(data) { data <- adaptation_level_small # bootstrap -data_bootstrap <- b_bootstrap(data, lm_statistic, n1=1000, n2=1000) +data_bootstrap <- b_bootstrap(data, lm_statistic, n1 = 1000, n2 = 1000) } \references{ -\url{http://www.sumsar.net/blog/2015/07/easy-bayesian-bootstrap-in-r/} +\url{https://www.sumsar.net/blog/2015/07/easy-bayesian-bootstrap-in-r/} Rubin, D. B. (1981). The Bayesian Bootstrap. The annals of statistics, 9(1), 130-134. } From 971eacdb145db76c7a731533e5baa56144030cc8 Mon Sep 17 00:00:00 2001 From: Jure Demsar Date: Fri, 24 Mar 2023 14:55:33 +0100 Subject: [PATCH 13/13] CRAN SUBMISSION --- CRAN-SUBMISSION | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/CRAN-SUBMISSION b/CRAN-SUBMISSION index 639277b..2b151cc 100644 --- a/CRAN-SUBMISSION +++ b/CRAN-SUBMISSION @@ -1,3 +1,3 @@ Version: 1.2.11 -Date: 2023-03-22 17:00:29 UTC -SHA: ef815e5c8ec8e6c2ea2d326184a17032f9771de8 +Date: 2023-03-22 19:58:11 UTC +SHA: ebd5258c187e9b08ce60238e50168d25c2b5e910