# tests/check_transf_thresholds.R In mvord: Multivariate Ordinal Regression Models

```library(mvord)
library(Matrix)
rho <- list()
rho\$ndim <- 5
rho\$threshold.values <- list(c(NA),
c(NA,NA,NA),
c(NA,NA,NA),
c(NA,NA,NA,NA,NA),
c(NA,NA))
rho\$ntheta <- sapply(rho\$threshold.values, length)
rho\$ind.thresholds <- list(1,2:4,5:7,8:12,13:14)

theta <- c(1,0,1,2,1,1,1,-5:-1, -10,-3)

mvord:::check(all.equal(mvord:::transf_thresholds_flexible(theta, rho), list(1,
c(0, exp(1), exp(1) + exp(2)),
c(1,1 + exp(1), 1+ exp(1) + exp(1)),
c(-5, -5 + exp(-4), -5 + exp(-4) + exp(-3), -5 + exp(-4) + exp(-3) +exp(-2), -5 + exp(-4) + exp(-3) + exp(-2) + exp(-1)),
c(-10, -10 + exp(-3)))))

rho\$threshold.constraints <- c(1,2,2,3,4)
rho\$ntheta <- sapply(rho\$threshold.values, length)
rho\$ind.thresholds <- list(1,2:4,2:4,5:9,10:11)

theta <- c(1,0,1,2,-5:-1, -10,-3)

mvord:::check(all.equal(mvord:::transf_thresholds_flexible(theta, rho), list(1,
c(0, exp(1), exp(1) + exp(2)),
c(0,exp(1), exp(1) + exp(2)),
c(-5, -5 + exp(-4), -5 + exp(-4) + exp(-3), -5 + exp(-4) + exp(-3) +exp(-2), -5 + exp(-4) + exp(-3) + exp(-2) + exp(-1)),
c(-10, -10 + exp(-3)))))
##################################################################################
rho\$threshold.values <- list(1,
c(2,NA,NA),
c(3,NA,NA),
c(4,NA,NA,NA,NA),
c(-9,NA))
rho\$ndim <- 5
rho\$ntheta <- sapply(rho\$threshold.values, length)
rho\$npar.theta <- sapply(1:rho\$ndim, function(j) sum(is.na(rho\$threshold.values[[j]])))
rho\$ind.thresholds <- list(c(),1:2,3:4,5:8,9)
rho\$threshold.values.fixed <- lapply(1:rho\$ndim, function(j){ind <- !is.na(rho\$threshold.values[[j]])
rho\$threshold.values[[j]][ind]})

rho\$ncat.first.ind <- cumsum(c(1,rho\$ntheta))[- (rho\$ndim + 1)]
rho\$inds.cat <- lapply(seq_len(rho\$ndim), function(j)
seq_len(rho\$ntheta[j]) +  rho\$ncat.first.ind[j] - 1)

rho\$constraints <- as.data.frame(c(1,1,1,1,1,1,1,1,1,1,1,1,1,1), ncol = 1)
rho\$contr_theta <- do.call("rbind", rep(list(diag(14)), 1))
## make function which compute correction factor for thresholds:
rho\$mat_center_scale <- bdiag(rho\$constraints) %*% 1

betatilde <- 0

rho\$thold_correction <-lapply(seq_len(rho\$ndim), mvord:::build_correction_thold_fun0, rho = rho)

theta <- c(1,0,1,2,1,1,-3,-10, -5)

mvord:::check(all.equal(mvord:::transf_thresholds_fix1_first(theta, rho, betatilde), list(1,
c(2, 2 + exp(1), 2 + exp(1) + exp(0)),
c(3, 3 + exp(1), 3 + exp(1)+ exp(2)),
c(4, 4+ exp(1), 4+ exp(1)+ exp(1), 4+ exp(1)+ exp(1)+ exp(-3),4+ exp(1)+ exp(1)+ exp(-3)+ exp(-10)),
c(-9, -9 + exp(-5)))))

rho\$threshold.values <- list(1,
c(2,NA,NA),
c(2,NA,NA),
c(4,NA,NA,NA,NA),
c(-9,NA))
rho\$ndim <- 5
rho\$ntheta <- sapply(rho\$threshold.values, length)
rho\$npar.theta <- sapply(1:rho\$ndim, function(j) sum(is.na(rho\$threshold.values[[j]])))
rho\$ind.thresholds <- list(c(),1:2,1:2,3:6,7)
rho\$threshold.values.fixed <- lapply(1:rho\$ndim, function(j){ind <- !is.na(rho\$threshold.values[[j]])
rho\$threshold.values[[j]][ind]})

rho\$ncat.first.ind <- cumsum(c(1,rho\$ntheta))[- (rho\$ndim + 1)]
rho\$inds.cat <- lapply(seq_len(rho\$ndim), function(j)
seq_len(rho\$ntheta[j]) +  rho\$ncat.first.ind[j] - 1)

rho\$constraints <- as.data.frame(c(1,1,1,1,1,1,1,1,1,1,1,1,1,1), ncol = 1)
rho\$contr_theta <- do.call("rbind", rep(list(diag(14)), 1))
## make function which compute correction factor for thresholds:
rho\$mat_center_scale <- bdiag(rho\$constraints) %*% 1

betatilde <- 0

theta <- c(1,0,1,1,-3,-10, -5)

mvord:::check(all.equal(mvord:::transf_thresholds_fix1_first(theta, rho, betatilde), list(1,
c(2, 2 + exp(1), 2 + exp(1) + exp(0)),
c(2, 2 + exp(1), 2 + exp(1) + exp(0)),
c(4, 4+ exp(1), 4+ exp(1)+ exp(1), 4+ exp(1)+ exp(1)+ exp(-3),4+ exp(1)+ exp(1)+ exp(-3)+ exp(-10)),
c(-9, -9 + exp(-5)))))

##################################################################################
rho\$threshold.values <- list(1,
c(2,4,NA),
c(3,4,NA),
c(4,10,NA,NA,NA),
c(-9,-6))
rho\$ndim <- 5
rho\$ntheta <- sapply(rho\$threshold.values, length)
rho\$npar.theta <- sapply(1:rho\$ndim, function(j) sum(is.na(rho\$threshold.values[[j]])))
rho\$ind.thresholds <- list(c(),1,2,3:5,c())
rho\$threshold.values.fixed <- lapply(1:rho\$ndim, function(j){ind <- !is.na(rho\$threshold.values[[j]])
rho\$threshold.values[[j]][ind]})

rho\$first.ind.theta <- c(NA,1,2,3,NA)

rho\$ncat.first.ind <- cumsum(c(1,rho\$ntheta))[- (rho\$ndim + 1)]
rho\$inds.cat <- lapply(seq_len(rho\$ndim), function(j)
seq_len(rho\$ntheta[j]) +  rho\$ncat.first.ind[j] - 1)

rho\$constraints <- as.data.frame(c(1,1,1,1,1,1,1,1,1,1,1,1,1,1), ncol = 1)
rho\$contr_theta <- do.call("rbind", rep(list(diag(14)), 1))
## make function which compute correction factor for thresholds:
rho\$mat_center_scale <- bdiag(rho\$constraints) %*% 1

betatilde <- 0

gamma <- c(1,0,1,2,-10)

mvord:::check(all.equal(mvord:::transf_thresholds_fix2_first(gamma, rho, betatilde), list(1,
c(2, 4, 4 + exp(1)),
c(3, 4, 4 + exp(0)),
c(4,10, 10 + exp(1), 10 + exp(1) + exp(2),10 + exp(1) + exp(2) + exp(-10)),
c(-9, -6))))

rho\$threshold.values <- list(1,
c(2,4,NA),
c(2,4,NA),
c(4,10,NA,NA,NA),
c(-9,-6))
rho\$ndim <- 5
rho\$ntheta <- sapply(rho\$threshold.values, length)
rho\$npar.theta <- sapply(1:rho\$ndim, function(j) sum(is.na(rho\$threshold.values[[j]])))
rho\$ind.thresholds <- list(c(),1,1,2:4,c())
rho\$threshold.values.fixed <- lapply(1:rho\$ndim, function(j){ind <- !is.na(rho\$threshold.values[[j]])
rho\$threshold.values[[j]][ind]})

rho\$first.ind.theta <- c(NA,1,1,2,NA)

rho\$ncat.first.ind <- cumsum(c(1,rho\$ntheta))[- (rho\$ndim + 1)]
rho\$inds.cat <- lapply(seq_len(rho\$ndim), function(j)
seq_len(rho\$ntheta[j]) +  rho\$ncat.first.ind[j] - 1)

rho\$constraints <- as.data.frame(c(1,1,1,1,1,1,1,1,1,1,1,1,1,1), ncol = 1)
rho\$contr_theta <- do.call("rbind", rep(list(diag(14)), 1))
## make function which compute correction factor for thresholds:
rho\$mat_center_scale <- bdiag(rho\$constraints) %*% 1

betatilde <- 0

gamma <- c(1,1,2,-10)

mvord:::check(all.equal(mvord:::transf_thresholds_fix2_first(gamma, rho, betatilde), list(1,
c(2, 4, 4 + exp(1)),
c(2, 4, 4 + exp(1)),
c(4,10, 10 + exp(1), 10 + exp(1) + exp(2),10 + exp(1) + exp(2) + exp(-10)),
c(-9, -6))))
##################################################################################
rho\$threshold.values <- list(1,
c(2,NA,4),
c(2,NA,4),
c(4,NA,NA,NA,10),
c(-9,-6))
rho\$ndim <- 5
rho\$ntheta <- sapply(rho\$threshold.values, length)
rho\$npar.theta <- sapply(1:rho\$ndim, function(j) sum(is.na(rho\$threshold.values[[j]])))
rho\$ind.thresholds <- list(c(),1,1,2:4,c())
rho\$threshold.values.fixed <- lapply(1:rho\$ndim, function(j){ind <- !is.na(rho\$threshold.values[[j]])
rho\$threshold.values[[j]][ind]})

gamma <- c(1,1,2,-10)

## (theta_j - theta_(j-1)/(1 - theta_j) = exp(gamma_j)/(1 + exp(gamma_j))
th41 <- exp(1)/(1 + exp(1))
th42 <- (exp(2) + th41)/(1 + exp(2))
th43 <- (exp(-10) + th42)/(1 + exp(-10))

mvord:::check(all.equal(mvord:::transf_thresholds_fix2_firstlast(gamma, rho, betatilde),
list(1,
c(2, c(exp(1)/(1 + exp(1))) *(4 - 2) + 2, 4),
c(2, c(exp(1)/(1 + exp(1))) *(4 - 2) + 2, 4),
c(4, c(th41, th42, th43) *(10 - 4) + 4, 10),
c(-9, -6))))

rho\$threshold.values <- list(1,
c(2,NA,4),
c(3,NA,4),
c(4,NA,NA,NA,10),
c(-9,-6))
rho\$ndim <- 5
rho\$ntheta <- sapply(rho\$threshold.values, length)
rho\$npar.theta <- sapply(1:rho\$ndim, function(j) sum(is.na(rho\$threshold.values[[j]])))
rho\$ind.thresholds <- list(c(),1,2,3:5,c())
rho\$threshold.values.fixed <- lapply(1:rho\$ndim, function(j){ind <- !is.na(rho\$threshold.values[[j]])
rho\$threshold.values[[j]][ind]})

gamma <- c(1,0,1,2,-10)

## (theta_j - theta_(j-1)/(1 - theta_j) = exp(gamma_j)/(1 + exp(gamma_j))
th41 <- exp(1)/(1 + exp(1))
th42 <- (exp(2) + th41)/(1 + exp(2))
th43 <- (exp(-10) + th42)/(1 + exp(-10))

mvord:::check(all.equal(mvord:::transf_thresholds_fix2_firstlast(gamma, rho, betatilde),
list(1,
c(2, c(exp(1)/(1 + exp(1))) *(4 - 2) + 2, 4),
c(3, c(exp(0)/(1 + exp(0))) *(4 - 3) + 3, 4),
c(4, c(th41, th42, th43) *(10 - 4) + 4, 10),
c(-9, -6))))
```

## Try the mvord package in your browser

Any scripts or data that you put into this service are public.

mvord documentation built on Oct. 31, 2018, 9:03 a.m.