Skip to content

Commit

Permalink
skpr 0.61.6: fix case in gen_design() where order of levels in the …
Browse files Browse the repository at this point in the history
…formula is not identical between the subplots and the whole plots
  • Loading branch information
tylermorganwall committed Oct 24, 2019
1 parent ded8821 commit a80104b
Show file tree
Hide file tree
Showing 6 changed files with 84 additions and 7 deletions.
2 changes: 1 addition & 1 deletion DESCRIPTION
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
Package: skpr
Title: Design of Experiments Suite: Generate and Evaluate Optimal Designs
Version: 0.61.5
Version: 0.61.6
Authors@R: c(person("Tyler", "Morgan-Wall", email = "[email protected]", role = c("aut", "cre")),
person("George", "Khoury", email = "[email protected]", role = c("aut")))
Description: Generates and evaluates D, I, A, Alias, E, T, and G optimal designs. Supports generation and evaluation of split/split-split/.../N-split plot designs. Includes parametric and Monte Carlo power evaluation functions, and supports calculating power for censored responses. Provides a framework to evaluate power using functions provided in other packages or written by the user. Includes a Shiny graphical user interface that displays the underlying code used to create and evaluate the design to improve ease-of-use and make analyses more reproducible.
Expand Down
8 changes: 4 additions & 4 deletions R/eval_design_mc.R
Original file line number Diff line number Diff line change
Expand Up @@ -748,18 +748,18 @@ eval_design_mc = function(design, model, alpha,
attr(retval, "I") = IOptimality(modelmatrix_cor, momentsMatrix = mm, blockedVar = diag(nrow(modelmatrix_cor)))
deffic = DOptimality(modelmatrix_cor)
if(!is.infinite(deffic)) {
attr(results, "D") = 100 * DOptimality(modelmatrix_cor) ^ (1 / ncol(modelmatrix_cor)) / nrow(modelmatrix_cor)
attr(retval, "D") = 100 * DOptimality(modelmatrix_cor) ^ (1 / ncol(modelmatrix_cor)) / nrow(modelmatrix_cor)
} else {
attr(results, "D") = 100 * DOptimalityLog(modelmatrix_cor) ^ (1 / ncol(modelmatrix_cor)) / nrow(modelmatrix_cor)
attr(retval, "D") = 100 * DOptimalityLog(modelmatrix_cor) ^ (1 / ncol(modelmatrix_cor)) / nrow(modelmatrix_cor)
}
} else {
attr(retval, "variance.matrix") = V
attr(retval, "I") = IOptimality(modelmatrix_cor, momentsMatrix = mm, blockedVar = V)
deffic = DOptimalityBlocked(modelmatrix_cor, blockedVar = V)
if(!is.infinite(deffic)) {
attr(results, "D") = 100 * DOptimalityBlocked(modelmatrix_cor, blockedVar = V) ^ (1 / ncol(modelmatrix_cor)) / nrow(modelmatrix_cor)
attr(retval, "D") = 100 * DOptimalityBlocked(modelmatrix_cor, blockedVar = V) ^ (1 / ncol(modelmatrix_cor)) / nrow(modelmatrix_cor)
} else {
attr(results, "D") = 100 * DOptimalityBlockedLog(modelmatrix_cor, blockedVar = V) ^ (1 / ncol(modelmatrix_cor)) / nrow(modelmatrix_cor)
attr(retval, "D") = 100 * DOptimalityBlockedLog(modelmatrix_cor, blockedVar = V) ^ (1 / ncol(modelmatrix_cor)) / nrow(modelmatrix_cor)
}
}
if(alpha_adjust) {
Expand Down
9 changes: 7 additions & 2 deletions R/gen_design.R
Original file line number Diff line number Diff line change
Expand Up @@ -424,7 +424,6 @@ gen_design = function(candidateset, model, trials,
interactioncounter = 1
interactionlist = list()
#get model matrix of everything except whole/interactions

for (interaction_col in interactionnames) {
term_vals = unlist(strsplit(interaction_col, split = "(\\s\\*\\s)|(:)", perl = TRUE))
if (any(term_vals %in% submm)) {
Expand Down Expand Up @@ -805,7 +804,13 @@ gen_design = function(candidateset, model, trials,
blockedFactors = c(colnames(blockedmodelmatrix), colnames(candidatesetmm)[-1])
blockedmm = gen_momentsmatrix(blockedFactors, levelvector, classvector)
} else {
interactionnames = interactionnames[!(interactionnames %in% colnames(blockedmodelmatrix))]
blocked_interactions = colnames(blockedmodelmatrix)[grepl(":",colnames(blockedmodelmatrix),fixed=TRUE)]
potential_blocked_interactions = c()
for(i in seq_len(length(blocked_interactions))) {
potential_blocked_interactions = c(potential_blocked_interactions,
potential_permuted_factors(unlist(strsplit(blocked_interactions[[i]],":"))))
}
interactionnames = interactionnames[!(interactionnames %in% potential_blocked_interactions)]
blockedFactors = c(colnames(blockedmodelmatrix), colnames(candidatesetmm)[-1], interactionnames)
blockedmm = gen_momentsmatrix(blockedFactors, levelvector, classvector)
}
Expand Down
36 changes: 36 additions & 0 deletions R/permutation_functions.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,36 @@
#'@title Permutations
#'
#'@description Return permutations
#'
#'@param n Number of elements to permute
#'@keywords internal
#'@return Matrix of permuted element ids
permutations = function(n){
if(n == 1){
return(matrix(1))
} else {
sp = permutations(n-1)
p = nrow(sp)
A = matrix(nrow=n*p,ncol=n)
for(i in 1:n){
A[(i-1)*p+1:p,] = cbind(i,sp+(sp>=i))
}
return(A)
}
}

#'@title Find potential permuted interactions
#'
#'@description Returns permuted interactions
#'
#'@param x character vector of interaction terms
#'@keywords internal
#'@return character vector of potential interaction terms
potential_permuted_factors = function(x) {
numberterms = length(x)
if(numberterms == 1) {
return(x)
}
return(unique(unlist(apply(matrix(x[permutations(numberterms)],ncol=numberterms),1, paste,collapse=":"))))
}

18 changes: 18 additions & 0 deletions man/permutations.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

18 changes: 18 additions & 0 deletions man/potential_permuted_factors.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

0 comments on commit a80104b

Please sign in to comment.