Testing bmd package

[Zhenglei-BCS]

I am trying to add more documentation to bmd while testing and validating the package.

bmd pull request

Suggestions for Improvement

Documentation Improvements

1. Comprehensive Vignettes:

   • Create step-by-step tutorials for common workflows
   • Include real-world examples for different response types
   • Provide guidance on model selection and interpretation

2. Function Documentation:

   • Add more detailed examples to each function's documentation
   • Clarify parameter choices and their implications
   • Include references to statistical literature for advanced methods

3. Theoretical Background:

   • Add more documentation explaining the statistical theory behind different BMD definitions
   • Explain when to use different confidence interval methods
   • Provide guidance on interpreting results

Code Structure and Functionality

1. Code Organization:

   • Refactor code to improve consistency in function naming and parameter order
   • Group related functions into coherent files rather than having one function per file
   • Implement S3 methods consistently across all object types

2. Error Handling:

   • Improve error messages to be more informative and actionable
   • Add input validation for major functions
   • Implement warning messages for potentially problematic parameter combinations

3. Performance Optimization:

   • Optimize computationally intensive functions (especially bootstrapping)
   • Implement parallel processing options for bootstrap procedures
   • Consider using Rcpp for performance-critical calculations

New Features

1. Expanded Model Types:

   • Add support for hormetic/non-monotonic dose-response models
   • Implement Bayesian dose-response models (brms ?)
   • Add time-to-event models for survival data

2. Model Diagnostics:

   • Add comprehensive model diagnostic functions
   • Implement goodness-of-fit tests specific to dose-response models
   • Add residual analysis tools for different response types

3. Export Capabilities:

   • Add functions to export results in formats compatible with regulatory submissions (maybe for drcHelper)

4. Other thoughts

   • Develop interfaces to toxicological databases

Technical Improvements

1. Testing Framework:

   • Implement comprehensive unit tests for all functions
   • Add integration tests for common workflows
   • Create benchmark tests for performance-critical functions

2. Documentation Generation:

   • Use roxygen2 consistently for all documentation
   • Include mathematical formulas in documentation using LaTeX
   • Generate a package website using pkgdown

[Zhenglei-BCS]

Issues when run tests:

• 1. require drc package and not automatically loaded.

• 2. encounter issues when using sandwich methods. $\rightarrow$ might be an issue caused by drc package version, this is dependent on 3.2.0

Backtrace:
    ▆
 1. └─base::library(drc) at test-bmd.R:33:1
 2.   └─base (local) testRversion(pkgInfo, package, pkgpath)

Error (test-bmd.R:162:3): bmd function computes correct sandwich confidence interval for a simple model
Error in `UseMethod("estfun")`: no applicable method for 'estfun' applied to an object of class "drc"
Backtrace:
    ▆
 1. └─bmd:::bmd(...) at test-bmd.R:162:3
 2.   └─sandwich::sandwich(object) at bmd/R/bmd.R:258:9
 3.     └─sandwich (local) meat.(x, ...)
 4.       └─sandwich::estfun(x, ...)

Error (test-bmd.R:179:3): bmd function computes BMD (point) correctly for ryegrass model
Error in `UseMethod("estfun")`: no applicable method for 'estfun' applied to an object of class "drc"
Backtrace:
    ▆
 1. └─bmd:::bmd(...) at test-bmd.R:179:3
 2.   └─sandwich::sandwich(object) at bmd/R/bmd.R:258:9
 3.     └─sandwich (local) meat.(x, ...)
 4.       └─sandwich::estfun(x, ...)

Error (test-bmd.R:206:3): bmd function computes BMD (extra) correctly for ryegrass model
Error in `UseMethod("estfun")`: no applicable method for 'estfun' applied to an object of class "drc"
Backtrace:
    ▆
 1. └─bmd:::bmd(...) at test-bmd.R:206:3
 2.   └─sandwich::sandwich(object) at bmd/R/bmd.R:258:9
 3.     └─sandwich (local) meat.(x, ...)
 4.       └─sandwich::estfun(x, ...)

Error (test-bmd.R:238:3): bmd function computes BMD (relative) correctly for ryegrass model
Error in `UseMethod("estfun")`: no applicable method for 'estfun' applied to an object of class "drc"
Backtrace:
    ▆
 1. └─bmd:::bmd(...) at test-bmd.R:238:3
 2.   └─sandwich::sandwich(object) at bmd/R/bmd.R:258:9
 3.     └─sandwich (local) meat.(x, ...)
 4.       └─sandwich::estfun(x, ...)

════

• 3. encounter issues with bmdBoot and bmdMA, tests have benn passed though. $\rightarrow$ make sure with Jens this is not an issue.

⠸ |        224 | bmdBoot                                                                                                Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt,  : 
  non-finite finite-difference value [1]
Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt,  : 
  non-finite finite-difference value [6]
Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt,  : 
  non-finite finite-difference value [7]
Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt,  : 
  non-finite finite-difference value [6]
Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt,  : 
  non-finite finite-difference value [7]
Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt,  : 
  non-finite finite-difference value [6]
Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt,  : 
  non-finite finite-difference value [1]
⠼ |        225 | bmdBoot                                                                                                Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt,  : 
  non-finite finite-difference value [1]
Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt,  : 
  non-finite finite-difference value [6]
Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt,  : 
  non-finite finite-difference value [7]
Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt,  : 
  non-finite finite-difference value [6]
Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt,  : 
  non-finite finite-difference value [7]
Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt,  : 
  non-finite finite-difference value [6]
Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt,  : 
  non-finite finite-difference value [1]
⠼ |        235 | bmdBoot                                                                                                Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt,  : 
  non-finite finite-difference value [1]
Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt,  : 
  non-finite finite-difference value [6]
Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt,  : 
  non-finite finite-difference value [7]
Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt,  : 
  non-finite finite-difference value [6]
✔ |        249 | bmdBoot [162.6s]

[Zhenglei-BCS]

• NAMESPACE issue:

Do we need to export all these functions? Just to confirm.

export(bmd, bmdBoot, bmdIso, bmdIsoBoot, PAV, bmdMA, bootDataGen, bmdMACurve, BCa, invBmd, expandBinomial,
getStackingWeights, drmOrdinal, bmdOrdinal, bmdOrdinalMA,
qplotDrc, qplotBmd, MACurve, monotonicityTest, trendTest, bmdHetVar, drmHetVar, drmMMRE, bmdHetVarMA)

[Zhenglei-BCS]

• Direct run examples (devtools::run_examples()) would result in missing package. This is solved.

> ## Model-averaged BMD using the Stacking Weights
> bmdMA(list(ryegrass.m1,ryegrass.m2,ryegrass.m3,ryegrass.m4), modelWeights="Stack", bmr=0.05, 
+       backgType = "hybridSD", def = "hybridAdd", type = "bootstrap",
+       bootstrapType = "nonparametric", bootInterval = "percentile", R = 50)
Loading required namespace: CVXR
Failed with error:  ‘there is no package called ‘CVXR’’
Error in getStackingWeights(modelList, stackingSplits) : 
  package "CVXR" must be installed to estimate stacking weights

[Zhenglei-BCS]

• Annoying warnings:

plot(bmd0, type = "confidence")
#> Warning in (tquan * sqrt(varVal + sumObjRV)) * c(-1, 1): Recycling array of length 1 in array-vector arithmetic is deprecated.
#>   Use c() or as.vector() instead.
#> Warning in (tquan * sqrt(varVal + sumObjRV)) * c(-1, 1): Recycling array of length 1 in array-vector arithmetic is deprecated.
#>   Use c() or as.vector() instead.
#> Warning in (tquan * sqrt(varVal + sumObjRV)) * c(-1, 1): Recycling array of length 1 in array-vector arithmetic is deprecated.
#>   Use c() or as.vector() instead.
#> Warning in (tquan * sqrt(varVal + sumObjRV)) * c(-1, 1): Recycling array of length 1 in array-vector arithmetic is deprecated.
#>   Use c() or as.vector() instead.

## Plot with confidence regions

qplotDrc(spinach.m1, col = TRUE, type = "confidence")
#> Warning in (tquan * sqrt(varVal + sumObjRV)) * c(-1, 1): Recycling array of length 1 in array-vector arithmetic is deprecated.
#>   Use c() or as.vector() instead.
#> Warning in (tquan * sqrt(varVal + sumObjRV)) * c(-1, 1): Recycling array of length 1 in array-vector arithmetic is deprecated.
#>   Use c() or as.vector() instead.
#> Warning in (tquan * sqrt(varVal + sumObjRV)) * c(-1, 1): Recycling array of length 1 in array-vector arithmetic is deprecated.
#>   Use c() or as.vector() instead.
#> Warning in (tquan * sqrt(varVal + sumObjRV)) * c(-1, 1): Recycling array of length 1 in array-vector arithmetic is deprecated.
#>   Use c() or as.vector() instead.
#> Warning in (tquan * sqrt(varVal + sumObjRV)) * c(-1, 1): Recycling array of length 1 in array-vector arithmetic is deprecated.
#>   Use c() or as.vector() instead.
#> Warning in (tquan * sqrt(varVal + sumObjRV)) * c(-1, 1): Recycling array of length 1 in array-vector arithmetic is deprecated.
#>   Use c() or as.vector() instead.
#> Warning in (tquan * sqrt(varVal + sumObjRV)) * c(-1, 1): Recycling array of length 1 in array-vector arithmetic is deprecated.
#>   Use c() or as.vector() instead.
#> Warning in (tquan * sqrt(varVal + sumObjRV)) * c(-1, 1): Recycling array of length 1 in array-vector arithmetic is deprecated.
#>   Use c() or as.vector() instead.
#> Warning in (tquan * sqrt(varVal + sumObjRV)) * c(-1, 1): Recycling array of length 1 in array-vector arithmetic is deprecated.
#>   Use c() or as.vector() instead.
#> Warning in (tquan * sqrt(varVal + sumObjRV)) * c(-1, 1): Recycling array of length 1 in array-vector arithmetic is deprecated.
#>   Use c() or as.vector() instead.
#> Warning in (tquan * sqrt(varVal + sumObjRV)) * c(-1, 1): Recycling array of length 1 in array-vector arithmetic is deprecated.
#>   Use c() or as.vector() instead.

[Zhenglei-BCS]

As a record, I have used these steps to replace manual rd with roxygen style and added vignettes and package github pages.

Rd2roxygen

# Install the package if not already installed
install.packages("Rd2roxygen")

# Load the package
library(Rd2roxygen)

# Convert all Rd files in a package to roxygen comments
Rd2roxygen("path/to/your/package")

Update DESCRIPTION file to inlcude

Roxygen: list(markdown = TRUE)
RoxygenNote: 7.2.3  # or whatever version you're using

Remove empty tage

# This script will remove empty tags from your R files
r_files <- list.files("R", pattern = "\\.R$", full.names = TRUE)
for (file in r_files) {
  lines <- readLines(file)
  # Remove lines with empty tags
  empty_tag_pattern <- "#'\\s+@(note|seealso|references)\\s*$"
  lines <- lines[!grepl(empty_tag_pattern, lines)]
  writeLines(lines, file)
}

Delete manual Rd and regenerate them

# Step 1: Remove all .Rd files
unlink("man/*.Rd")

# Step 2: Fix empty tags
r_files <- list.files("R", pattern = "\\.R$", full.names = TRUE)
for (file in r_files) {
  lines <- readLines(file)
  # Remove lines with empty tags
  empty_tag_pattern <- "#'\\s+@(note|seealso|references|examples)\\s*$"
  lines <- lines[!grepl(empty_tag_pattern, lines)]
  writeLines(lines, file)
}

# Step 3: Fix S3 method exports
s3_methods <- c(
  "print.drcMMRE", "print.drcHetVar", "logLik.drcOrdinal", 
  "print.bmdHetVar", "logLik.drcHetVar", "plot.drcHetVar",
  "AIC.drcOrdinal", "vcov.drcMMRE", "BIC.drcHetVar", 
  "print.bmdOrdinal", "coef.drcHetVar", "AIC.drcHetVar",
  "BIC.drcOrdinal", "print.drcOrdinal", "print.bmd", 
  "plot.drcOrdinal"
)

for (method in s3_methods) {
  # Extract function name and class
  parts <- strsplit(method, "\\.")[[1]]
  func_name <- parts[1]
  class_name <- parts[2]
  
  # Find the R file containing this method
  r_file <- paste0("R/", method, ".R")
  if (file.exists(r_file)) {
    lines <- readLines(r_file)
    # Add @export tag before the function definition
    func_line <- grep(paste0(method, "\\s*(<-|=)\\s*function"), lines)[1]
    if (!is.na(func_line)) {
      # Check if there are already roxygen comments
      if (func_line > 1 && grepl("#'", lines[func_line-1])) {
        # Add @export tag to existing roxygen block
        export_line <- max(grep("#'", lines[1:(func_line-1)]))
        lines <- c(
          lines[1:export_line],
          "#' @export",
          lines[(export_line+1):length(lines)]
        )
      } else {
        # Create new roxygen block with @export
        lines <- c(
          lines[1:(func_line-1)],
          "#' @title S3 method",
          "#' @export",
          lines[func_line:length(lines)]
        )
      }
      writeLines(lines, r_file)
      cat("Added @export tag to", method, "in", r_file, "\n")
    }
  }
}

# Step 4: Update NAMESPACE
file.copy("NAMESPACE", "NAMESPACE.backup", overwrite = TRUE)
writeLines("# Generated by roxygen2: do not edit by hand", "NAMESPACE")

# Step 5: Update DESCRIPTION file
desc <- readLines("DESCRIPTION")
if (!any(grepl("^Roxygen:", desc))) {
  desc <- c(desc, "", "Roxygen: list(markdown = TRUE)", 
            paste0("RoxygenNote: ", packageVersion("roxygen2")))
  writeLines(desc, "DESCRIPTION")
}

# Step 6: Regenerate documentation
devtools::document()

Note that all inline documentation need to be checked to add at export tag to be consistent with the manual written NAMESPACE

Moving detailed examples into vignettes

#' @examples
#' # See vignette("function-examples") for detailed examples
#' 
#' # Simple example:
#' x <- c(0, 1, 2, 4, 8)
#' y <- c(1.0, 0.8, 0.6, 0.4, 0.2)
#' if (requireNamespace("drc", quietly = TRUE)) {
#'   model <- drc::drm(y ~ x, fct = drc::LL.4())
#'   result <- bmd(model, bmr = 0.1, def = "relative")
#'   # Print result summary
#'   print(result)
#' }

### Function categories

Define function categories

function_categories <- list(
"Basic BMD Functions" = c("bmd", "bmdBoot", "invBmd"),
"Model Averaging" = c("bmdMA", "bmdMACurve", "MACurve", "getStackingWeights"),
"Specialized Models" = c("drmOrdinal", "bmdOrdinal", "drmHetVar", "bmdHetVar", "drmMMRE"),
"Visualization" = c("qplotDrc", "qplotBmd", "plot.bmd"),
"Statistical Tests" = c("monotonicityTest", "trendTest"),
"Utilities" = c("PAV", "expandBinomial", "bootDataGen", "BCa")
)

[Zhenglei-BCS]

Functions I used to extact all examples and make them into .Rmd

# Step 1: Extract examples from your package
examples <- extract_package_examples(package = ".", dir = "extracted_examples")
# OR
examples <- extract_examples_using_example(package = ".", dir = "extracted_examples")

# Step 2: Create a vignette with all examples
create_examples_vignette(examples_dir = "extracted_examples", vignette_name = "function-examples")

extract_package_examples <- function(package = ".", dir = "extracted_examples") {
  # Create directory for extracted examples if it doesn't exist
  if (!dir.exists(dir)) {
    dir.create(dir)
  }
  
  # Get package path
  pkg_path <- devtools::as.package(package)$path
  
  # Get all Rd files
  rd_files <- list.files(file.path(pkg_path, "man"), pattern = "\\.Rd$", full.names = TRUE)
  
  # Process each Rd file
  examples_list <- list()
  
  for (rd_file in rd_files) {
    # Parse the Rd file
    rd <- tools::parse_Rd(rd_file)
    
    # Extract function name
    name_idx <- which(sapply(rd, attr, "Rd_tag") == "\\name")
    if (length(name_idx) == 0) next
    func_name <- paste(rd[[name_idx]][[1]], collapse = "")
    
    # Extract examples
    example_idx <- which(sapply(rd, attr, "Rd_tag") == "\\examples")
    if (length(example_idx) == 0) next
    
    # Get example text
    example_text <- paste(sapply(rd[[example_idx]], as.character), collapse = "\n")
    
    # Clean up the example text
    example_text <- gsub("\\\\dont(run|test)\\{", "", example_text)
    example_text <- gsub("\\}", "", example_text)
    
    # Store in list
    if (nchar(example_text) > 0) {
      examples_list[[func_name]] <- example_text
      
      # Write to file
      r_file <- file.path(dir, paste0(func_name, "_example.R"))
      writeLines(example_text, r_file)
      
      # Create Rmd snippet
      rmd_content <- c(
        paste0("## Example for `", func_name, "` function"),
        "",
        "```r",
        example_text,
        "```",
        ""
      )
      rmd_file <- file.path(dir, paste0(func_name, "_example.Rmd"))
      writeLines(rmd_content, rmd_file)
    }
  }
  
  # Create a combined file with all examples
  if (length(examples_list) > 0) {
    all_examples <- character()
    for (func_name in names(examples_list)) {
      all_examples <- c(
        all_examples,
        paste0("### Examples for ", func_name),
        "",
        examples_list[[func_name]],
        "",
        "# -------------------------------------------------------------",
        ""
      )
    }
    writeLines(all_examples, file.path(dir, "all_examples.R"))
  }
  
  return(examples_list)
}

create_examples_vignette <- function(examples_dir = "extracted_examples", 
                                    vignette_name = "function-examples") {
  # Create vignette
  usethis::use_vignette(vignette_name)
  
  # Get all Rmd files
  example_files <- list.files(examples_dir, pattern = "_example\\.Rmd$", full.names = TRUE)
  
  if (length(example_files) == 0) {
    stop("No example Rmd files found in ", examples_dir)
  }
  
  # Read vignette template
  vignette_file <- file.path("vignettes", paste0(vignette_name, ".Rmd"))
  vignette_content <- readLines(vignette_file)
  
  # Find where to insert examples (after setup chunk)
  setup_chunks <- grep("```\\\\s*\$", vignette_content)
  if (length(setup_chunks) >= 2) {
    setup_end <- setup_chunks[2]
  } else {
    setup_end <- length(vignette_content)
  }
  
  # Prepare combined examples content
  combined_examples <- character()
  for (file in example_files) {
    file_content <- readLines(file)
    combined_examples <- c(combined_examples, file_content, "", "---", "")
  }
  
  # Insert examples
  new_content <- c(
    vignette_content[1:setup_end],
    "",
    "# Function Examples",
    "",
    "This vignette provides examples for all functions in the package.",
    "",
    combined_examples,
    ""
  )
  
  # Write updated vignette
  writeLines(new_content, vignette_file)
  
  message("Created vignette: ", vignette_file)
}

[Zhenglei-BCS]

Summary of Repository Transformation to Roxygen Style with Devtools Workflow

1. Documentation Conversion

Converting .Rd Files to Roxygen Comments

• Used Rd2roxygen package to convert existing .Rd documentation files to roxygen2-style comments
• Placed roxygen comments directly above corresponding function definitions in R source files
• Fixed issues with empty tags and incomplete conversions manually

Enhancing Roxygen Documentation

• Added proper @param, @return, @examples, and other tags to all functions
• Fixed S3 method documentation by adding proper @export tags

Package-Level Documentation

• Created proper package documentation using the "_PACKAGE" approach
• Fixed duplicate alias issues between package documentation and function documentation
• Added comprehensive package description, imports, and metadata

2. NAMESPACE Management

Converting Static NAMESPACE to Roxygen-Generated

• Replaced the manually maintained NAMESPACE file with roxygen-generated one
• Added proper @export tags to all exported functions
• Added @importFrom and @import tags for package dependencies
• Ensured S3 methods were properly exported with @export tags

3. Vignettes Management

Converting Examples to Vignettes

• Created vignettes to showcase package functionality by extracting examples from existing documentation
• Potentially we could move complex examples from function documentation to vignettes
• Fixed example execution issues in documentation

4. Package Structure Modernization

Implementing Devtools Workflow

• Adopted devtools::document() for documentation generation

• Used devtools::check() for package validation

• Used devtools::test() for running unit tests

• Set up devtools::build_vignettes() for vignette building

• R folder can be organized by group functions with similar application scenario together.

5. GitHub Pages Integration

Setting Up Package Website

• Used pkgdown to generate a package website
• In the future we could configure GitHub Actions to automatically build and deploy the website

This transformation has modernized the package development workflow, making it easier to maintain, extend, and document the package in the future. The roxygen2-based documentation system ensures that code and documentation stay in sync, while the devtools workflow simplifies package development tasks. The GitHub Pages integration provides accessible documentation for users of the package.

[Zhenglei-BCS]

This has been merged by Jens.