Errors in running code in vignettes:
when running code in ‘free1way.Rnw’
  ...

> yr <- format(dt <- as.Date(packageDescription("free1way.docreg")$Date), 
+     "%Y")

> vs <- packageDescription("free1way.docreg")$Version

> title <- "Semiparametrically Efficient Population and Permutation Inference in \n       Distribution-free Stratified $K$-sample Oneway Layouts"

> DOI <- paste0("10.32614/CRAN.package.", packageDescription("free1way.docreg")$Package)

> Nsim <- 100

> options(digits = 5)

> .table2list <- function(x) {
+     dx <- dim(x)
+     if (length(dx) == 1) 
+         stop("incorrect dimensions")
+     if (length(dx) == 2) 
+     .... [TRUNCATED] 

> .nll <- function(parm, x, mu = 0, rightcensored = FALSE) {
+     bidx <- seq_len(ncol(x) - 1)
+     delta <- c(0, mu + parm[bidx])
+     intercepts  .... [TRUNCATED] 

> .nsc <- function(parm, x, mu = 0, rightcensored = FALSE) {
+     bidx <- seq_len(ncol(x) - 1)
+     delta <- c(0, mu + parm[bidx])
+     intercepts  .... [TRUNCATED] 

> .nsr <- function(parm, x, mu = 0, rightcensored = FALSE) {
+     bidx <- seq_len(ncol(x) - 1)
+     delta <- c(0, mu + parm[bidx])
+     intercepts  .... [TRUNCATED] 

> .hes <- function(parm, x, mu = 0, rightcensored = FALSE, 
+     full = FALSE) {
+     bidx <- seq_len(ncol(x) - 1)
+     delta <- c(0, mu + parm[bid .... [TRUNCATED] 

> .snll <- function(parm, x, mu = 0, rightcensored = FALSE) {
+     C <- vapply(x, NROW, 0)
+     K <- unique(do.call("c", lapply(x, ncol)))
+     B < .... [TRUNCATED] 

> .snsc <- function(parm, x, mu = 0, rightcensored = FALSE) {
+     C <- vapply(x, NROW, 0)
+     K <- unique(do.call("c", lapply(x, ncol)))
+     B < .... [TRUNCATED] 

> .shes <- function(parm, x, mu = 0, xrc = NULL, full = FALSE, 
+     retMatrix = FALSE) {
+     C <- vapply(x, NROW, 0)
+     K <- unique(do.call("c" .... [TRUNCATED] 

> .snsr <- function(parm, x, mu = 0, rightcensored = FALSE) {
+     C <- vapply(x, NROW, 0)
+     K <- unique(do.call("c", lapply(x, ncol)))
+     B < .... [TRUNCATED] 

> .free1wayML <- function(x, link, mu = 0, start = NULL, 
+     fix = NULL, residuals = TRUE, score = TRUE, hessian = TRUE, 
+     MPL_Jeffreys = FALS .... [TRUNCATED] 

> .SW <- function(res, xt) {
+     if (length(dim(xt)) == 3) {
+         res <- matrix(res, nrow = dim(xt)[1], ncol = dim(xt)[3])
+         STAT <- Ex .... [TRUNCATED] 

> .resample <- function(res, xt, B = 10000) {
+     if (length(dim(xt)) == 2) 
+         xt <- as.table(array(xt, dim = c(dim(xt), 1)))
+     res <- m .... [TRUNCATED] 

> .p <- function(link, q, ...) link$linkinv(q = q, ...)

> .q <- function(link, p, ...) link$link(p = p, ...)

> .d <- function(link, x, ...) link$dlinkinv(x = x, 
+     ...)

> .dd <- function(link, x, ...) link$ddlinkinv(x = x, 
+     ...)

> .ddd <- function(link, x, ...) link$dddlinkinv(x = x, 
+     ...)

> .dd2d <- function(link, x, ...) link$dd2dlinkinv(x = x, 
+     ...)

> linkfun <- function(name, alias, model, parm, link, 
+     linkinv, dlinkinv, ddlinkinv, ...) {
+     ret <- list(name = name, alias = alias, model  .... [TRUNCATED] 

> logit <- function() linkfun(name = "Logit", alias = c("Wilcoxon", 
+     "Kruskal-Wallis"), model = "proportional odds", parm = "log-odds ratio", 
+ .... [TRUNCATED] 

> probit <- function() linkfun(name = "Probit", alias = "van der Waerden normal scores", 
+     model = "latent normal shift", parm = "generalised Coh ..." ... [TRUNCATED] 

> cloglog <- function() linkfun(name = "Complementary Log-log", 
+     alias = "Savage", model = "proportional hazards", parm = "log-hazard ratio", 
+ .... [TRUNCATED] 

> loglog <- function() linkfun(name = "Log-log", alias = "Lehmann", 
+     model = "Lehmann", parm = "log-reverse time hazard ratio", 
+     link = fu .... [TRUNCATED] 

> .NewtonRaphson <- function(start, objective, gradient, 
+     hessian, control = list(iter.max = 150, trace = trace, objtol = 5e-04, 
+         grad .... [TRUNCATED] 

> library("free1way.docreg")

> (x <- matrix(c(10, 5, 7, 11, 8, 9), nrow = 2))
     [,1] [,2] [,3]
[1,]   10    7    8
[2,]    5   11    9

> d <- expand.grid(y = relevel(gl(2, 1), "2"), t = gl(3, 
+     1))

> d$x <- c(x)

> m <- glm(y ~ t, data = d, weights = x, family = binomial())

> (cf <- coef(m))
(Intercept)          t2          t3 
    0.69315    -1.14513    -0.81093 

> F <- plogis

> f <- dlogis

> op <- optim(par = c(mt2 = 0, mt3 = 0, `(Intercept)` = 0), 
+     fn = .nll, gr = .nsc, x = x, method = "BFGS", hessian = TRUE)

> cbind(glm = c(cf[-1] * -1, cf[1]), free1way = op$par)
                glm free1way
t2          1.14513  1.14513
t3          0.81093  0.81093
(Intercept) 0.69315  0.69315

> logLik(m)
'log Lik.' -33.33 (df=3)

> -op$value
[1] -33.33

> fp <- function(x) {
+     p <- plogis(x)
+     p * (1 - p)^2 - p^2 * (1 - p)
+ }

> H <- .hes(op$par, x)

> solve(H$Z - crossprod(H$X, Matrix::solve(H$A, H$X)))
        [,1]    [,2]
[1,] 0.53377 0.30000
[2,] 0.30000 0.53611

> solve(op$hessian)[1:2, 1:2]
        mt2     mt3
mt2 0.53377 0.30000
mt3 0.30000 0.53611

> vcov(m)[-1, -1]
        t2      t3
t2 0.53377 0.30000
t3 0.30000 0.53611

> obj <- .free1wayML(as.table(x), link = logit())

> obj$coefficients
      B       C 
1.14513 0.81093 

> -obj$value
[1] -33.33

> obj$vcov
        B       C
B 0.53377 0.30000
C 0.30000 0.53611

> (x <- as.table(array(c(10, 5, 7, 11, 8, 9, 9, 4, 8, 
+     15, 5, 4), dim = c(2, 3, 2))))
, , A

   A  B  C
A 10  7  8
B  5 11  9

, , B

   A  B  C
A  9  8  5
B  4 15  4


> d <- expand.grid(y = relevel(gl(2, 1), "2"), t = gl(3, 
+     1), s = gl(2, 1))

> d$x <- c(x)

> m <- glm(y ~ 0 + s + t, data = d, weights = x, family = binomial())

> logLik(m)
'log Lik.' -62.523 (df=4)

> (cf <- coef(m))
      s1       s2       t2       t3 
 0.72433  0.77385 -1.30224 -0.74147 

> xl <- .table2list(x)$xlist

> op <- optim(par = c(mt2 = 0, mt3 = 0, `(Intercept 1)` = 0, 
+     `(Intercept 2)` = 0), fn = .snll, gr = .snsc, x = xl, method = "BFGS", 
+     hess .... [TRUNCATED] 

> cbind(glm = c(cf[-(1:2)] * -1, cf[1:2]), free1way = op$par)
       glm free1way
t2 1.30224  1.30225
t3 0.74147  0.74149
s1 0.72433  0.72435
s2 0.77385  0.77386

> logLik(m)
'log Lik.' -62.523 (df=4)

> -op$value
[1] -62.523

> solve(.shes(op$par, xl))
        [,1]    [,2]
[1,] 0.27095 0.16154
[2,] 0.16154 0.32006

> solve(op$hessian)[1:2, 1:2]
        mt2     mt3
mt2 0.27095 0.16154
mt3 0.16154 0.32006

> vcov(m)[-(1:2), -(1:2)]
        t2      t3
t2 0.27095 0.16153
t3 0.16153 0.32006

> obj <- .free1wayML(as.table(x), link = logit())

> obj$coefficients
      B       C 
1.30224 0.74147 

> -obj$value
[1] -62.523

> obj$vcov
        B       C
B 0.27095 0.16153
C 0.16153 0.32006

> N <- 10000

> P <- 30

> X <- matrix(rnorm(N * P), ncol = P)

> y <- X %*% runif(P) + rnorm(nrow(X))

> f <- function(par) sum((y - X %*% par)^2)

> g <- function(par) colSums(-2 * c(y - X %*% par) * 
+     X)

> h <- function(par) 2 * crossprod(X)

> start <- runif(P)

> cf <- .NewtonRaphson(start = start, objective = f, 
+     gradient = g, hessian = h)

> cf2 <- coef(m <- lm(y ~ 0 + X))

> all.equal(sum((y - fitted(m))^2), cf$objective)
[1] TRUE

> all.equal(unname(cf$par), unname(cf2))
[1] TRUE

> N <- 10

> a <- matrix(c(5, 6, 4, 3, 5, 7, 3, 4, 5, 3, 5, 6, 
+     0, 0, 0, 4, 6, 5), ncol = 3, byrow = TRUE)

> x <- as.table(array(c(a[1:3, ], a[-(1:3), ]), dim = c(3, 
+     3, 2)))

> x
, , A

  A B C
A 5 6 4
B 3 5 7
C 3 4 5

, , B

  A B C
A 3 5 6
B 0 0 0
C 4 6 5


> ret <- .free1wayML(x, logit())

> ret[c("value", "par")]
$value
[1] 65.932

$par
[1]  0.052627  0.204750 -0.489670  1.020466  0.028488


> cf <- ret$par

> cf[1:2] <- cf[1:2] + 0.5

> c(cf[1:2], cf[3], log(cf[4] - cf[3]), cf[5])
[1]  0.552627  0.704750 -0.489670  0.412200  0.028488

> .free1wayML(x, logit(), start = cf, fix = 1:2)[c("value", 
+     "par")]
$value
[1] 66.381

$par
[1]  0.55263  0.70475 -0.11885  1.40796  0.40841


> .free1wayML(x, logit(), start = cf, fix = 2)[c("value", 
+     "par")]
$value
[1] 66.299

$par
[1]  0.36095  0.70475 -0.18771  1.34153  0.33443


> .free1wayML(x, logit(), start = cf, fix = seq_along(ret$par))[c("value", 
+     "par")]
$value
[1] 67.771

$par
[1]  0.552627  0.704750 -0.489670  1.020466  0.028488


> set.seed(29)

> w <- gl(2, 15)

> (s <- .SW(r <- rank(u <- runif(length(w))), model.matrix(~0 + 
+     w)))
$Statistic
   [,1]
x2  287

$Expectation
[1] 232.5

$Covariance
       [,1]
[1,] 581.25


> ps <- .resample(r, model.matrix(~0 + w), B = 1e+05)

> ps$testStat^2
[1] 5.1101

> mean(abs(ps$permStat) > abs(ps$testStat) - .Machine$double.eps)
[1] 0.02435

> pchisq(ps$testStat^ifelse(ps$DF == 1, 2, 1), df = ps$DF, 
+     lower.tail = FALSE)
[1] 0.023787

> kruskal.test(u ~ w)

	Kruskal-Wallis rank sum test

data:  u by w
Kruskal-Wallis chi-squared = 5.11, df = 1, p-value = 0.024


> library("coin")
Loading required package: survival

> kruskal_test(u ~ w, distribution = approximate(1e+05))

	Approximative Kruskal-Wallis Test

data:  u by w (1, 2)
chi-squared = 5.11, p-value = 0.024


> wilcox_test(u ~ w, distribution = "exact")

	Exact Wilcoxon-Mann-Whitney Test

data:  u by w (1, 2)
Z = -2.26, p-value = 0.023
alternative hypothesis: true mu is not equal to 0


> free1way(u ~ w, exact = TRUE)

	2-sample Wilcoxon test against proportional odds alternatives

data:  u by w (1, 2)
W = 2.38, p-value = 0.023
alternative hypothesis: true log-odds ratio is not equal to 0


> wilcox_test(u ~ w, distribution = "exact", alternative = "less")

	Exact Wilcoxon-Mann-Whitney Test

data:  u by w (1, 2)
Z = -2.26, p-value = 0.012
alternative hypothesis: true mu is less than 0


> print(free1way(u ~ w, exact = TRUE), alternative = "greater")

	2-sample Wilcoxon test against proportional odds alternatives

data:  u by w (1, 2)
W = 2.38, p-value = 0.012
alternative hypothesis: true log-odds ratio is greater than 0


> wilcox_test(u ~ w, distribution = "exact", alternative = "greater")

	Exact Wilcoxon-Mann-Whitney Test

data:  u by w (1, 2)
Z = -2.26, p-value = 0.99
alternative hypothesis: true mu is greater than 0


> print(free1way(u ~ w, exact = TRUE), alternative = "less")

	2-sample Wilcoxon test against proportional odds alternatives

data:  u by w (1, 2)
W = 2.38, p-value = 0.99
alternative hypothesis: true log-odds ratio is less than 0


> x
, , A

  A B C
A 5 6 4
B 3 5 7
C 3 4 5

, , B

  A B C
A 3 5 6
B 0 0 0
C 4 6 5


> (ft <- free1way(x))

	Stratified 3-sample Kruskal-Wallis test against proportional odds
	alternatives

data:   by  (c("A", "B", "C")) 
	 stratified by 
Perm chi-squared = 0.145, df = 2, p-value = 0.93
alternative hypothesis: true log-odds ratio is not equal to 0


> coef(ft)
       B        C 
0.052627 0.204750 

> vcov(ft)
        B       C
B 0.34916 0.20834
C 0.20834 0.33850

> summary(ft)

Call:
free1way.table(y = x)

Coefficients:
  log-odds ratio Std. Error z value P(>|z|)
B         0.0526     0.5909  0.0891    0.93
C         0.2047     0.5818  0.3519    0.72

> library("multcomp")
Loading required package: mvtnorm
Loading required package: TH.data
Loading required package: MASS

Attaching package: ‘TH.data’

The following object is masked from ‘package:MASS’:

    geyser


> summary(glht(ft), test = univariate())

	 Simultaneous Tests for General Linear Hypotheses

Fit: free1way.table(y = x)

Linear Hypotheses:
       Estimate Std. Error z value Pr(>|z|)
B == 0   0.0526     0.5909    0.09     0.93
C == 0   0.2047     0.5818    0.35     0.72
(Univariate p values reported)


> summary(ft, test = "Wald")

	Stratified 3-sample Kruskal-Wallis test against proportional odds
	alternatives

data:   by  (c("A", "B", "C")) 
	 stratified by 
Wald chi-squared = 0.148, df = 2, p-value = 0.93
alternative hypothesis: true log-odds ratio is not equal to 0


> summary(glht(ft), test = Chisqtest())

	 General Linear Hypotheses

Linear Hypotheses:
       Estimate
B == 0   0.0526
C == 0   0.2047

Global Test:
  Chisq DF Pr(>Chisq)
1 0.148  2      0.929

> summary(ft, test = "Rao")

	Stratified 3-sample Kruskal-Wallis test against proportional odds
	alternatives

data:   by  (c("A", "B", "C")) 
	 stratified by 
Rao chi-squared = 0.148, df = 2, p-value = 0.93
alternative hypothesis: true log-odds ratio is not equal to 0


> summary(ft, test = "Permutation")

	Stratified 3-sample Kruskal-Wallis test against proportional odds
	alternatives

data:   by  (c("A", "B", "C")) 
	 stratified by 
Perm chi-squared = 0.145, df = 2, p-value = 0.93
alternative hypothesis: true log-odds ratio is not equal to 0


> summary(ft, test = "LRT")

	Stratified 3-sample Kruskal-Wallis test against proportional odds
	alternatives

data:   by  (c("A", "B", "C")) 
	 stratified by 
logLR chi-squared = 0.148, df = 2, p-value = 0.93
alternative hypothesis: true log-odds ratio is not equal to 0


> confint(glht(ft), calpha = univariate_calpha())

	 Simultaneous Confidence Intervals

Fit: free1way.table(y = x)

Quantile = 1.96
95% confidence level
 

Linear Hypotheses:
       Estimate lwr     upr    
B == 0  0.0526  -1.1055  1.2108
C == 0  0.2047  -0.9356  1.3451


> confint(ft, test = "Wald")
     2.5 % 97.5 %
B -1.10552 1.2108
C -0.93558 1.3451

> confint(ft, test = "Rao")
     2.5 % 97.5 %
B -1.09075 1.1954
C -0.92072 1.3301

> confint(ft, test = "LRT")
    2.5 % 97.5 %
B -1.1113 1.2205
C -0.9387 1.3573

> N <- 25

> w <- gl(2, N)

> y <- rlogis(length(w), location = c(0, 1)[w])

> ft <- free1way(y ~ w)

> summary(ft)

Call:
free1way.formula(formula = y ~ w)

Coefficients:
   log-odds ratio Std. Error z value P(>|z|)
w2          1.290      0.524   2.462    0.01

> wilcox.test(y ~ w, alternative = "greater", correct = FALSE, 
+     exact = FALSE)$p.value
[1] 0.99215

> pvalue(wilcox_test(y ~ w, alternative = "greater"))
[1] 0.99215

> summary(ft, test = "Permutation", alternative = "less")$p.value
 Perm Z 
0.99215 

> wilcox.test(y ~ w, alternative = "less", correct = FALSE, 
+     exact = FALSE)$p.value
[1] 0.0078535

> pvalue(wilcox_test(y ~ w, alternative = "less"))
[1] 0.0078535

> summary(ft, test = "Permutation", alternative = "greater")$p.value
   Perm Z 
0.0078535 

> wilcox.test(y ~ w, correct = FALSE, exact = FALSE)$p.value
[1] 0.015707

> kruskal.test(y ~ w)$p.value
[1] 0.015707

> pvalue(wilcox_test(y ~ w))
[1] 0.015707

> summary(ft, test = "Permutation")$p.value
  Perm Z 
0.015707 

> summary(ft, test = "Wald", alternative = "less")

	2-sample Wilcoxon test against proportional odds alternatives

data:  y by w (1, 2)
Wald Z = 2.46, p-value = 0.99
alternative hypothesis: true log-odds ratio is less than 0


> summary(ft, test = "Wald", alternative = "greater")

	2-sample Wilcoxon test against proportional odds alternatives

data:  y by w (1, 2)
Wald Z = 2.46, p-value = 0.0069
alternative hypothesis: true log-odds ratio is greater than 0


> summary(ft, test = "Wald")

	2-sample Wilcoxon test against proportional odds alternatives

data:  y by w (1, 2)
Wald chi-squared = 6.06, df = 1, p-value = 0.014
alternative hypothesis: true log-odds ratio is not equal to 0


> summary(ft, test = "Rao", alternative = "less")

	2-sample Wilcoxon test against proportional odds alternatives

data:  y by w (1, 2)
Rao Z = 2.5, p-value = 0.99
alternative hypothesis: true log-odds ratio is less than 0


> summary(ft, test = "Rao", alternative = "greater")

	2-sample Wilcoxon test against proportional odds alternatives

data:  y by w (1, 2)
Rao Z = 2.5, p-value = 0.0063
alternative hypothesis: true log-odds ratio is greater than 0


> summary(ft, test = "Rao")

	2-sample Wilcoxon test against proportional odds alternatives

data:  y by w (1, 2)
Rao chi-squared = 6.24, df = 1, p-value = 0.013
alternative hypothesis: true log-odds ratio is not equal to 0


> summary(ft, test = "LRT")

	2-sample Wilcoxon test against proportional odds alternatives

data:  y by w (1, 2)
logLR chi-squared = 6.33, df = 1, p-value = 0.012
alternative hypothesis: true log-odds ratio is not equal to 0


> confint(ft, test = "Permutation")
     2.5 % 97.5 %
w2 0.23829 2.3423

> confint(ft, test = "LRT")
     2.5 % 97.5 %
w2 0.28151 2.3446

> confint(ft, test = "Wald")
     2.5 % 97.5 %
w2 0.26308 2.3168

> confint(ft, test = "Rao")
     2.5 % 97.5 %
w2 0.27214 2.3084

> confint(ft, test = "Rao", what = "AUC")
     2.5 %  97.5 %
w2 0.54524 0.82742

> library("rms")

  When sourcing ‘free1way.R’:
Error: there is no package called ‘rms’
Execution halted

  ‘free1way.Rnw’ using ‘UTF-8’... failed