CRAN Package Check Results for Package mvord

Last updated on 2020-05-29 10:47:28 CEST.

Flavor Version Tinstall Tcheck Ttotal Status Flags
r-devel-linux-x86_64-debian-clang 1.0.0 14.20 171.56 185.76 ERROR
r-devel-linux-x86_64-debian-gcc 1.0.0 11.16 130.68 141.84 ERROR
r-devel-linux-x86_64-fedora-clang 1.0.1 235.66 OK
r-devel-linux-x86_64-fedora-gcc 1.0.1 230.69 OK
r-devel-windows-ix86+x86_64 1.0.0 32.00 315.00 347.00 ERROR
r-patched-linux-x86_64 1.0.0 13.98 176.40 190.38 OK
r-patched-solaris-x86 1.0.1 306.80 OK
r-release-linux-x86_64 1.0.0 11.77 175.53 187.30 OK
r-release-osx-x86_64 1.0.0 OK
r-release-windows-ix86+x86_64 1.0.0 23.00 420.00 443.00 OK
r-oldrel-osx-x86_64 1.0.0 NOTE
r-oldrel-windows-ix86+x86_64 1.0.0 30.00 253.00 283.00 OK

Check Details

Version: 1.0.0
Check: tests
Result: ERROR
     Running 'check_methods.R' [6s/7s]
     Running 'check_mvord_data.R' [2s/3s]
     Running 'check_mvord_errors.R' [3s/3s]
     Running 'check_set_threshold_type.R' [2s/3s]
     Running 'check_toy_example.R' [32s/39s]
     Running 'check_transf_sigmas.R' [2s/3s]
     Running 'check_transf_thresholds.R' [2s/3s]
    Running the tests in 'tests/check_toy_example.R' failed.
    Complete output:
     > library(mvord)
     Loading required package: minqa
     Loading required package: BB
     Loading required package: ucminf
     Loading required package: dfoptim
     > library(MASS)
     >
     >
     > #data(data_toy_example)
     > tolerance <- 1e-6
     >
     > mult.obs <- 2
     > sigma <- matrix(c(1,0.8,0.8,1), ncol = 2)
     > betas <- list(c(0.8,-0.5),
     + c(0.8,-0.5))
     > thresholds <- list(c(-1,1),c(-1,1))
     > nobs <- 100
     > suppressWarnings(RNGversion("3.5.0"))
     > set.seed(2017)
     > errors <- mvrnorm(n = nobs, mu = rep(0, mult.obs), Sigma = sigma)
     >
     > X1 <- rnorm(nobs, 0, 1)
     > X2 <- rnorm(nobs, 0, 1)
     >
     > pred <- cbind(X1, X2)
     >
     > y <- sapply(1:mult.obs, function(j) pred %*% betas[[j]] + errors[, j], simplify = "array")
     > y.ord <- sapply(1:mult.obs, function(j) cut(y[, , j], c(min(y[, , j]) - 1, c(thresholds[[j]]), max(y[, , j]) + 1),
     + labels = FALSE), simplify = "array")
     > predictors.fixed <- lapply(1:mult.obs, function(j) pred)
     > y <- as.data.frame(y.ord)
     >
     > for(i in 1:mult.obs){
     + y[, i] <- factor(y[, i], levels = sort(unique(y[, i])),
     + ordered = TRUE)
     + }
     >
     >
     >
     >
     > data_toy_example <- cbind.data.frame(y, predictors.fixed[[1]])
     > colnames(data_toy_example) <- c("Y1","Y2", "X1", "X2")
     > w <- c(rep(1/20, 20), rep(1/30, 30), rep(1/20, 20), rep(1/30, 30))
     >
     >
     >
     > # convert data_toy_example into long format
     > df <- cbind.data.frame("i" = rep(1:100,2),
     + "j" = rep(1:2,each = 100),
     + "Y" = c(data_toy_example$Y1,data_toy_example$Y2),
     + "X1" = rep(data_toy_example$X1,2),
     + "X2" = rep(data_toy_example$X2,2),
     + "f1" = factor(sample(rep(data_toy_example$Y2,2)), ordered =F),
     + "f2" = factor(rep(data_toy_example$Y1,2), ordered=F),
     + w=rep(w,2))
     >
     >
     >
     > # library(ROI)
     > # ROI_solver <- function(starting.values, objFun, control){
     > # n <- length(starting.values)
     > # op <- OP(objective = F_objective(objFun, n = n),
     > # bounds = V_bound(li = seq_len(n), lb = rep.int(-Inf, n)))
     > # optRes <- ROI_solve(op, solver = "nlminb",
     > # control = c(list(start = starting.values),
     > # control))
     > # list(optpar = optRes$solution,
     > # objective = optRes$objval) # a vector of length equal to number of parameters to optimize
     > # }
     > #
     > #
     > #
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + link = mvprobit(),
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = c(1,1),
     + contrasts = list(f1 = function(x)
     + contr.treatment(nlevels(df$f1), base = 1),
     + f2 = "contr.sum"),
     + control= mvord.control(solver="BFGS",se=T))
     Warning message:
     In mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, link = mvprobit(), :
     Variables f1 and f2 are absent, the contrasts will be ignored.
     > options(digits = 22)
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cor_general(~1),
     link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), contrasts = list(f1 = function(x) contr.treatment(nlevels(df$f1),
     base = 1), f2 = "contr.sum"), control = mvord.control(solver = "BFGS",
     se = T))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -134.91 280.34 294.06 32
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -0.96257386663519672876 0.16613952738605439197 -5.7937700000000003087
     1 2|3 1.03347036873223707687 0.15004482537617938598 6.8877399999999999736
     2 1|2 -0.96257386663519672876 0.16613952738605439197 -5.7937700000000003087
     2 2|3 1.03347036873223707687 0.15004482537617938598 6.8877399999999999736
     Pr(>|z|)
     1 1|2 6.8825e-09 ***
     1 2|3 5.6684e-12 ***
     2 1|2 6.8825e-09 ***
     2 2|3 5.6684e-12 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     X1 1 0.63801035062404309883 0.13576747301357891540 4.6992900000000004113
     X2 1 -0.42672524816263474046 0.13643665802446264257 -3.1276399999999999757
     Pr(>|z|)
     X1 1 2.6107e-06 ***
     X2 1 0.0017621 **
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.854266728221226845363 0.062440889990360744222 13.681210000000000093
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.96257386663519672876, 1.03347036873223707687, -0.96257386663519672876, 1.03347036873223707687), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.63801035062404309883, -0.42672524816263474046), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.85426672822122684536), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.16613952738605441972, 0.15004482537617935822, 0.16613952738605441972, 0.15004482537617935822), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.13576747301357894315, 0.13643665802446261481), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.062440889990360744222), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -134.90867383086322207, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 280.3436634512001433, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 294.05508548271637892, tolerance = tolerance))
     >
     >
     > res2 <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = list(matrix(rep(1,4), ncol = 1), matrix(rep(1,4), ncol = 1)))
     >
     > mvord:::check(all.equal(res$beta, res2$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res2$sebeta, tolerance = tolerance))
     > ########################################################################
     > ## No coefficients
     > res <- mvord:::mvord(formula = MMO(Y) ~ -1,
     + data = df,
     + link = mvprobit(),
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + control= mvord.control(solver="BFGS",se=T))
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ -1, data = df, error.structure = cor_general(~1),
     link = mvprobit(), threshold.constraints = c(1, 1), control = mvord.control(solver = "BFGS",
     se = T))
    
     Formula: MMO(Y) ~ -1
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -153.66 313.51 321.57 26
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -0.75479706538110091785 0.13105603927553963195 -5.7593500000000004135
     1 2|3 0.86086304364935783973 0.13360642334434200129 6.4432799999999996743
     2 1|2 -0.75479706538110091785 0.13105603927553963195 -5.7593500000000004135
     2 2|3 0.86086304364935783973 0.13360642334434200129 6.4432799999999996743
     Pr(>|z|)
     1 1|2 8.4440e-09 ***
     1 2|3 1.1692e-10 ***
     2 1|2 8.4440e-09 ***
     2 2|3 1.1692e-10 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.905795171446422409112 0.038855435949034601573 23.311930000000000263
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.75479706538110091785, 0.86086304364935783973, -0.75479706538110091785, 0.86086304364935783973), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(numeric(0)), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.90579517144642240911), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.13105603927553965971, 0.13360642334434202905, 0.13105603927553965971, 0.13360642334434202905), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(numeric(0)), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.038855435949034601573), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -153.66397119528727444, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 313.51350940088383368, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 321.57073678022845797, tolerance = tolerance))
     >
     > #polychor
     > res <- mvord:::mvord(formula = MMO(Y) ~ 1,
     + data = df,
     + link = mvprobit(),
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + control= mvord.control(solver="BFGS",se=T))
     Note: First threshold for each response is fixed to 0 in order to ensure identifiability!
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 1, data = df, error.structure = cor_general(~1),
     link = mvprobit(), threshold.constraints = c(1, 1), control = mvord.control(solver = "BFGS",
     se = T))
    
     Formula: MMO(Y) ~ 1
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit fix1first 100 2 -153.56 315.46 326.32 52
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 0.00000000000000000000 0.00000000000000000000 NA
     1 2|3 1.61576075469784496974 0.15607113684416781818 10.3527199999999997
     2 1|2 0.00000000000000000000 0.00000000000000000000 NA
     2 2|3 1.61576075469784496974 0.15607113684416781818 10.3527199999999997
     Pr(>|z|)
     1 1|2 NA
     1 2|3 < 2.22e-16 ***
     2 1|2 NA
     2 2|3 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error
     (Intercept) 1 0.73786225358709867095 0.13500906208906057748
     (Intercept) 2 0.77255634837212872057 0.14063945843173736305
     z value Pr(>|z|)
     (Intercept) 1 5.4652799999999999159 4.6218e-08 ***
     (Intercept) 2 5.4931700000000001083 3.9478e-08 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.906373932662642656233 0.038916052675185490439 23.290489999999998361
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(0.0000000000000000000, 1.6157607546978449697, 0.0000000000000000000, 1.6157607546978449697), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.73786225358709867095, 0.77255634837212872057), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.90637393266264265623), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.00000000000000000000, 0.15607113684416776267, 0.00000000000000000000, 0.15607113684416776267), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.13500906208906049422, 0.14063945843173727979), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.038916052675185497378), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -153.5640565887688922, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 315.46144651087109878, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 326.31632228582151356, tolerance = tolerance))
     > #######################################################################
     > ## cor_general(~factor)
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + error.structure = cor_general(~f2),
     + threshold.constraints = c(1,1),
     + coef.constraints = c(1,1),
     + contrasts = list(f2 = "contr.sum"))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cor_general(~f2),
     link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), contrasts = list(f2 = "contr.sum"), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -134.17 283.39 303 38
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -0.90552506790519915469 0.18265647646915181279 -4.9575300000000002143
     1 2|3 1.00420547521745429087 0.17467164336594764862 5.7491000000000003212
     2 1|2 -0.90552506790519915469 0.18265647646915181279 -4.9575300000000002143
     2 2|3 1.00420547521745429087 0.17467164336594764862 5.7491000000000003212
     Pr(>|z|)
     1 1|2 7.1395e-07 ***
     1 2|3 8.9718e-09 ***
     2 1|2 7.1395e-07 ***
     2 2|3 8.9718e-09 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     X1 1 0.64793126338493944871 0.14153585035881874332 4.5778600000000002623
     X2 1 -0.42893128334048874484 0.13989943737309903926 -3.0659999999999998366
     Pr(>|z|)
     X1 1 4.6976e-06 ***
     X2 1 0.0021695 **
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error
     corr f21 1 2 0.73667859723415496376 0.17823839430849072740
     corr f22 1 2 0.91829826305234418804 0.06885865619995423792
     corr f23 1 2 0.83710611971260706632 0.13103608299872890330
     z value Pr(>|z|)
     corr f21 1 2 4.1331100000000002836 3.5789e-05 ***
     corr f22 1 2 13.3359900000000006770 < 2.22e-16 ***
     corr f23 1 2 6.3883599999999995944 1.6767e-10 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate,
     + c(-0.90552506790519915469, 1.00420547521745429087, -0.90552506790519915469, 1.00420547521745429087), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate,
     + c(0.64793126338493944871, -0.42893128334048874484), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate,
     + c(0.73667859723415496376, 0.91829826305234418804, 0.83710611971260706632), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`,
     + c(0.18265647646915181279, 0.17467164336594759311, 0.18265647646915181279, 0.17467164336594759311), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`,
     + c(0.14153585035881874332, 0.13989943737309898375), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`,
     + c(0.17823839430849069965, 0.06885865619995423792, 0.13103608299872890330), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -134.17046176709035876, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 283.39468697504094052, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 303.00349482656417877, tolerance = tolerance))
     >
     > ##################################################################################
     > res <- mvord:::mvord(formula = MMO(Y, i, j) ~ 0 + X1 + X2,
     + data = df,
     + link = mvlogit(df = 8L),
     + error.structure = cor_general(~1))
     > #threshold.constraints = c(1,1),
     > #coef.constraints = c(1,1))
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y, i, j) ~ 0 + X1 + X2, data = df,
     error.structure = cor_general(~1), link = mvlogit(df = 8L))
    
     Formula: MMO(Y, i, j) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvlogit flexible 100 2 -132.87 285.51 311.28 541
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -1.61992816122470717488 0.33738752096507018674 -4.8013899999999996027
     1 2|3 1.83247502787750837783 0.36742758401929348500 4.9873099999999999099
     2 1|2 -1.61590236494067829298 0.33963203103205336086 -4.7577999999999995850
     2 2|3 1.74613466597926825230 0.30112279474893460796 5.7987500000000000711
     Pr(>|z|)
     1 1|2 1.5757e-06 ***
     1 2|3 6.1226e-07 ***
     2 1|2 1.9571e-06 ***
     2 2|3 6.6813e-09 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     X1 1 1.28924245583443863872 0.30395056815485199664 4.2416200000000001680
     X1 2 0.88420412270574577640 0.29740838248021878032 2.9730300000000000615
     X2 1 -0.84290979913631414178 0.30080366448985135230 -2.8021899999999999586
     X2 2 -0.71568734501553299410 0.25697668783955229799 -2.7850299999999998946
     Pr(>|z|)
     X1 1 2.2191e-05 ***
     X1 2 0.0029488 **
     X2 1 0.0050757 **
     X2 2 0.0053523 **
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.883893242218400709831 0.071174772540211109217 12.41863000000000028
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > # mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-1.6170817306633420, 1.7855897338762188, -1.6170817306633420, 1.7855897338762188), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$coefficients$Estimate, c(1.07242200717115987, 1.07242200717115987, -0.76715925377701732, -0.76715925377701732), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.85317690560688531), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.28997405649500174, 0.27427389826231802, 0.28997405649500174, 0.27427389826231802), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.24111402270822993, 0.24111402270822993, 0.24156664773225886, 0.24156664773225886), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.063316529381183581), tolerance = tolerance))
     > # mvord:::check(all.equal(logLik(res), -135.41665313840898, tolerance = tolerance))
     > # mvord:::check(all.equal(AIC(res), 281.35962206629165, tolerance = tolerance))
     > # mvord:::check(all.equal(BIC(res), 295.07104409780789, tolerance = tolerance))
     >
     > ##################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cov_general(~1),
     + threshold.constraints = c(1,1),
     + threshold.values = list(c(-1,NA),
     + c(-1,NA)),
     + coef.constraints = c(1,1))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cov_general(~1),
     link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), threshold.values = list(c(-1, NA), c(-1, NA)), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit fix1first 100 2 -134.64 282.04 298.67 37
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -1.00000000000000000000 0.00000000000000000000 NA
     1 2|3 1.08266798008347264748 0.24926043418153945308 4.3435199999999998255
     2 1|2 -1.00000000000000000000 0.00000000000000000000 NA
     2 2|3 1.08266798008347264748 0.24926043418153945308 4.3435199999999998255
     Pr(>|z|)
     1 1|2 NA
     1 2|3 1.4022e-05 ***
     2 1|2 NA
     2 2|3 1.4022e-05 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     X1 1 0.68834946142116071321 0.15128036741855890868 4.5501599999999999824
     X2 1 -0.45406517655044392745 0.15704134440910905157 -2.8913700000000002177
     Pr(>|z|)
     X1 1 5.3606e-06 ***
     X2 1 0.0038356 **
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.856153637803035327813 0.061822221469868911259 13.8486399999999996169
     sigma 1 1.001934583868995254363 0.185187960683942831608 5.4103700000000003456
     sigma 2 1.089828199315897139243 0.197593289723884180109 5.5155099999999999127
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     sigma 1 6.2896e-08 ***
     sigma 2 3.4777e-08 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-1.000000000000000000000, 1.082667980083374503764, -1.000000000000000000000, 1.082667980083374503764), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.6883494614209317852271, -0.4540651765505133163892), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.8561536378031743277361, 1.0019345838689188710191, 1.0898281993157663549709), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.0000000000000000000000, 0.2492604341815170820862, 0.0000000000000000000000, 0.2492604341815170820862), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1512803674185349001036, 0.1570413444091093568833), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.06182222146981963123435, 0.18518796068392143205905, 0.19759328972385112321852), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -134.6391112878561102661, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 282.0441800225207202857, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 298.6729258905298252103, tolerance = tolerance))
     > ########################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_equi(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = cbind(c(1,1),c(1,2)))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cor_equi(~1),
     link = mvprobit(), coef.constraints = cbind(c(1, 1), c(1,
     2)), threshold.constraints = c(1, 1), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -134.84 282.45 299.08 35
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -0.96275665235846763501 0.16703865743956475276 -5.7636799999999999145
     1 2|3 1.03388029135962500327 0.15203091606549479220 6.8004600000000001714
     2 1|2 -0.96275665235846763501 0.16703865743956475276 -5.7636799999999999145
     2 2|3 1.03388029135962478122 0.15203091606549479220 6.8004600000000001714
     Pr(>|z|)
     1 1|2 8.2301e-09 ***
     1 2|3 1.0429e-11 ***
     2 1|2 8.2301e-09 ***
     2 2|3 1.0429e-11 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     X1 1 0.63820467410958414689 0.13670492809263015688 4.6684799999999997411
     X2 1 -0.44676561119553748203 0.15900879662664607617 -2.8096899999999997988
     X2 2 -0.40750155519132619242 0.13850979189699536009 -2.9420399999999999885
     Pr(>|z|)
     X1 1 3.0343e-06 ***
     X2 1 0.0049589 **
     X2 2 0.0032606 **
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     (Intercept) 1.27290057125462996446 0.23221611847930476169 5.4815300000000002356
     Pr(>|z|)
     (Intercept) 4.2165e-08 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.9627566523584676350112, 1.0338802913596250032668, -0.9627566523584676350112, 1.0338802913596247812222), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.6382046741095841468905, -0.4467656111955374820255, -0.4075015551913261924177), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$error.structure$Estimate, c(1.272900571254629964457), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.1670386574395647250046, 0.1520309160654948199554, 0.1670386574395647250046, 0.1520309160654948199554), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1367049280926301846328, 0.1590087966266460206555, 0.1385097918969951935608), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.2322161184793045674013), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -134.8432321319771745038, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 282.4524217107628487611, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 299.0811675787719536856, tolerance = tolerance))
     >
     > res2 <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_equi(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = list(X2 = cbind(c(1,1,0,0), c(0,0,1,1)),
     + X1 = matrix(rep(1,4), ncol = 1)))
     >
     > mvord:::check(all.equal(res$beta, res2$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res2$sebeta, tolerance = tolerance))
     >
     > ########################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_ar1(~ 1 + X1),
     + threshold.constraints = c(1,1),
     + coef.constraints = c(1,1))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cor_ar1(~1 +
     X1), link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -133.94 280.64 297.27 35
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -0.95722951152617552495 0.16558288168119875139 -5.7809699999999999420
     1 2|3 1.03746793968331618707 0.15306952616698107916 6.7777599999999997848
     2 1|2 -0.95722951152617552495 0.16558288168119875139 -5.7809699999999999420
     2 2|3 1.03746793968331618707 0.15306952616698107916 6.7777599999999997848
     Pr(>|z|)
     1 1|2 7.4272e-09 ***
     1 2|3 1.2206e-11 ***
     2 1|2 7.4272e-09 ***
     2 2|3 1.2206e-11 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     X1 1 0.65304208674282193670 0.13584016202004611795 4.8074300000000000921
     X2 1 -0.42273136793544752177 0.13820217716271032682 -3.0587900000000001199
     Pr(>|z|)
     X1 1 1.5288e-06 ***
     X2 1 0.0022223 **
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error
     (Intercept) 1.29891554197132030879 0.24517466007553267993
     X1 0.29350892362932679003 0.29416988646865149803
     z value Pr(>|z|)
     (Intercept) 5.29792000000000040671 1.1713e-07 ***
     X1 0.99775000000000002576 0.3184
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.9572295115261755249492, 1.0374679396833161870717, -0.9572295115261755249492, 1.0374679396833161870717), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.6530420867428219366957, -0.4227313679354475217664), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(1.298915541971320308789, 0.293508923629326790028), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.1655828816811988069002, 0.1530695261669812456962, 0.1655828816811988069002, 0.1530695261669812456962), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1358401620200463122412, 0.1382021771627104100855), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.2451746600755327076815, 0.2941698864686516090572), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -133.938827675498004055, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 280.6436127978045078635, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 297.272358665813612788, tolerance = tolerance))
     >
     > res2 <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_ar1(~1 + X1),
     + threshold.constraints = c(1,1),
     + coef.constraints = list(matrix(rep(1,4), ncol = 1), matrix(rep(1,4),
     + ncol = 1)))
     >
     > mvord:::check(all.equal(res$beta, res2$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res2$sebeta, tolerance = tolerance))
     >
     > ########################################################################################
     > res <- mvord(formula = MMO2(Y1,Y2) ~ 0 + X1 + X2,
     + data = data_toy_example,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = cbind(c(1,2),c(NA,1)))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord(formula = MMO2(Y1, Y2) ~ 0 + X1 + X2, data = data_toy_example,
     error.structure = cor_general(~1), link = mvprobit(), coef.constraints = cbind(c(1,
     2), c(NA, 1)), threshold.constraints = c(1, 1), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO2(Y1, Y2) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -137.65 288.06 304.69 63
    
     Thresholds:
     Estimate Std. Error z value
     Y1 1|2 -0.89892942527340835568 0.16246077762271224354 -5.5332100000000004059
     Y1 2|3 0.98413365940583297231 0.15981657213961053543 6.1578900000000000858
     Y2 1|2 -0.89892942527340835568 0.16246077762271224354 -5.5332100000000004059
     Y2 2|3 0.98413365940583286129 0.15981657213961053543 6.1578900000000000858
     Pr(>|z|)
     Y1 1|2 3.1442e-08 ***
     Y1 2|3 7.3718e-10 ***
     Y2 1|2 3.1442e-08 ***
     Y2 2|3 7.3718e-10 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     X1 1 0.681944859288649340989 0.174961892922603001743 3.89767999999999981142
     X1 2 0.468376583108446209458 0.157106006693623234671 2.98127999999999993008
     X2 1 -0.052179801987792852336 0.089798215416079740780 -0.58108000000000004093
     Pr(>|z|)
     X1 1 9.712e-05 ***
     X1 2 0.0028705 **
     X2 1 0.5611876
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr Y1 Y2 0.890390657144670694123 0.063175533024021041095 14.09392000000000067
     Pr(>|z|)
     corr Y1 Y2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.8989294252736504953205, 0.9841336594057152886705, -0.8989294252736503842982, 0.9841336594057153996928), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.68194485928859438494953, 0.46837658310856356003171, -0.05217980198794509860694), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.8903906571446290607597), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.1624607776227008359982, 0.1598165721396140603883, 0.1624607776227008359982, 0.1598165721396140603883), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.17496189292260000414103, 0.15710600669362786985239, 0.08979821541610534529898), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.0631755330240382911855), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -137.6494615446075613363, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 288.064880536023622426, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 304.6936264040327273506, tolerance = tolerance))
     >
     >
     > res2 <- mvord:::mvord(formula = MMO2(Y1,Y2) ~ 0 + X1 + X2,
     + data = data_toy_example,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = list(X1 = cbind(c(1,1,0,0), c(0,0,1,1)),
     + X2 = matrix(c(rep(0,2),rep(1,2)), ncol = 1)))
     >
     > mvord:::check(all.equal(res$beta, res2$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res2$sebeta, tolerance = tolerance))
     >
     > ########################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + offset(X2),
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,2),
     + coef.constraints = list(X1 = cbind(c(1,1,0,0), c(0,0,1,1))))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + offset(X2), data = df,
     error.structure = cor_general(~1), link = mvprobit(), coef.constraints = list(X1 = cbind(c(1,
     1, 0, 0), c(0, 0, 1, 1))), threshold.constraints = c(1,
     2), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + offset(X2)
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -194.33 403.71 423.31 78
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -0.76726526230273051077 0.12023341892343347714 -6.3814599999999996882
     1 2|3 0.89190483010004850684 0.13058942975202986192 6.8298399999999999110
     2 1|2 -0.79675344626870869824 0.11172347359918323451 -7.1314799999999998192
     2 2|3 0.87294080476149371606 0.13662238116424227363 6.3894399999999995643
     Pr(>|z|)
     1 1|2 1.7540e-10 ***
     1 2|3 8.5010e-12 ***
     2 1|2 9.9297e-13 ***
     2 2|3 1.6649e-10 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     X1 1 0.53370081959015414075 0.13636430505433183580 3.9137900000000001022
     X1 2 0.34385576067174200565 0.10634167794017361508 3.2334999999999998188
     Pr(>|z|)
     X1 1 9.086e-05 ***
     X1 2 0.0012228 **
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.940808527820466200531 0.033697206769499193912 27.919480000000000075
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.7672652623027305107684, 0.8919048301000485068357, -0.7967534462687086982413, 0.8729408047614937160574), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.5337008195901541407480, 0.3438557606717420056519), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.9408085278204662005308), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.1202334189234334493879, 0.1305894297520299174309, 0.1117234735991832761393, 0.1366223811642423013879), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1363643050543318635537, 0.1063416779401736012023), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.03369720676949919391241), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -194.3258468555947047207, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 403.705457152049632441, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 423.314265003572927526, tolerance = tolerance))
     >
     > res2 <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1,
     + offset = list(df$X2[1:100], df$X2[101:200]),
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,2),
     + coef.constraints = list(X1 = cbind(c(1,1,0,0), c(0,0,1,1))))
     >
     > mvord:::check(all.equal(res$beta, res2$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res2$sebeta, tolerance = tolerance))
     >
     > res3 <- mvord:::mvord(formula = MMO2(Y1,Y2) ~ 0 + X1 + offset(X2),
     + data = data_toy_example,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,2),
     + coef.constraints = list(X1 = cbind(c(1,1,0,0), c(0,0,1,1))))
     > mvord:::check(all.equal(res$beta, res3$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res3$sebeta, tolerance = tolerance))
     >
     > ########################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvlogit(),
     + #solver = "newuoa",
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = cbind(c(NA, NA), c(1,2)),
     + coef.values = cbind(c(1, 1), c(NA,NA)))
     > #rho <- res$rho
     >
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cor_general(~1),
     link = mvlogit(), coef.constraints = cbind(c(NA, NA), c(1,
     2)), coef.values = cbind(c(1, 1), c(NA, NA)), threshold.constraints = c(1,
     1))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvlogit flexible 100 2 -135.52 281.57 295.28 166
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -1.58328891554802453356 0.25851084829875303761 -6.1246499999999999275
     1 2|3 1.75864250142195510662 0.24251593708509780467 7.2516600000000002169
     2 1|2 -1.58328891554802453356 0.25851084829875303761 -6.1246499999999999275
     2 2|3 1.75864250142195510662 0.24251593708509780467 7.2516600000000002169
     Pr(>|z|)
     1 1|2 9.0882e-10 ***
     1 2|3 4.1170e-13 ***
     2 1|2 9.0882e-10 ***
     2 2|3 4.1170e-13 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     X2 1 -0.77306084234178251702 0.26732809287377440333 -2.8918099999999999916
     X2 2 -0.72430326757362628598 0.24403984199061926064 -2.9679700000000002191
     Pr(>|z|)
     X2 1 0.0038304 **
     X2 2 0.0029977 **
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.855894446593468805062 0.060599195266431746254 14.123860000000000525
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     >
     > options(digits = 22)
     >
     > tolerance2 <- 1e-4
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-1.583288915548024533564, 1.758642501421955106622, -1.583288915548024533564, 1.758642501421955106622), tolerance = tolerance2))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(-0.7730608423417825170176, -0.7243032675736262859800), tolerance = tolerance2))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.8558944465934688050623), tolerance = tolerance2))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.2585108482987530376107, 0.2425159370850978601819, 0.2585108482987530376107, 0.2425159370850978601819), tolerance = tolerance2))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.2673280928737744588375, 0.2440398419906193439033), tolerance = tolerance2))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.06059919526643173237623), tolerance = tolerance2))
     > mvord:::check(all.equal(logLik(res)[[1]], -135.5210996495467554723, tolerance = tolerance2))
     > mvord:::check(all.equal(AIC(res), 281.568515088567210114, tolerance = tolerance2))
     > mvord:::check(all.equal(BIC(res), 295.2799371200834457341, tolerance = tolerance2))
     >
     > # res1 <- mvord:::mvord(formula = MMO(Y) ~ 1 + X1 + X2,
     > # data = df,
     > # index = c("i", "j"),
     > # link = mvprobit(),
     > # solver = "BFGS",
     > # se = TRUE,
     > # error.structure = cor_equi(~1),
     > # threshold.constraints = c(1,1),
     > # coef.constraints = list(Intercept = matrix(rep(1,6), ncol = 1),
     > # X2 = cbind(c(1,1,1,0,0,0), c(0,0,0,1,1,1)), X1 = matrix(rep(1,6), ncol = 1)))
     >
     >
     > ########################################################################################
     > df$X3 <- cut(df$X2, c(-Inf, -0.2, 0.2, Inf))
     >
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X3,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = T,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = c(1,1))
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X3, data = df, error.structure = cor_general(~1),
     link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X3
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -136.05 284.87 301.5 34
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -1.21270698725061820689 0.21131407873814661569 -5.7388799999999999812
     1 2|3 0.75845153901163941956 0.18860790313556302644 4.0213099999999997181
     2 1|2 -1.21270698725061820689 0.21131407873814661569 -5.7388799999999999812
     2 2|3 0.75845153901163941956 0.18860790313556302644 4.0213099999999997181
     Pr(>|z|)
     1 1|2 9.5302e-09 ***
     1 2|3 5.7874e-05 ***
     2 1|2 9.5302e-09 ***
     2 2|3 5.7874e-05 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error
     X1 1 0.582346304092615452142 0.133874274950941285489
     X3(-0.2,0.2] 1 0.065949564659511766829 0.394250784463482328857
     X3(0.2, Inf] 1 -0.639867558058318963710 0.239596299509496207802
     z value Pr(>|z|)
     X1 1 4.34994999999999976126 1.3617e-05 ***
     X3(-0.2,0.2] 1 0.16728000000000001202 0.8671512
     X3(0.2, Inf] 1 -2.67060999999999992838 0.0075714 **
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.857975765661140976626 0.062391829886069061217 13.751409999999999911
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     >
     > options(digits = 22)
     >
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-1.212706987250618206886, 0.758451539011639419563, -1.212706987250618206886, 0.758451539011639419563), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.58234630409261545214150, 0.06594956465951176682871, -0.63986755805831896370961), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.8579757656611409766256), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.2113140787381464491546, 0.1886079031355630819533, 0.2113140787381464491546, 0.1886079031355630819533), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1338742749509412854891, 0.3942507844634823288565, 0.2395962995094961522913), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.06239182988606906121731), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -136.0526219854373835005, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 284.8712014176832667545, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 301.499947285692371679, tolerance = tolerance))
     > ########################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 1 + X1 * X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = T,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + threshold.values = list(c(-1,NA),
     + c(-1,NA)),
     + coef.constraints = c(1,1))
     > res.summary <- summary(res, short = TRUE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 1 + X1 * X2, data = df, error.structure = cor_general(~1),
     link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), threshold.values = list(c(-1, NA), c(-1, NA)), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 1 + X1 * X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit fix1first 100 2 -134.63 282.02 298.65 34
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -1.00000000000000000000 0.00000000000000000000 NA
     1 2|3 1.00269513578732416548 0.22781659168150072969 4.4013299999999997425
     Pr(>|z|)
     1 1|2 NA
     1 2|3 1.0759e-05 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error
     (Intercept) 1 -0.021913582738025384755 0.167100624889530591233
     X1 1 0.621621046245967145971 0.138540865251019801319
     X2 1 -0.427624354284268537452 0.140667605271754703189
     X1:X2 1 -0.102567255999721287929 0.150381573286269054623
     z value Pr(>|z|)
     (Intercept) 1 -0.13114000000000000656 0.8956645
     X1 1 4.48690999999999995396 7.2262e-06 ***
     X2 1 -3.03996000000000021757 0.0023661 **
     X1:X2 1 -0.68205000000000004512 0.4952094
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.853574805224538435411 0.063090989221389948138 13.529270000000000351
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     >
     > options(digits = 22)
     >
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-1.000000000000000000000, 1.002695135787324165477), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(-0.02191358273802538475516, 0.62162104624596714597118, -0.42762435428426853745165, -0.10256725599972128792903), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.8535748052245384354109), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.0000000000000000000000, 0.2278165916815008962271), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1671006248895306467439, 0.1385408652510198845853, 0.1406676052717547309445, 0.1503815732862690546234), tolerance = tolerance))
     > #mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.063097735777156036), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.06309098922138993426056), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -134.6255603757719541136, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 282.0170781983524079806, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 298.6458240663615129051, tolerance = tolerance))
     > ########################################################################################
     > df_NA <- df[-c(1,90:110),]
     >
     >
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df_NA,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = T,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,2)
     + #coef.constraints = list(X1 = cbind(c(1,1,1,0,0,0), c(0,0,0,1,1,1)))
     + )
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df_NA, error.structure = cor_general(~1),
     link = mvprobit(), threshold.constraints = c(1, 2), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 99 2 -119.45 258.71 284.4 51
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -1.02012637737423839113 0.19996641229613928981 -5.1014900000000000801
     1 2|3 1.14180364925241861762 0.22379090693938255563 5.1021000000000000796
     2 1|2 -0.90662822300388590246 0.19059651595705001670 -4.7567899999999996297
     2 2|3 0.99945117269142824679 0.17403069740199303417 5.7429600000000000648
     Pr(>|z|)
     1 1|2 3.3699e-07 ***
     1 2|3 3.3590e-07 ***
     2 1|2 1.9669e-06 ***
     2 2|3 9.3036e-09 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     X1 1 0.83733616646585673493 0.18898470050447499502 4.4307100000000003703
     X1 2 0.49821075302580530852 0.16630384260856759249 2.9957899999999999530
     X2 1 -0.44740279449498343567 0.18763575165139206868 -2.3844199999999999839
     X2 2 -0.35390389963398127815 0.14489624593604019664 -2.4424600000000000755
     Pr(>|z|)
     X1 1 9.3924e-06 ***
     X1 2 0.0027374 **
     X2 1 0.0171060 *
     X2 2 0.0145874 *
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.910650307820911830703 0.070553566078614848855 12.907220000000000582
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate,
     + c(-1.0201263773742383911269, 1.1418036492524186176212, -0.9066282230038859024646, 0.9994511726914282467860), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.8373361664658567349306, 0.4982107530258053085248, -0.4474027944949834356692, -0.3539038996339812781500),
     + tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.9106503078209118307029), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.1999664122961394840949, 0.2237909069393826111405, 0.1905965159570501277209, 0.1740306974019930341679), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1889847005044749950198, 0.1663038426085674814647, 0.1876357516513920686840, 0.1448962459360401411335), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.07055356607861483497768), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -119.4526280916837492896, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 258.7052561833675099479, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 284.3969426996999345647, tolerance = tolerance))
     > ##########################################################################################
     > df_23 <- df
     >
     > df_23$Y[which(df_23$Y[1:100] == 3)] <- 2
     >
     >
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df_23,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,2))
     Error in check_response_missings(rho) :
     Model is not identifiable. For at least one dimension, not all response categories are observed.
     Either remove the non-observed category or set constraints on the threshold parameters to ensure identifiability.
     Calls: <Anonymous> -> check_response_missings
     Execution halted
Flavor: r-devel-linux-x86_64-debian-clang

Version: 1.0.0
Check: tests
Result: ERROR
     Running ‘check_methods.R’ [4s/6s]
     Running ‘check_mvord_data.R’ [2s/3s]
     Running ‘check_mvord_errors.R’ [2s/3s]
     Running ‘check_set_threshold_type.R’ [2s/3s]
     Running ‘check_toy_example.R’ [24s/35s]
     Running ‘check_transf_sigmas.R’ [2s/3s]
     Running ‘check_transf_thresholds.R’ [2s/3s]
    Running the tests in ‘tests/check_toy_example.R’ failed.
    Complete output:
     > library(mvord)
     Loading required package: minqa
     Loading required package: BB
     Loading required package: ucminf
     Loading required package: dfoptim
     > library(MASS)
     >
     >
     > #data(data_toy_example)
     > tolerance <- 1e-6
     >
     > mult.obs <- 2
     > sigma <- matrix(c(1,0.8,0.8,1), ncol = 2)
     > betas <- list(c(0.8,-0.5),
     + c(0.8,-0.5))
     > thresholds <- list(c(-1,1),c(-1,1))
     > nobs <- 100
     > suppressWarnings(RNGversion("3.5.0"))
     > set.seed(2017)
     > errors <- mvrnorm(n = nobs, mu = rep(0, mult.obs), Sigma = sigma)
     >
     > X1 <- rnorm(nobs, 0, 1)
     > X2 <- rnorm(nobs, 0, 1)
     >
     > pred <- cbind(X1, X2)
     >
     > y <- sapply(1:mult.obs, function(j) pred %*% betas[[j]] + errors[, j], simplify = "array")
     > y.ord <- sapply(1:mult.obs, function(j) cut(y[, , j], c(min(y[, , j]) - 1, c(thresholds[[j]]), max(y[, , j]) + 1),
     + labels = FALSE), simplify = "array")
     > predictors.fixed <- lapply(1:mult.obs, function(j) pred)
     > y <- as.data.frame(y.ord)
     >
     > for(i in 1:mult.obs){
     + y[, i] <- factor(y[, i], levels = sort(unique(y[, i])),
     + ordered = TRUE)
     + }
     >
     >
     >
     >
     > data_toy_example <- cbind.data.frame(y, predictors.fixed[[1]])
     > colnames(data_toy_example) <- c("Y1","Y2", "X1", "X2")
     > w <- c(rep(1/20, 20), rep(1/30, 30), rep(1/20, 20), rep(1/30, 30))
     >
     >
     >
     > # convert data_toy_example into long format
     > df <- cbind.data.frame("i" = rep(1:100,2),
     + "j" = rep(1:2,each = 100),
     + "Y" = c(data_toy_example$Y1,data_toy_example$Y2),
     + "X1" = rep(data_toy_example$X1,2),
     + "X2" = rep(data_toy_example$X2,2),
     + "f1" = factor(sample(rep(data_toy_example$Y2,2)), ordered =F),
     + "f2" = factor(rep(data_toy_example$Y1,2), ordered=F),
     + w=rep(w,2))
     >
     >
     >
     > # library(ROI)
     > # ROI_solver <- function(starting.values, objFun, control){
     > # n <- length(starting.values)
     > # op <- OP(objective = F_objective(objFun, n = n),
     > # bounds = V_bound(li = seq_len(n), lb = rep.int(-Inf, n)))
     > # optRes <- ROI_solve(op, solver = "nlminb",
     > # control = c(list(start = starting.values),
     > # control))
     > # list(optpar = optRes$solution,
     > # objective = optRes$objval) # a vector of length equal to number of parameters to optimize
     > # }
     > #
     > #
     > #
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + link = mvprobit(),
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = c(1,1),
     + contrasts = list(f1 = function(x)
     + contr.treatment(nlevels(df$f1), base = 1),
     + f2 = "contr.sum"),
     + control= mvord.control(solver="BFGS",se=T))
     Warning message:
     In mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, link = mvprobit(), :
     Variables f1 and f2 are absent, the contrasts will be ignored.
     > options(digits = 22)
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cor_general(~1),
     link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), contrasts = list(f1 = function(x) contr.treatment(nlevels(df$f1),
     base = 1), f2 = "contr.sum"), control = mvord.control(solver = "BFGS",
     se = T))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -134.91 280.34 294.06 32
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -0.96257386663519672876 0.16613952738605439197 -5.7937700000000003087
     1 2|3 1.03347036873223707687 0.15004482537617938598 6.8877399999999999736
     2 1|2 -0.96257386663519672876 0.16613952738605439197 -5.7937700000000003087
     2 2|3 1.03347036873223707687 0.15004482537617938598 6.8877399999999999736
     Pr(>|z|)
     1 1|2 6.8825e-09 ***
     1 2|3 5.6684e-12 ***
     2 1|2 6.8825e-09 ***
     2 2|3 5.6684e-12 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     X1 1 0.63801035062404309883 0.13576747301357891540 4.6992900000000004113
     X2 1 -0.42672524816263474046 0.13643665802446264257 -3.1276399999999999757
     Pr(>|z|)
     X1 1 2.6107e-06 ***
     X2 1 0.0017621 **
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.854266728221226845363 0.062440889990360744222 13.681210000000000093
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.96257386663519672876, 1.03347036873223707687, -0.96257386663519672876, 1.03347036873223707687), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.63801035062404309883, -0.42672524816263474046), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.85426672822122684536), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.16613952738605441972, 0.15004482537617935822, 0.16613952738605441972, 0.15004482537617935822), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.13576747301357894315, 0.13643665802446261481), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.062440889990360744222), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -134.90867383086322207, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 280.3436634512001433, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 294.05508548271637892, tolerance = tolerance))
     >
     >
     > res2 <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = list(matrix(rep(1,4), ncol = 1), matrix(rep(1,4), ncol = 1)))
     >
     > mvord:::check(all.equal(res$beta, res2$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res2$sebeta, tolerance = tolerance))
     > ########################################################################
     > ## No coefficients
     > res <- mvord:::mvord(formula = MMO(Y) ~ -1,
     + data = df,
     + link = mvprobit(),
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + control= mvord.control(solver="BFGS",se=T))
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ -1, data = df, error.structure = cor_general(~1),
     link = mvprobit(), threshold.constraints = c(1, 1), control = mvord.control(solver = "BFGS",
     se = T))
    
     Formula: MMO(Y) ~ -1
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -153.66 313.51 321.57 26
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -0.75479706538110091785 0.13105603927553963195 -5.7593500000000004135
     1 2|3 0.86086304364935783973 0.13360642334434200129 6.4432799999999996743
     2 1|2 -0.75479706538110091785 0.13105603927553963195 -5.7593500000000004135
     2 2|3 0.86086304364935783973 0.13360642334434200129 6.4432799999999996743
     Pr(>|z|)
     1 1|2 8.4440e-09 ***
     1 2|3 1.1692e-10 ***
     2 1|2 8.4440e-09 ***
     2 2|3 1.1692e-10 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.905795171446422409112 0.038855435949034601573 23.311930000000000263
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.75479706538110091785, 0.86086304364935783973, -0.75479706538110091785, 0.86086304364935783973), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(numeric(0)), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.90579517144642240911), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.13105603927553965971, 0.13360642334434202905, 0.13105603927553965971, 0.13360642334434202905), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(numeric(0)), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.038855435949034601573), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -153.66397119528727444, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 313.51350940088383368, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 321.57073678022845797, tolerance = tolerance))
     >
     > #polychor
     > res <- mvord:::mvord(formula = MMO(Y) ~ 1,
     + data = df,
     + link = mvprobit(),
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + control= mvord.control(solver="BFGS",se=T))
     Note: First threshold for each response is fixed to 0 in order to ensure identifiability!
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 1, data = df, error.structure = cor_general(~1),
     link = mvprobit(), threshold.constraints = c(1, 1), control = mvord.control(solver = "BFGS",
     se = T))
    
     Formula: MMO(Y) ~ 1
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit fix1first 100 2 -153.56 315.46 326.32 52
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 0.00000000000000000000 0.00000000000000000000 NA
     1 2|3 1.61576075469784496974 0.15607113684416781818 10.3527199999999997
     2 1|2 0.00000000000000000000 0.00000000000000000000 NA
     2 2|3 1.61576075469784496974 0.15607113684416781818 10.3527199999999997
     Pr(>|z|)
     1 1|2 NA
     1 2|3 < 2.22e-16 ***
     2 1|2 NA
     2 2|3 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error
     (Intercept) 1 0.73786225358709867095 0.13500906208906057748
     (Intercept) 2 0.77255634837212872057 0.14063945843173736305
     z value Pr(>|z|)
     (Intercept) 1 5.4652799999999999159 4.6218e-08 ***
     (Intercept) 2 5.4931700000000001083 3.9478e-08 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.906373932662642656233 0.038916052675185490439 23.290489999999998361
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(0.0000000000000000000, 1.6157607546978449697, 0.0000000000000000000, 1.6157607546978449697), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.73786225358709867095, 0.77255634837212872057), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.90637393266264265623), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.00000000000000000000, 0.15607113684416776267, 0.00000000000000000000, 0.15607113684416776267), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.13500906208906049422, 0.14063945843173727979), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.038916052675185497378), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -153.5640565887688922, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 315.46144651087109878, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 326.31632228582151356, tolerance = tolerance))
     > #######################################################################
     > ## cor_general(~factor)
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + error.structure = cor_general(~f2),
     + threshold.constraints = c(1,1),
     + coef.constraints = c(1,1),
     + contrasts = list(f2 = "contr.sum"))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cor_general(~f2),
     link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), contrasts = list(f2 = "contr.sum"), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -134.17 283.39 303 38
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -0.90552506790519915469 0.18265647646915181279 -4.9575300000000002143
     1 2|3 1.00420547521745429087 0.17467164336594764862 5.7491000000000003212
     2 1|2 -0.90552506790519915469 0.18265647646915181279 -4.9575300000000002143
     2 2|3 1.00420547521745429087 0.17467164336594764862 5.7491000000000003212
     Pr(>|z|)
     1 1|2 7.1395e-07 ***
     1 2|3 8.9718e-09 ***
     2 1|2 7.1395e-07 ***
     2 2|3 8.9718e-09 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     X1 1 0.64793126338493944871 0.14153585035881874332 4.5778600000000002623
     X2 1 -0.42893128334048874484 0.13989943737309903926 -3.0659999999999998366
     Pr(>|z|)
     X1 1 4.6976e-06 ***
     X2 1 0.0021695 **
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error
     corr f21 1 2 0.73667859723415496376 0.17823839430849072740
     corr f22 1 2 0.91829826305234418804 0.06885865619995423792
     corr f23 1 2 0.83710611971260706632 0.13103608299872890330
     z value Pr(>|z|)
     corr f21 1 2 4.1331100000000002836 3.5789e-05 ***
     corr f22 1 2 13.3359900000000006770 < 2.22e-16 ***
     corr f23 1 2 6.3883599999999995944 1.6767e-10 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate,
     + c(-0.90552506790519915469, 1.00420547521745429087, -0.90552506790519915469, 1.00420547521745429087), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate,
     + c(0.64793126338493944871, -0.42893128334048874484), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate,
     + c(0.73667859723415496376, 0.91829826305234418804, 0.83710611971260706632), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`,
     + c(0.18265647646915181279, 0.17467164336594759311, 0.18265647646915181279, 0.17467164336594759311), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`,
     + c(0.14153585035881874332, 0.13989943737309898375), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`,
     + c(0.17823839430849069965, 0.06885865619995423792, 0.13103608299872890330), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -134.17046176709035876, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 283.39468697504094052, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 303.00349482656417877, tolerance = tolerance))
     >
     > ##################################################################################
     > res <- mvord:::mvord(formula = MMO(Y, i, j) ~ 0 + X1 + X2,
     + data = df,
     + link = mvlogit(df = 8L),
     + error.structure = cor_general(~1))
     > #threshold.constraints = c(1,1),
     > #coef.constraints = c(1,1))
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y, i, j) ~ 0 + X1 + X2, data = df,
     error.structure = cor_general(~1), link = mvlogit(df = 8L))
    
     Formula: MMO(Y, i, j) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvlogit flexible 100 2 -132.87 285.51 311.28 597
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -1.61993043593746333109 0.33738759349733610238 -4.8013899999999996027
     1 2|3 1.83247603989284058024 0.36742763311918208435 4.9873099999999999099
     2 1|2 -1.61590382450539626902 0.33963203525352064771 -4.7578100000000000946
     2 2|3 1.74613531250955356100 0.30112281186577066761 5.7987500000000000711
     Pr(>|z|)
     1 1|2 1.5757e-06 ***
     1 2|3 6.1225e-07 ***
     2 1|2 1.9571e-06 ***
     2 2|3 6.6812e-09 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     X1 1 1.28924362568011741104 0.30395050451031768723 4.2416200000000001680
     X1 2 0.88420512353335356526 0.29740834083412664990 2.9730300000000000615
     X2 1 -0.84291068182789685714 0.30080361920504455897 -2.8022000000000000242
     X2 2 -0.71568802554483457179 0.25697665338227909659 -2.7850299999999998946
     Pr(>|z|)
     X1 1 2.2191e-05 ***
     X1 2 0.0029487 **
     X2 1 0.0050756 **
     X2 2 0.0053523 **
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.883893395936291792481 0.071174687242751624727 12.418649999999999523
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > # mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-1.6170817306633420, 1.7855897338762188, -1.6170817306633420, 1.7855897338762188), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$coefficients$Estimate, c(1.07242200717115987, 1.07242200717115987, -0.76715925377701732, -0.76715925377701732), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.85317690560688531), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.28997405649500174, 0.27427389826231802, 0.28997405649500174, 0.27427389826231802), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.24111402270822993, 0.24111402270822993, 0.24156664773225886, 0.24156664773225886), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.063316529381183581), tolerance = tolerance))
     > # mvord:::check(all.equal(logLik(res), -135.41665313840898, tolerance = tolerance))
     > # mvord:::check(all.equal(AIC(res), 281.35962206629165, tolerance = tolerance))
     > # mvord:::check(all.equal(BIC(res), 295.07104409780789, tolerance = tolerance))
     >
     > ##################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cov_general(~1),
     + threshold.constraints = c(1,1),
     + threshold.values = list(c(-1,NA),
     + c(-1,NA)),
     + coef.constraints = c(1,1))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cov_general(~1),
     link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), threshold.values = list(c(-1, NA), c(-1, NA)), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit fix1first 100 2 -134.64 282.04 298.67 37
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -1.00000000000000000000 0.00000000000000000000 NA
     1 2|3 1.08266798008332409964 0.24926043418148299824 4.3435199999999998255
     2 1|2 -1.00000000000000000000 0.00000000000000000000 NA
     2 2|3 1.08266798008332409964 0.24926043418148299824 4.3435199999999998255
     Pr(>|z|)
     1 1|2 NA
     1 2|3 1.4022e-05 ***
     2 1|2 NA
     2 2|3 1.4022e-05 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     X1 1 0.68834946142083863752 0.15128036741851591529 4.5501599999999999824
     X2 1 -0.45406517655023831415 0.15704134440906966641 -2.8913700000000002177
     Pr(>|z|)
     X1 1 5.3606e-06 ***
     X2 1 0.0038356 **
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.856153637803131029038 0.061822221469835146601 13.8486399999999996169
     sigma 1 1.001934583868862471689 0.185187960683899310865 5.4103700000000003456
     sigma 2 1.089828199315680645753 0.197593289723817816528 5.5155099999999999127
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     sigma 1 6.2896e-08 ***
     sigma 2 3.4777e-08 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-1.000000000000000000000, 1.082667980083374503764, -1.000000000000000000000, 1.082667980083374503764), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.6883494614209317852271, -0.4540651765505133163892), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.8561536378031743277361, 1.0019345838689188710191, 1.0898281993157663549709), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.0000000000000000000000, 0.2492604341815170820862, 0.0000000000000000000000, 0.2492604341815170820862), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1512803674185349001036, 0.1570413444091093568833), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.06182222146981963123435, 0.18518796068392143205905, 0.19759328972385112321852), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -134.6391112878561102661, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 282.0441800225207202857, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 298.6729258905298252103, tolerance = tolerance))
     > ########################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_equi(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = cbind(c(1,1),c(1,2)))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cor_equi(~1),
     link = mvprobit(), coef.constraints = cbind(c(1, 1), c(1,
     2)), threshold.constraints = c(1, 1), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -134.84 282.45 299.08 35
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -0.96275665235846763501 0.16703865743956475276 -5.7636799999999999145
     1 2|3 1.03388029135962500327 0.15203091606549479220 6.8004600000000001714
     2 1|2 -0.96275665235846763501 0.16703865743956475276 -5.7636799999999999145
     2 2|3 1.03388029135962478122 0.15203091606549479220 6.8004600000000001714
     Pr(>|z|)
     1 1|2 8.2301e-09 ***
     1 2|3 1.0429e-11 ***
     2 1|2 8.2301e-09 ***
     2 2|3 1.0429e-11 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     X1 1 0.63820467410958414689 0.13670492809263015688 4.6684799999999997411
     X2 1 -0.44676561119553748203 0.15900879662664607617 -2.8096899999999997988
     X2 2 -0.40750155519132619242 0.13850979189699536009 -2.9420399999999999885
     Pr(>|z|)
     X1 1 3.0343e-06 ***
     X2 1 0.0049589 **
     X2 2 0.0032606 **
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     (Intercept) 1.27290057125462996446 0.23221611847930476169 5.4815300000000002356
     Pr(>|z|)
     (Intercept) 4.2165e-08 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.9627566523584676350112, 1.0338802913596250032668, -0.9627566523584676350112, 1.0338802913596247812222), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.6382046741095841468905, -0.4467656111955374820255, -0.4075015551913261924177), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$error.structure$Estimate, c(1.272900571254629964457), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.1670386574395647250046, 0.1520309160654948199554, 0.1670386574395647250046, 0.1520309160654948199554), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1367049280926301846328, 0.1590087966266460206555, 0.1385097918969951935608), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.2322161184793045674013), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -134.8432321319771745038, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 282.4524217107628487611, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 299.0811675787719536856, tolerance = tolerance))
     >
     > res2 <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_equi(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = list(X2 = cbind(c(1,1,0,0), c(0,0,1,1)),
     + X1 = matrix(rep(1,4), ncol = 1)))
     >
     > mvord:::check(all.equal(res$beta, res2$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res2$sebeta, tolerance = tolerance))
     >
     > ########################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_ar1(~ 1 + X1),
     + threshold.constraints = c(1,1),
     + coef.constraints = c(1,1))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cor_ar1(~1 +
     X1), link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -133.94 280.64 297.27 35
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -0.95722951152617552495 0.16558288168119875139 -5.7809699999999999420
     1 2|3 1.03746793968331618707 0.15306952616698107916 6.7777599999999997848
     2 1|2 -0.95722951152617552495 0.16558288168119875139 -5.7809699999999999420
     2 2|3 1.03746793968331618707 0.15306952616698107916 6.7777599999999997848
     Pr(>|z|)
     1 1|2 7.4272e-09 ***
     1 2|3 1.2206e-11 ***
     2 1|2 7.4272e-09 ***
     2 2|3 1.2206e-11 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     X1 1 0.65304208674282193670 0.13584016202004611795 4.8074300000000000921
     X2 1 -0.42273136793544752177 0.13820217716271032682 -3.0587900000000001199
     Pr(>|z|)
     X1 1 1.5288e-06 ***
     X2 1 0.0022223 **
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error
     (Intercept) 1.29891554197132030879 0.24517466007553267993
     X1 0.29350892362932679003 0.29416988646865149803
     z value Pr(>|z|)
     (Intercept) 5.29792000000000040671 1.1713e-07 ***
     X1 0.99775000000000002576 0.3184
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.9572295115261755249492, 1.0374679396833161870717, -0.9572295115261755249492, 1.0374679396833161870717), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.6530420867428219366957, -0.4227313679354475217664), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(1.298915541971320308789, 0.293508923629326790028), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.1655828816811988069002, 0.1530695261669812456962, 0.1655828816811988069002, 0.1530695261669812456962), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1358401620200463122412, 0.1382021771627104100855), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.2451746600755327076815, 0.2941698864686516090572), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -133.938827675498004055, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 280.6436127978045078635, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 297.272358665813612788, tolerance = tolerance))
     >
     > res2 <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_ar1(~1 + X1),
     + threshold.constraints = c(1,1),
     + coef.constraints = list(matrix(rep(1,4), ncol = 1), matrix(rep(1,4),
     + ncol = 1)))
     >
     > mvord:::check(all.equal(res$beta, res2$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res2$sebeta, tolerance = tolerance))
     >
     > ########################################################################################
     > res <- mvord(formula = MMO2(Y1,Y2) ~ 0 + X1 + X2,
     + data = data_toy_example,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = cbind(c(1,2),c(NA,1)))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord(formula = MMO2(Y1, Y2) ~ 0 + X1 + X2, data = data_toy_example,
     error.structure = cor_general(~1), link = mvprobit(), coef.constraints = cbind(c(1,
     2), c(NA, 1)), threshold.constraints = c(1, 1), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO2(Y1, Y2) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -137.65 288.06 304.69 63
    
     Thresholds:
     Estimate Std. Error z value
     Y1 1|2 -0.89892942527339636527 0.16246077762273250511 -5.5332100000000004059
     Y1 2|3 0.98413365940575758817 0.15981657213962935371 6.1578900000000000858
     Y2 1|2 -0.89892942527339636527 0.16246077762273250511 -5.5332100000000004059
     Y2 2|3 0.98413365940575769919 0.15981657213962935371 6.1578900000000000858
     Pr(>|z|)
     Y1 1|2 3.1442e-08 ***
     Y1 2|3 7.3718e-10 ***
     Y2 1|2 3.1442e-08 ***
     Y2 2|3 7.3718e-10 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     X1 1 0.681944859288474369841 0.174961892922601003342 3.89767999999999981142
     X1 2 0.468376583108382760212 0.157106006693620736669 2.98127999999999993008
     X2 1 -0.052179801987817860109 0.089798215416046725523 -0.58108000000000004093
     Pr(>|z|)
     X1 1 9.712e-05 ***
     X1 2 0.0028705 **
     X2 1 0.5611876
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr Y1 Y2 0.890390657144797814659 0.063175533023959645762 14.09392000000000067
     Pr(>|z|)
     corr Y1 Y2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.8989294252736504953205, 0.9841336594057152886705, -0.8989294252736503842982, 0.9841336594057153996928), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.68194485928859438494953, 0.46837658310856356003171, -0.05217980198794509860694), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.8903906571446290607597), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.1624607776227008359982, 0.1598165721396140603883, 0.1624607776227008359982, 0.1598165721396140603883), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.17496189292260000414103, 0.15710600669362786985239, 0.08979821541610534529898), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.0631755330240382911855), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -137.6494615446075613363, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 288.064880536023622426, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 304.6936264040327273506, tolerance = tolerance))
     >
     >
     > res2 <- mvord:::mvord(formula = MMO2(Y1,Y2) ~ 0 + X1 + X2,
     + data = data_toy_example,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = list(X1 = cbind(c(1,1,0,0), c(0,0,1,1)),
     + X2 = matrix(c(rep(0,2),rep(1,2)), ncol = 1)))
     >
     > mvord:::check(all.equal(res$beta, res2$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res2$sebeta, tolerance = tolerance))
     >
     > ########################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + offset(X2),
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,2),
     + coef.constraints = list(X1 = cbind(c(1,1,0,0), c(0,0,1,1))))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + offset(X2), data = df,
     error.structure = cor_general(~1), link = mvprobit(), coef.constraints = list(X1 = cbind(c(1,
     1, 0, 0), c(0, 0, 1, 1))), threshold.constraints = c(1,
     2), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + offset(X2)
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -194.33 403.71 423.31 78
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -0.76726526230273051077 0.12023341892343347714 -6.3814599999999996882
     1 2|3 0.89190483010004850684 0.13058942975202986192 6.8298399999999999110
     2 1|2 -0.79675344626870869824 0.11172347359918323451 -7.1314799999999998192
     2 2|3 0.87294080476149371606 0.13662238116424227363 6.3894399999999995643
     Pr(>|z|)
     1 1|2 1.7540e-10 ***
     1 2|3 8.5010e-12 ***
     2 1|2 9.9297e-13 ***
     2 2|3 1.6649e-10 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     X1 1 0.53370081959015414075 0.13636430505433183580 3.9137900000000001022
     X1 2 0.34385576067174200565 0.10634167794017361508 3.2334999999999998188
     Pr(>|z|)
     X1 1 9.086e-05 ***
     X1 2 0.0012228 **
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.940808527820466200531 0.033697206769499193912 27.919480000000000075
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.7672652623027305107684, 0.8919048301000485068357, -0.7967534462687086982413, 0.8729408047614937160574), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.5337008195901541407480, 0.3438557606717420056519), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.9408085278204662005308), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.1202334189234334493879, 0.1305894297520299174309, 0.1117234735991832761393, 0.1366223811642423013879), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1363643050543318635537, 0.1063416779401736012023), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.03369720676949919391241), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -194.3258468555947047207, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 403.705457152049632441, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 423.314265003572927526, tolerance = tolerance))
     >
     > res2 <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1,
     + offset = list(df$X2[1:100], df$X2[101:200]),
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,2),
     + coef.constraints = list(X1 = cbind(c(1,1,0,0), c(0,0,1,1))))
     >
     > mvord:::check(all.equal(res$beta, res2$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res2$sebeta, tolerance = tolerance))
     >
     > res3 <- mvord:::mvord(formula = MMO2(Y1,Y2) ~ 0 + X1 + offset(X2),
     + data = data_toy_example,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,2),
     + coef.constraints = list(X1 = cbind(c(1,1,0,0), c(0,0,1,1))))
     > mvord:::check(all.equal(res$beta, res3$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res3$sebeta, tolerance = tolerance))
     >
     > ########################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvlogit(),
     + #solver = "newuoa",
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = cbind(c(NA, NA), c(1,2)),
     + coef.values = cbind(c(1, 1), c(NA,NA)))
     > #rho <- res$rho
     >
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cor_general(~1),
     link = mvlogit(), coef.constraints = cbind(c(NA, NA), c(1,
     2)), coef.values = cbind(c(1, 1), c(NA, NA)), threshold.constraints = c(1,
     1))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvlogit flexible 100 2 -135.52 281.57 295.28 166
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -1.58328891554802453356 0.25851084829875303761 -6.1246499999999999275
     1 2|3 1.75864250142195510662 0.24251593708509780467 7.2516600000000002169
     2 1|2 -1.58328891554802453356 0.25851084829875303761 -6.1246499999999999275
     2 2|3 1.75864250142195510662 0.24251593708509780467 7.2516600000000002169
     Pr(>|z|)
     1 1|2 9.0882e-10 ***
     1 2|3 4.1170e-13 ***
     2 1|2 9.0882e-10 ***
     2 2|3 4.1170e-13 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     X2 1 -0.77306084234178251702 0.26732809287377440333 -2.8918099999999999916
     X2 2 -0.72430326757362628598 0.24403984199061926064 -2.9679700000000002191
     Pr(>|z|)
     X2 1 0.0038304 **
     X2 2 0.0029977 **
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.855894446593468805062 0.060599195266431746254 14.123860000000000525
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     >
     > options(digits = 22)
     >
     > tolerance2 <- 1e-4
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-1.583288915548024533564, 1.758642501421955106622, -1.583288915548024533564, 1.758642501421955106622), tolerance = tolerance2))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(-0.7730608423417825170176, -0.7243032675736262859800), tolerance = tolerance2))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.8558944465934688050623), tolerance = tolerance2))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.2585108482987530376107, 0.2425159370850978601819, 0.2585108482987530376107, 0.2425159370850978601819), tolerance = tolerance2))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.2673280928737744588375, 0.2440398419906193439033), tolerance = tolerance2))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.06059919526643173237623), tolerance = tolerance2))
     > mvord:::check(all.equal(logLik(res)[[1]], -135.5210996495467554723, tolerance = tolerance2))
     > mvord:::check(all.equal(AIC(res), 281.568515088567210114, tolerance = tolerance2))
     > mvord:::check(all.equal(BIC(res), 295.2799371200834457341, tolerance = tolerance2))
     >
     > # res1 <- mvord:::mvord(formula = MMO(Y) ~ 1 + X1 + X2,
     > # data = df,
     > # index = c("i", "j"),
     > # link = mvprobit(),
     > # solver = "BFGS",
     > # se = TRUE,
     > # error.structure = cor_equi(~1),
     > # threshold.constraints = c(1,1),
     > # coef.constraints = list(Intercept = matrix(rep(1,6), ncol = 1),
     > # X2 = cbind(c(1,1,1,0,0,0), c(0,0,0,1,1,1)), X1 = matrix(rep(1,6), ncol = 1)))
     >
     >
     > ########################################################################################
     > df$X3 <- cut(df$X2, c(-Inf, -0.2, 0.2, Inf))
     >
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X3,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = T,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = c(1,1))
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X3, data = df, error.structure = cor_general(~1),
     link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X3
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -136.05 284.87 301.5 34
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -1.21270698725061820689 0.21131407873814661569 -5.7388799999999999812
     1 2|3 0.75845153901163941956 0.18860790313556302644 4.0213099999999997181
     2 1|2 -1.21270698725061820689 0.21131407873814661569 -5.7388799999999999812
     2 2|3 0.75845153901163941956 0.18860790313556302644 4.0213099999999997181
     Pr(>|z|)
     1 1|2 9.5302e-09 ***
     1 2|3 5.7874e-05 ***
     2 1|2 9.5302e-09 ***
     2 2|3 5.7874e-05 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error
     X1 1 0.582346304092615452142 0.133874274950941285489
     X3(-0.2,0.2] 1 0.065949564659511766829 0.394250784463482328857
     X3(0.2, Inf] 1 -0.639867558058318963710 0.239596299509496207802
     z value Pr(>|z|)
     X1 1 4.34994999999999976126 1.3617e-05 ***
     X3(-0.2,0.2] 1 0.16728000000000001202 0.8671512
     X3(0.2, Inf] 1 -2.67060999999999992838 0.0075714 **
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.857975765661140976626 0.062391829886069061217 13.751409999999999911
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     >
     > options(digits = 22)
     >
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-1.212706987250618206886, 0.758451539011639419563, -1.212706987250618206886, 0.758451539011639419563), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.58234630409261545214150, 0.06594956465951176682871, -0.63986755805831896370961), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.8579757656611409766256), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.2113140787381464491546, 0.1886079031355630819533, 0.2113140787381464491546, 0.1886079031355630819533), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1338742749509412854891, 0.3942507844634823288565, 0.2395962995094961522913), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.06239182988606906121731), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -136.0526219854373835005, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 284.8712014176832667545, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 301.499947285692371679, tolerance = tolerance))
     > ########################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 1 + X1 * X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = T,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + threshold.values = list(c(-1,NA),
     + c(-1,NA)),
     + coef.constraints = c(1,1))
     > res.summary <- summary(res, short = TRUE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 1 + X1 * X2, data = df, error.structure = cor_general(~1),
     link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), threshold.values = list(c(-1, NA), c(-1, NA)), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 1 + X1 * X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit fix1first 100 2 -134.63 282.02 298.65 34
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -1.00000000000000000000 0.00000000000000000000 NA
     1 2|3 1.00269513578732416548 0.22781659168150072969 4.4013299999999997425
     Pr(>|z|)
     1 1|2 NA
     1 2|3 1.0759e-05 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error
     (Intercept) 1 -0.021913582738025384755 0.167100624889530591233
     X1 1 0.621621046245967145971 0.138540865251019801319
     X2 1 -0.427624354284268537452 0.140667605271754703189
     X1:X2 1 -0.102567255999721287929 0.150381573286269054623
     z value Pr(>|z|)
     (Intercept) 1 -0.13114000000000000656 0.8956645
     X1 1 4.48690999999999995396 7.2262e-06 ***
     X2 1 -3.03996000000000021757 0.0023661 **
     X1:X2 1 -0.68205000000000004512 0.4952094
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.853574805224538435411 0.063090989221389948138 13.529270000000000351
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     >
     > options(digits = 22)
     >
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-1.000000000000000000000, 1.002695135787324165477), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(-0.02191358273802538475516, 0.62162104624596714597118, -0.42762435428426853745165, -0.10256725599972128792903), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.8535748052245384354109), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.0000000000000000000000, 0.2278165916815008962271), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1671006248895306467439, 0.1385408652510198845853, 0.1406676052717547309445, 0.1503815732862690546234), tolerance = tolerance))
     > #mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.063097735777156036), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.06309098922138993426056), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -134.6255603757719541136, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 282.0170781983524079806, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 298.6458240663615129051, tolerance = tolerance))
     > ########################################################################################
     > df_NA <- df[-c(1,90:110),]
     >
     >
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df_NA,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = T,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,2)
     + #coef.constraints = list(X1 = cbind(c(1,1,1,0,0,0), c(0,0,0,1,1,1)))
     + )
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df_NA, error.structure = cor_general(~1),
     link = mvprobit(), threshold.constraints = c(1, 2), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 99 2 -119.45 258.71 284.4 51
    
     Thresholds:
     Estimate Std. Error z value
     1 1|2 -1.02012637737423839113 0.19996641229613928981 -5.1014900000000000801
     1 2|3 1.14180364925241861762 0.22379090693938255563 5.1021000000000000796
     2 1|2 -0.90662822300388590246 0.19059651595705001670 -4.7567899999999996297
     2 2|3 0.99945117269142824679 0.17403069740199303417 5.7429600000000000648
     Pr(>|z|)
     1 1|2 3.3699e-07 ***
     1 2|3 3.3590e-07 ***
     2 1|2 1.9669e-06 ***
     2 2|3 9.3036e-09 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     X1 1 0.83733616646585673493 0.18898470050447499502 4.4307100000000003703
     X1 2 0.49821075302580530852 0.16630384260856759249 2.9957899999999999530
     X2 1 -0.44740279449498343567 0.18763575165139206868 -2.3844199999999999839
     X2 2 -0.35390389963398127815 0.14489624593604019664 -2.4424600000000000755
     Pr(>|z|)
     X1 1 9.3924e-06 ***
     X1 2 0.0027374 **
     X2 1 0.0171060 *
     X2 2 0.0145874 *
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.910650307820911830703 0.070553566078614848855 12.907220000000000582
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001000000000000000020817 '**' 0.01000000000000000020817 '*'
     0.05000000000000000277556 '.' 0.1000000000000000055511 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate,
     + c(-1.0201263773742383911269, 1.1418036492524186176212, -0.9066282230038859024646, 0.9994511726914282467860), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.8373361664658567349306, 0.4982107530258053085248, -0.4474027944949834356692, -0.3539038996339812781500),
     + tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.9106503078209118307029), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.1999664122961394840949, 0.2237909069393826111405, 0.1905965159570501277209, 0.1740306974019930341679), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1889847005044749950198, 0.1663038426085674814647, 0.1876357516513920686840, 0.1448962459360401411335), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.07055356607861483497768), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -119.4526280916837492896, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 258.7052561833675099479, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 284.3969426996999345647, tolerance = tolerance))
     > ##########################################################################################
     > df_23 <- df
     >
     > df_23$Y[which(df_23$Y[1:100] == 3)] <- 2
     >
     >
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df_23,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,2))
     Error in check_response_missings(rho) :
     Model is not identifiable. For at least one dimension, not all response categories are observed.
     Either remove the non-observed category or set constraints on the threshold parameters to ensure identifiability.
     Calls: <Anonymous> -> check_response_missings
     Execution halted
Flavor: r-devel-linux-x86_64-debian-gcc

Version: 1.0.0
Check: running tests for arch ‘i386’
Result: ERROR
     Running 'check_methods.R' [6s]
     Running 'check_mvord_data.R' [3s]
     Running 'check_mvord_errors.R' [3s]
     Running 'check_set_threshold_type.R' [3s]
     Running 'check_toy_example.R' [35s]
     Running 'check_transf_sigmas.R' [3s]
     Running 'check_transf_thresholds.R' [3s]
    Running the tests in 'tests/check_toy_example.R' failed.
    Complete output:
     > library(mvord)
     Loading required package: minqa
     Loading required package: BB
     Loading required package: ucminf
     Loading required package: dfoptim
     > library(MASS)
     >
     >
     > #data(data_toy_example)
     > tolerance <- 1e-6
     >
     > mult.obs <- 2
     > sigma <- matrix(c(1,0.8,0.8,1), ncol = 2)
     > betas <- list(c(0.8,-0.5),
     + c(0.8,-0.5))
     > thresholds <- list(c(-1,1),c(-1,1))
     > nobs <- 100
     > suppressWarnings(RNGversion("3.5.0"))
     > set.seed(2017)
     > errors <- mvrnorm(n = nobs, mu = rep(0, mult.obs), Sigma = sigma)
     >
     > X1 <- rnorm(nobs, 0, 1)
     > X2 <- rnorm(nobs, 0, 1)
     >
     > pred <- cbind(X1, X2)
     >
     > y <- sapply(1:mult.obs, function(j) pred %*% betas[[j]] + errors[, j], simplify = "array")
     > y.ord <- sapply(1:mult.obs, function(j) cut(y[, , j], c(min(y[, , j]) - 1, c(thresholds[[j]]), max(y[, , j]) + 1),
     + labels = FALSE), simplify = "array")
     > predictors.fixed <- lapply(1:mult.obs, function(j) pred)
     > y <- as.data.frame(y.ord)
     >
     > for(i in 1:mult.obs){
     + y[, i] <- factor(y[, i], levels = sort(unique(y[, i])),
     + ordered = TRUE)
     + }
     >
     >
     >
     >
     > data_toy_example <- cbind.data.frame(y, predictors.fixed[[1]])
     > colnames(data_toy_example) <- c("Y1","Y2", "X1", "X2")
     > w <- c(rep(1/20, 20), rep(1/30, 30), rep(1/20, 20), rep(1/30, 30))
     >
     >
     >
     > # convert data_toy_example into long format
     > df <- cbind.data.frame("i" = rep(1:100,2),
     + "j" = rep(1:2,each = 100),
     + "Y" = c(data_toy_example$Y1,data_toy_example$Y2),
     + "X1" = rep(data_toy_example$X1,2),
     + "X2" = rep(data_toy_example$X2,2),
     + "f1" = factor(sample(rep(data_toy_example$Y2,2)), ordered =F),
     + "f2" = factor(rep(data_toy_example$Y1,2), ordered=F),
     + w=rep(w,2))
     >
     >
     >
     > # library(ROI)
     > # ROI_solver <- function(starting.values, objFun, control){
     > # n <- length(starting.values)
     > # op <- OP(objective = F_objective(objFun, n = n),
     > # bounds = V_bound(li = seq_len(n), lb = rep.int(-Inf, n)))
     > # optRes <- ROI_solve(op, solver = "nlminb",
     > # control = c(list(start = starting.values),
     > # control))
     > # list(optpar = optRes$solution,
     > # objective = optRes$objval) # a vector of length equal to number of parameters to optimize
     > # }
     > #
     > #
     > #
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + link = mvprobit(),
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = c(1,1),
     + contrasts = list(f1 = function(x)
     + contr.treatment(nlevels(df$f1), base = 1),
     + f2 = "contr.sum"),
     + control= mvord.control(solver="BFGS",se=T))
     Warning message:
     In mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, link = mvprobit(), :
     Variables f1 and f2 are absent, the contrasts will be ignored.
     > options(digits = 22)
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cor_general(~1),
     link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), contrasts = list(f1 = function(x) contr.treatment(nlevels(df$f1),
     base = 1), f2 = "contr.sum"), control = mvord.control(solver = "BFGS",
     se = T))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -134.91 280.34 294.06 32
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -0.96257386663519673 0.16613952738605439 -5.7937700000000003 6.8825e-09
     1 2|3 1.03347036873223708 0.15004482537617939 6.8877400000000000 5.6684e-12
     2 1|2 -0.96257386663519673 0.16613952738605439 -5.7937700000000003 6.8825e-09
     2 2|3 1.03347036873223708 0.15004482537617939 6.8877400000000000 5.6684e-12
    
     1 1|2 ***
     1 2|3 ***
     2 1|2 ***
     2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
     X1 1 0.63801035062404310 0.13576747301357892 4.6992900000000004 2.6107e-06
     X2 1 -0.42672524816263474 0.13643665802446264 -3.1276400000000000 0.0017621
    
     X1 1 ***
     X2 1 **
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value Pr(>|z|)
     corr 1 2 0.854266728221226845 0.062440889990360744 13.68121 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.96257386663519672876, 1.03347036873223707687, -0.96257386663519672876, 1.03347036873223707687), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.63801035062404309883, -0.42672524816263474046), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.85426672822122684536), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.16613952738605441972, 0.15004482537617935822, 0.16613952738605441972, 0.15004482537617935822), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.13576747301357894315, 0.13643665802446261481), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.062440889990360744222), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -134.90867383086322207, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 280.3436634512001433, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 294.05508548271637892, tolerance = tolerance))
     >
     >
     > res2 <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = list(matrix(rep(1,4), ncol = 1), matrix(rep(1,4), ncol = 1)))
     >
     > mvord:::check(all.equal(res$beta, res2$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res2$sebeta, tolerance = tolerance))
     > ########################################################################
     > ## No coefficients
     > res <- mvord:::mvord(formula = MMO(Y) ~ -1,
     + data = df,
     + link = mvprobit(),
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + control= mvord.control(solver="BFGS",se=T))
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ -1, data = df, error.structure = cor_general(~1),
     link = mvprobit(), threshold.constraints = c(1, 1), control = mvord.control(solver = "BFGS",
     se = T))
    
     Formula: MMO(Y) ~ -1
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -153.66 313.51 321.57 26
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -0.75479706538110092 0.13105603927553963 -5.7593500000000004 8.4440e-09
     1 2|3 0.86086304364935784 0.13360642334434200 6.4432799999999997 1.1692e-10
     2 1|2 -0.75479706538110092 0.13105603927553963 -5.7593500000000004 8.4440e-09
     2 2|3 0.86086304364935784 0.13360642334434200 6.4432799999999997 1.1692e-10
    
     1 1|2 ***
     1 2|3 ***
     2 1|2 ***
     2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
    
     Error Structure:
     Estimate Std. Error z value Pr(>|z|)
     corr 1 2 0.905795171446422409 0.038855435949034602 23.31193 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.75479706538110091785, 0.86086304364935783973, -0.75479706538110091785, 0.86086304364935783973), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(numeric(0)), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.90579517144642240911), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.13105603927553965971, 0.13360642334434202905, 0.13105603927553965971, 0.13360642334434202905), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(numeric(0)), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.038855435949034601573), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -153.66397119528727444, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 313.51350940088383368, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 321.57073678022845797, tolerance = tolerance))
     >
     > #polychor
     > res <- mvord:::mvord(formula = MMO(Y) ~ 1,
     + data = df,
     + link = mvprobit(),
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + control= mvord.control(solver="BFGS",se=T))
     Note: First threshold for each response is fixed to 0 in order to ensure identifiability!
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 1, data = df, error.structure = cor_general(~1),
     link = mvprobit(), threshold.constraints = c(1, 1), control = mvord.control(solver = "BFGS",
     se = T))
    
     Formula: MMO(Y) ~ 1
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit fix1first 100 2 -153.56 315.46 326.32 52
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 0.00000000000000000 0.00000000000000000 NA NA
     1 2|3 1.61576075469784497 0.15607113684416782 10.35272 < 2.22e-16 ***
     2 1|2 0.00000000000000000 0.00000000000000000 NA NA
     2 2|3 1.61576075469784497 0.15607113684416782 10.35272 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     (Intercept) 1 0.73786225358709867 0.13500906208906058 5.4652799999999999
     (Intercept) 2 0.77255634837212872 0.14063945843173736 5.4931700000000001
     Pr(>|z|)
     (Intercept) 1 4.6218e-08 ***
     (Intercept) 2 3.9478e-08 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value Pr(>|z|)
     corr 1 2 0.90637393266264266 0.03891605267518549 23.290489999999998 < 2.22e-16
    
     corr 1 2 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(0.0000000000000000000, 1.6157607546978449697, 0.0000000000000000000, 1.6157607546978449697), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.73786225358709867095, 0.77255634837212872057), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.90637393266264265623), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.00000000000000000000, 0.15607113684416776267, 0.00000000000000000000, 0.15607113684416776267), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.13500906208906049422, 0.14063945843173727979), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.038916052675185497378), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -153.5640565887688922, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 315.46144651087109878, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 326.31632228582151356, tolerance = tolerance))
     > #######################################################################
     > ## cor_general(~factor)
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + error.structure = cor_general(~f2),
     + threshold.constraints = c(1,1),
     + coef.constraints = c(1,1),
     + contrasts = list(f2 = "contr.sum"))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cor_general(~f2),
     link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), contrasts = list(f2 = "contr.sum"), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -134.17 283.39 303 38
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -0.90552506790519915 0.18265647646915181 -4.9575300000000002 7.1395e-07
     1 2|3 1.00420547521745429 0.17467164336594765 5.7491000000000003 8.9718e-09
     2 1|2 -0.90552506790519915 0.18265647646915181 -4.9575300000000002 7.1395e-07
     2 2|3 1.00420547521745429 0.17467164336594765 5.7491000000000003 8.9718e-09
    
     1 1|2 ***
     1 2|3 ***
     2 1|2 ***
     2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
     X1 1 0.64793126338493945 0.14153585035881874 4.5778600000000003 4.6976e-06
     X2 1 -0.42893128334048874 0.13989943737309904 -3.0659999999999998 0.0021695
    
     X1 1 ***
     X2 1 **
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr f21 1 2 0.736678597234154964 0.178238394308490727 4.1331100000000003
     corr f22 1 2 0.918298263052344188 0.068858656199954238 13.3359900000000007
     corr f23 1 2 0.837106119712607066 0.131036082998728903 6.3883599999999996
     Pr(>|z|)
     corr f21 1 2 3.5789e-05 ***
     corr f22 1 2 < 2.22e-16 ***
     corr f23 1 2 1.6767e-10 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate,
     + c(-0.90552506790519915469, 1.00420547521745429087, -0.90552506790519915469, 1.00420547521745429087), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate,
     + c(0.64793126338493944871, -0.42893128334048874484), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate,
     + c(0.73667859723415496376, 0.91829826305234418804, 0.83710611971260706632), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`,
     + c(0.18265647646915181279, 0.17467164336594759311, 0.18265647646915181279, 0.17467164336594759311), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`,
     + c(0.14153585035881874332, 0.13989943737309898375), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`,
     + c(0.17823839430849069965, 0.06885865619995423792, 0.13103608299872890330), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -134.17046176709035876, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 283.39468697504094052, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 303.00349482656417877, tolerance = tolerance))
     >
     > ##################################################################################
     > res <- mvord:::mvord(formula = MMO(Y, i, j) ~ 0 + X1 + X2,
     + data = df,
     + link = mvlogit(df = 8L),
     + error.structure = cor_general(~1))
     > #threshold.constraints = c(1,1),
     > #coef.constraints = c(1,1))
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y, i, j) ~ 0 + X1 + X2, data = df,
     error.structure = cor_general(~1), link = mvlogit(df = 8L))
    
     Formula: MMO(Y, i, j) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvlogit flexible 100 2 -132.87 285.51 311.28 562
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -1.61992769101566214 0.33738754131593512 -4.8013899999999996 1.5757e-06
     1 2|3 1.83247577798092465 0.36742756957437139 4.9873099999999999 6.1225e-07
     2 1|2 -1.61590165396073338 0.33963202935670594 -4.7577999999999996 1.9571e-06
     2 2|3 1.74613517856085321 0.30112278942370385 5.7987500000000001 6.6812e-09
    
     1 1|2 ***
     1 2|3 ***
     2 1|2 ***
     2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
     X1 1 1.28924211373739284 0.30395058430878913 4.2416200000000002 2.2191e-05
     X1 2 0.88420369072675176 0.29740837377318347 2.9730300000000001 0.0029488
     X2 1 -0.84291006112195865 0.30080364905710810 -2.8021900000000000 0.0050756
     X2 2 -0.71568747719099213 0.25697669461943834 -2.7850299999999999 0.0053523
    
     X1 1 ***
     X1 2 **
     X2 1 **
     X2 2 **
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value Pr(>|z|)
     corr 1 2 0.883893221236931481 0.071174779980675715 12.41863 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > # mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-1.6170817306633420, 1.7855897338762188, -1.6170817306633420, 1.7855897338762188), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$coefficients$Estimate, c(1.07242200717115987, 1.07242200717115987, -0.76715925377701732, -0.76715925377701732), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.85317690560688531), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.28997405649500174, 0.27427389826231802, 0.28997405649500174, 0.27427389826231802), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.24111402270822993, 0.24111402270822993, 0.24156664773225886, 0.24156664773225886), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.063316529381183581), tolerance = tolerance))
     > # mvord:::check(all.equal(logLik(res), -135.41665313840898, tolerance = tolerance))
     > # mvord:::check(all.equal(AIC(res), 281.35962206629165, tolerance = tolerance))
     > # mvord:::check(all.equal(BIC(res), 295.07104409780789, tolerance = tolerance))
     >
     > ##################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cov_general(~1),
     + threshold.constraints = c(1,1),
     + threshold.values = list(c(-1,NA),
     + c(-1,NA)),
     + coef.constraints = c(1,1))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cov_general(~1),
     link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), threshold.values = list(c(-1, NA), c(-1, NA)), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit fix1first 100 2 -134.64 282.04 298.67 37
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -1.00000000000000000 0.00000000000000000 NA NA
     1 2|3 1.08266798008337450 0.24926043418151714 4.3435199999999998 1.4022e-05
     2 1|2 -1.00000000000000000 0.00000000000000000 NA NA
     2 2|3 1.08266798008337450 0.24926043418151714 4.3435199999999998 1.4022e-05
    
     1 1|2
     1 2|3 ***
     2 1|2
     2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
     X1 1 0.68834946142093179 0.15128036741853496 4.5501600000000000 5.3606e-06
     X2 1 -0.45406517655051332 0.15704134440910930 -2.8913700000000002 0.0038356
    
     X1 1 ***
     X2 1 **
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.856153637803174328 0.061822221469819631 13.8486399999999996
     sigma 1 1.001934583868918871 0.185187960683921543 5.4103700000000003
     sigma 2 1.089828199315766355 0.197593289723851151 5.5155099999999999
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     sigma 1 6.2896e-08 ***
     sigma 2 3.4777e-08 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-1.000000000000000000000, 1.082667980083374503764, -1.000000000000000000000, 1.082667980083374503764), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.6883494614209317852271, -0.4540651765505133163892), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.8561536378031743277361, 1.0019345838689188710191, 1.0898281993157663549709), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.0000000000000000000000, 0.2492604341815170820862, 0.0000000000000000000000, 0.2492604341815170820862), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1512803674185349001036, 0.1570413444091093568833), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.06182222146981963123435, 0.18518796068392143205905, 0.19759328972385112321852), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -134.6391112878561102661, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 282.0441800225207202857, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 298.6729258905298252103, tolerance = tolerance))
     > ########################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_equi(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = cbind(c(1,1),c(1,2)))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cor_equi(~1),
     link = mvprobit(), coef.constraints = cbind(c(1, 1), c(1,
     2)), threshold.constraints = c(1, 1), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -134.84 282.45 299.08 35
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -0.96275665235833874 0.16703865743957075 -5.7636799999999999 8.2301e-09
     1 2|3 1.03388029135981419 0.15203091606548902 6.8004600000000002 1.0429e-11
     2 1|2 -0.96275665235833885 0.16703865743957075 -5.7636799999999999 8.2301e-09
     2 2|3 1.03388029135981419 0.15203091606548902 6.8004600000000002 1.0429e-11
    
     1 1|2 ***
     1 2|3 ***
     2 1|2 ***
     2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
     X1 1 0.63820467410973414 0.13670492809262946 4.6684799999999997 3.0343e-06
     X2 1 -0.44676561119532299 0.15900879662665307 -2.8096899999999998 0.0049589
     X2 2 -0.40750155519121800 0.13850979189699270 -2.9420400000000000 0.0032606
    
     X1 1 ***
     X2 1 **
     X2 2 **
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     (Intercept) 1.27290057125464440 0.23221611847930698 5.4815300000000002
     Pr(>|z|)
     (Intercept) 4.2165e-08 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.9627566523584676350112, 1.0338802913596250032668, -0.9627566523584676350112, 1.0338802913596247812222), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.6382046741095841468905, -0.4467656111955374820255, -0.4075015551913261924177), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$error.structure$Estimate, c(1.272900571254629964457), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.1670386574395647250046, 0.1520309160654948199554, 0.1670386574395647250046, 0.1520309160654948199554), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1367049280926301846328, 0.1590087966266460206555, 0.1385097918969951935608), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.2322161184793045674013), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -134.8432321319771745038, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 282.4524217107628487611, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 299.0811675787719536856, tolerance = tolerance))
     >
     > res2 <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_equi(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = list(X2 = cbind(c(1,1,0,0), c(0,0,1,1)),
     + X1 = matrix(rep(1,4), ncol = 1)))
     >
     > mvord:::check(all.equal(res$beta, res2$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res2$sebeta, tolerance = tolerance))
     >
     > ########################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_ar1(~ 1 + X1),
     + threshold.constraints = c(1,1),
     + coef.constraints = c(1,1))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cor_ar1(~1 +
     X1), link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -133.94 280.64 297.27 35
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -0.95722951152617552 0.16558288168119875 -5.7809699999999999 7.4272e-09
     1 2|3 1.03746793968331619 0.15306952616698108 6.7777599999999998 1.2206e-11
     2 1|2 -0.95722951152617552 0.16558288168119875 -5.7809699999999999 7.4272e-09
     2 2|3 1.03746793968331619 0.15306952616698108 6.7777599999999998 1.2206e-11
    
     1 1|2 ***
     1 2|3 ***
     2 1|2 ***
     2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
     X1 1 0.65304208674282194 0.13584016202004612 4.8074300000000001 1.5288e-06
     X2 1 -0.42273136793544752 0.13820217716271033 -3.0587900000000001 0.0022223
    
     X1 1 ***
     X2 1 **
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     (Intercept) 1.29891554197132031 0.24517466007553268 5.29792000000000041
     X1 0.29350892362932679 0.29416988646865150 0.99775000000000003
     Pr(>|z|)
     (Intercept) 1.1713e-07 ***
     X1 0.3184
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.9572295115261755249492, 1.0374679396833161870717, -0.9572295115261755249492, 1.0374679396833161870717), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.6530420867428219366957, -0.4227313679354475217664), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(1.298915541971320308789, 0.293508923629326790028), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.1655828816811988069002, 0.1530695261669812456962, 0.1655828816811988069002, 0.1530695261669812456962), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1358401620200463122412, 0.1382021771627104100855), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.2451746600755327076815, 0.2941698864686516090572), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -133.938827675498004055, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 280.6436127978045078635, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 297.272358665813612788, tolerance = tolerance))
     >
     > res2 <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_ar1(~1 + X1),
     + threshold.constraints = c(1,1),
     + coef.constraints = list(matrix(rep(1,4), ncol = 1), matrix(rep(1,4),
     + ncol = 1)))
     >
     > mvord:::check(all.equal(res$beta, res2$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res2$sebeta, tolerance = tolerance))
     >
     > ########################################################################################
     > res <- mvord(formula = MMO2(Y1,Y2) ~ 0 + X1 + X2,
     + data = data_toy_example,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = cbind(c(1,2),c(NA,1)))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord(formula = MMO2(Y1, Y2) ~ 0 + X1 + X2, data = data_toy_example,
     error.structure = cor_general(~1), link = mvprobit(), coef.constraints = cbind(c(1,
     2), c(NA, 1)), threshold.constraints = c(1, 1), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO2(Y1, Y2) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -137.65 288.06 304.69 63
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     Y1 1|2 -0.89892942527339637 0.16246077762273251 -5.5332100000000004 3.1442e-08
     Y1 2|3 0.98413365940575759 0.15981657213962935 6.1578900000000001 7.3718e-10
     Y2 1|2 -0.89892942527339637 0.16246077762273251 -5.5332100000000004 3.1442e-08
     Y2 2|3 0.98413365940575770 0.15981657213962935 6.1578900000000001 7.3718e-10
    
     Y1 1|2 ***
     Y1 2|3 ***
     Y2 1|2 ***
     Y2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
     X1 1 0.681944859288474370 0.174961892922601003 3.89767999999999981 9.712e-05
     X1 2 0.468376583108382760 0.157106006693620737 2.98127999999999993 0.0028705
     X2 1 -0.052179801987817860 0.089798215416046726 -0.58108000000000004 0.5611876
    
     X1 1 ***
     X1 2 **
     X2 1
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr Y1 Y2 0.890390657144797815 0.063175533023959646 14.093920000000001
     Pr(>|z|)
     corr Y1 Y2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.8989294252736504953205, 0.9841336594057152886705, -0.8989294252736503842982, 0.9841336594057153996928), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.68194485928859438494953, 0.46837658310856356003171, -0.05217980198794509860694), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.8903906571446290607597), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.1624607776227008359982, 0.1598165721396140603883, 0.1624607776227008359982, 0.1598165721396140603883), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.17496189292260000414103, 0.15710600669362786985239, 0.08979821541610534529898), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.0631755330240382911855), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -137.6494615446075613363, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 288.064880536023622426, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 304.6936264040327273506, tolerance = tolerance))
     >
     >
     > res2 <- mvord:::mvord(formula = MMO2(Y1,Y2) ~ 0 + X1 + X2,
     + data = data_toy_example,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = list(X1 = cbind(c(1,1,0,0), c(0,0,1,1)),
     + X2 = matrix(c(rep(0,2),rep(1,2)), ncol = 1)))
     >
     > mvord:::check(all.equal(res$beta, res2$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res2$sebeta, tolerance = tolerance))
     >
     > ########################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + offset(X2),
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,2),
     + coef.constraints = list(X1 = cbind(c(1,1,0,0), c(0,0,1,1))))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + offset(X2), data = df,
     error.structure = cor_general(~1), link = mvprobit(), coef.constraints = list(X1 = cbind(c(1,
     1, 0, 0), c(0, 0, 1, 1))), threshold.constraints = c(1,
     2), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + offset(X2)
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -194.33 403.71 423.31 78
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -0.76726526230273051 0.12023341892343348 -6.3814599999999997 1.7540e-10
     1 2|3 0.89190483010004851 0.13058942975202986 6.8298399999999999 8.5010e-12
     2 1|2 -0.79675344626870870 0.11172347359918323 -7.1314799999999998 9.9297e-13
     2 2|3 0.87294080476149372 0.13662238116424227 6.3894399999999996 1.6649e-10
    
     1 1|2 ***
     1 2|3 ***
     2 1|2 ***
     2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
     X1 1 0.53370081959015414 0.13636430505433184 3.9137900000000001 9.086e-05 ***
     X1 2 0.34385576067174201 0.10634167794017362 3.2334999999999998 0.0012228 **
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value Pr(>|z|)
     corr 1 2 0.940808527820466201 0.033697206769499194 27.91948 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.7672652623027305107684, 0.8919048301000485068357, -0.7967534462687086982413, 0.8729408047614937160574), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.5337008195901541407480, 0.3438557606717420056519), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.9408085278204662005308), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.1202334189234334493879, 0.1305894297520299174309, 0.1117234735991832761393, 0.1366223811642423013879), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1363643050543318635537, 0.1063416779401736012023), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.03369720676949919391241), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -194.3258468555947047207, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 403.705457152049632441, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 423.314265003572927526, tolerance = tolerance))
     >
     > res2 <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1,
     + offset = list(df$X2[1:100], df$X2[101:200]),
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,2),
     + coef.constraints = list(X1 = cbind(c(1,1,0,0), c(0,0,1,1))))
     >
     > mvord:::check(all.equal(res$beta, res2$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res2$sebeta, tolerance = tolerance))
     >
     > res3 <- mvord:::mvord(formula = MMO2(Y1,Y2) ~ 0 + X1 + offset(X2),
     + data = data_toy_example,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,2),
     + coef.constraints = list(X1 = cbind(c(1,1,0,0), c(0,0,1,1))))
     > mvord:::check(all.equal(res$beta, res3$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res3$sebeta, tolerance = tolerance))
     >
     > ########################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvlogit(),
     + #solver = "newuoa",
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = cbind(c(NA, NA), c(1,2)),
     + coef.values = cbind(c(1, 1), c(NA,NA)))
     > #rho <- res$rho
     >
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cor_general(~1),
     link = mvlogit(), coef.constraints = cbind(c(NA, NA), c(1,
     2)), coef.values = cbind(c(1, 1), c(NA, NA)), threshold.constraints = c(1,
     1))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvlogit flexible 100 2 -135.52 281.57 295.28 157
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -1.58328882997919007 0.25851083789750051 -6.1246499999999999 9.0882e-10
     1 2|3 1.75864247509507332 0.24251592784905907 7.2516600000000002 4.1170e-13
     2 1|2 -1.58328882997919029 0.25851083789750051 -6.1246499999999999 9.0882e-10
     2 2|3 1.75864247509507310 0.24251592784905907 7.2516600000000002 4.1170e-13
    
     1 1|2 ***
     1 2|3 ***
     2 1|2 ***
     2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
     X2 1 -0.77306127413840775 0.26732808332585767 -2.8918100000000000 0.0038303 **
     X2 2 -0.72430349613959055 0.24403985107405990 -2.9679700000000002 0.0029977 **
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.855894394025347838 0.060599214162206316 14.123849999999999
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     >
     > options(digits = 22)
     >
     > tolerance2 <- 1e-4
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-1.583288915548024533564, 1.758642501421955106622, -1.583288915548024533564, 1.758642501421955106622), tolerance = tolerance2))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(-0.7730608423417825170176, -0.7243032675736262859800), tolerance = tolerance2))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.8558944465934688050623), tolerance = tolerance2))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.2585108482987530376107, 0.2425159370850978601819, 0.2585108482987530376107, 0.2425159370850978601819), tolerance = tolerance2))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.2673280928737744588375, 0.2440398419906193439033), tolerance = tolerance2))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.06059919526643173237623), tolerance = tolerance2))
     > mvord:::check(all.equal(logLik(res)[[1]], -135.5210996495467554723, tolerance = tolerance2))
     > mvord:::check(all.equal(AIC(res), 281.568515088567210114, tolerance = tolerance2))
     > mvord:::check(all.equal(BIC(res), 295.2799371200834457341, tolerance = tolerance2))
     >
     > # res1 <- mvord:::mvord(formula = MMO(Y) ~ 1 + X1 + X2,
     > # data = df,
     > # index = c("i", "j"),
     > # link = mvprobit(),
     > # solver = "BFGS",
     > # se = TRUE,
     > # error.structure = cor_equi(~1),
     > # threshold.constraints = c(1,1),
     > # coef.constraints = list(Intercept = matrix(rep(1,6), ncol = 1),
     > # X2 = cbind(c(1,1,1,0,0,0), c(0,0,0,1,1,1)), X1 = matrix(rep(1,6), ncol = 1)))
     >
     >
     > ########################################################################################
     > df$X3 <- cut(df$X2, c(-Inf, -0.2, 0.2, Inf))
     >
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X3,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = T,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = c(1,1))
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X3, data = df, error.structure = cor_general(~1),
     link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X3
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -136.05 284.87 301.5 34
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -1.21270698725061821 0.21131407873814662 -5.7388800000000000 9.5302e-09
     1 2|3 0.75845153901163942 0.18860790313556303 4.0213099999999997 5.7874e-05
     2 1|2 -1.21270698725061821 0.21131407873814662 -5.7388800000000000 9.5302e-09
     2 2|3 0.75845153901163942 0.18860790313556303 4.0213099999999997 5.7874e-05
    
     1 1|2 ***
     1 2|3 ***
     2 1|2 ***
     2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     X1 1 0.582346304092615452 0.133874274950941285 4.34994999999999976
     X3(-0.2,0.2] 1 0.065949564659511767 0.394250784463482329 0.16728000000000001
     X3(0.2, Inf] 1 -0.639867558058318964 0.239596299509496208 -2.67060999999999993
     Pr(>|z|)
     X1 1 1.3617e-05 ***
     X3(-0.2,0.2] 1 0.8671512
     X3(0.2, Inf] 1 0.0075714 **
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value Pr(>|z|)
     corr 1 2 0.857975765661140977 0.062391829886069061 13.75141 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     >
     > options(digits = 22)
     >
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-1.212706987250618206886, 0.758451539011639419563, -1.212706987250618206886, 0.758451539011639419563), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.58234630409261545214150, 0.06594956465951176682871, -0.63986755805831896370961), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.8579757656611409766256), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.2113140787381464491546, 0.1886079031355630819533, 0.2113140787381464491546, 0.1886079031355630819533), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1338742749509412854891, 0.3942507844634823288565, 0.2395962995094961522913), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.06239182988606906121731), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -136.0526219854373835005, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 284.8712014176832667545, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 301.499947285692371679, tolerance = tolerance))
     > ########################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 1 + X1 * X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = T,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + threshold.values = list(c(-1,NA),
     + c(-1,NA)),
     + coef.constraints = c(1,1))
     > res.summary <- summary(res, short = TRUE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 1 + X1 * X2, data = df, error.structure = cor_general(~1),
     link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), threshold.values = list(c(-1, NA), c(-1, NA)), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 1 + X1 * X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit fix1first 100 2 -134.63 282.02 298.65 34
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -1.00000000000000000 0.00000000000000000 NA NA
     1 2|3 1.00269513578732417 0.22781659168150073 4.4013299999999997 1.0759e-05
    
     1 1|2
     1 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     (Intercept) 1 -0.021913582738025385 0.167100624889530508 -0.13114000000000001
     X1 1 0.621621046245967146 0.138540865251019801 4.48690999999999995
     X2 1 -0.427624354284268537 0.140667605271754703 -3.03996000000000022
     X1:X2 1 -0.102567255999721288 0.150381573286269055 -0.68205000000000005
     Pr(>|z|)
     (Intercept) 1 0.8956645
     X1 1 7.2262e-06 ***
     X2 1 0.0023661 **
     X1:X2 1 0.4952094
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value Pr(>|z|)
     corr 1 2 0.853574805224538435 0.063090989221389948 13.52927 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     >
     > options(digits = 22)
     >
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-1.000000000000000000000, 1.002695135787324165477), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(-0.02191358273802538475516, 0.62162104624596714597118, -0.42762435428426853745165, -0.10256725599972128792903), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.8535748052245384354109), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.0000000000000000000000, 0.2278165916815008962271), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1671006248895306467439, 0.1385408652510198845853, 0.1406676052717547309445, 0.1503815732862690546234), tolerance = tolerance))
     > #mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.063097735777156036), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.06309098922138993426056), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -134.6255603757719541136, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 282.0170781983524079806, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 298.6458240663615129051, tolerance = tolerance))
     > ########################################################################################
     > df_NA <- df[-c(1,90:110),]
     >
     >
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df_NA,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = T,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,2)
     + #coef.constraints = list(X1 = cbind(c(1,1,1,0,0,0), c(0,0,0,1,1,1)))
     + )
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df_NA, error.structure = cor_general(~1),
     link = mvprobit(), threshold.constraints = c(1, 2), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 99 2 -119.45 258.71 284.4 51
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -1.02012637737374523 0.19996641229615930 -5.1014900000000001 3.3699e-07
     1 2|3 1.14180364925261624 0.22379090693938217 5.1021000000000001 3.3590e-07
     2 1|2 -0.90662822300335544 0.19059651595707131 -4.7567899999999996 1.9669e-06
     2 2|3 0.99945117269154393 0.17403069740199875 5.7429600000000001 9.3036e-09
    
     1 1|2 ***
     1 2|3 ***
     2 1|2 ***
     2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
     X1 1 0.83733616646601672 0.18898470050450009 4.4307100000000004 9.3924e-06
     X1 2 0.49821075302585316 0.16630384260859424 2.9957900000000000 0.0027374
     X2 1 -0.44740279449439818 0.18763575165142432 -2.3844200000000000 0.0171060
     X2 2 -0.35390389963344654 0.14489624593605005 -2.4424600000000001 0.0145874
    
     X1 1 ***
     X1 2 **
     X2 1 *
     X2 2 *
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.910650307820846217 0.070553566078659979 12.907220000000001
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate,
     + c(-1.0201263773742383911269, 1.1418036492524186176212, -0.9066282230038859024646, 0.9994511726914282467860), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.8373361664658567349306, 0.4982107530258053085248, -0.4474027944949834356692, -0.3539038996339812781500),
     + tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.9106503078209118307029), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.1999664122961394840949, 0.2237909069393826111405, 0.1905965159570501277209, 0.1740306974019930341679), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1889847005044749950198, 0.1663038426085674814647, 0.1876357516513920686840, 0.1448962459360401411335), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.07055356607861483497768), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -119.4526280916837492896, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 258.7052561833675099479, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 284.3969426996999345647, tolerance = tolerance))
     > ##########################################################################################
     > df_23 <- df
     >
     > df_23$Y[which(df_23$Y[1:100] == 3)] <- 2
     >
     >
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df_23,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,2))
     Error in check_response_missings(rho) :
     Model is not identifiable. For at least one dimension, not all response categories are observed.
     Either remove the non-observed category or set constraints on the threshold parameters to ensure identifiability.
     Calls: <Anonymous> -> check_response_missings
     Execution halted
Flavor: r-devel-windows-ix86+x86_64

Version: 1.0.0
Check: running tests for arch ‘x64’
Result: ERROR
     Running 'check_methods.R' [7s]
     Running 'check_mvord_data.R' [3s]
     Running 'check_mvord_errors.R' [3s]
     Running 'check_set_threshold_type.R' [3s]
     Running 'check_toy_example.R' [38s]
     Running 'check_transf_sigmas.R' [3s]
     Running 'check_transf_thresholds.R' [3s]
    Running the tests in 'tests/check_toy_example.R' failed.
    Complete output:
     > library(mvord)
     Loading required package: minqa
     Loading required package: BB
     Loading required package: ucminf
     Loading required package: dfoptim
     > library(MASS)
     >
     >
     > #data(data_toy_example)
     > tolerance <- 1e-6
     >
     > mult.obs <- 2
     > sigma <- matrix(c(1,0.8,0.8,1), ncol = 2)
     > betas <- list(c(0.8,-0.5),
     + c(0.8,-0.5))
     > thresholds <- list(c(-1,1),c(-1,1))
     > nobs <- 100
     > suppressWarnings(RNGversion("3.5.0"))
     > set.seed(2017)
     > errors <- mvrnorm(n = nobs, mu = rep(0, mult.obs), Sigma = sigma)
     >
     > X1 <- rnorm(nobs, 0, 1)
     > X2 <- rnorm(nobs, 0, 1)
     >
     > pred <- cbind(X1, X2)
     >
     > y <- sapply(1:mult.obs, function(j) pred %*% betas[[j]] + errors[, j], simplify = "array")
     > y.ord <- sapply(1:mult.obs, function(j) cut(y[, , j], c(min(y[, , j]) - 1, c(thresholds[[j]]), max(y[, , j]) + 1),
     + labels = FALSE), simplify = "array")
     > predictors.fixed <- lapply(1:mult.obs, function(j) pred)
     > y <- as.data.frame(y.ord)
     >
     > for(i in 1:mult.obs){
     + y[, i] <- factor(y[, i], levels = sort(unique(y[, i])),
     + ordered = TRUE)
     + }
     >
     >
     >
     >
     > data_toy_example <- cbind.data.frame(y, predictors.fixed[[1]])
     > colnames(data_toy_example) <- c("Y1","Y2", "X1", "X2")
     > w <- c(rep(1/20, 20), rep(1/30, 30), rep(1/20, 20), rep(1/30, 30))
     >
     >
     >
     > # convert data_toy_example into long format
     > df <- cbind.data.frame("i" = rep(1:100,2),
     + "j" = rep(1:2,each = 100),
     + "Y" = c(data_toy_example$Y1,data_toy_example$Y2),
     + "X1" = rep(data_toy_example$X1,2),
     + "X2" = rep(data_toy_example$X2,2),
     + "f1" = factor(sample(rep(data_toy_example$Y2,2)), ordered =F),
     + "f2" = factor(rep(data_toy_example$Y1,2), ordered=F),
     + w=rep(w,2))
     >
     >
     >
     > # library(ROI)
     > # ROI_solver <- function(starting.values, objFun, control){
     > # n <- length(starting.values)
     > # op <- OP(objective = F_objective(objFun, n = n),
     > # bounds = V_bound(li = seq_len(n), lb = rep.int(-Inf, n)))
     > # optRes <- ROI_solve(op, solver = "nlminb",
     > # control = c(list(start = starting.values),
     > # control))
     > # list(optpar = optRes$solution,
     > # objective = optRes$objval) # a vector of length equal to number of parameters to optimize
     > # }
     > #
     > #
     > #
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + link = mvprobit(),
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = c(1,1),
     + contrasts = list(f1 = function(x)
     + contr.treatment(nlevels(df$f1), base = 1),
     + f2 = "contr.sum"),
     + control= mvord.control(solver="BFGS",se=T))
     Warning message:
     In mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, link = mvprobit(), :
     Variables f1 and f2 are absent, the contrasts will be ignored.
     > options(digits = 22)
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cor_general(~1),
     link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), contrasts = list(f1 = function(x) contr.treatment(nlevels(df$f1),
     base = 1), f2 = "contr.sum"), control = mvord.control(solver = "BFGS",
     se = T))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -134.91 280.34 294.06 32
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -0.96257386663519673 0.16613952738605439 -5.7937700000000003 6.8825e-09
     1 2|3 1.03347036873223708 0.15004482537617939 6.8877400000000000 5.6684e-12
     2 1|2 -0.96257386663519673 0.16613952738605439 -5.7937700000000003 6.8825e-09
     2 2|3 1.03347036873223708 0.15004482537617939 6.8877400000000000 5.6684e-12
    
     1 1|2 ***
     1 2|3 ***
     2 1|2 ***
     2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
     X1 1 0.63801035062404310 0.13576747301357892 4.6992900000000004 2.6107e-06
     X2 1 -0.42672524816263474 0.13643665802446264 -3.1276400000000000 0.0017621
    
     X1 1 ***
     X2 1 **
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value Pr(>|z|)
     corr 1 2 0.854266728221226845 0.062440889990360744 13.68121 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.96257386663519672876, 1.03347036873223707687, -0.96257386663519672876, 1.03347036873223707687), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.63801035062404309883, -0.42672524816263474046), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.85426672822122684536), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.16613952738605441972, 0.15004482537617935822, 0.16613952738605441972, 0.15004482537617935822), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.13576747301357894315, 0.13643665802446261481), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.062440889990360744222), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -134.90867383086322207, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 280.3436634512001433, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 294.05508548271637892, tolerance = tolerance))
     >
     >
     > res2 <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = list(matrix(rep(1,4), ncol = 1), matrix(rep(1,4), ncol = 1)))
     >
     > mvord:::check(all.equal(res$beta, res2$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res2$sebeta, tolerance = tolerance))
     > ########################################################################
     > ## No coefficients
     > res <- mvord:::mvord(formula = MMO(Y) ~ -1,
     + data = df,
     + link = mvprobit(),
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + control= mvord.control(solver="BFGS",se=T))
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ -1, data = df, error.structure = cor_general(~1),
     link = mvprobit(), threshold.constraints = c(1, 1), control = mvord.control(solver = "BFGS",
     se = T))
    
     Formula: MMO(Y) ~ -1
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -153.66 313.51 321.57 26
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -0.75479706538110092 0.13105603927553963 -5.7593500000000004 8.4440e-09
     1 2|3 0.86086304364935784 0.13360642334434200 6.4432799999999997 1.1692e-10
     2 1|2 -0.75479706538110092 0.13105603927553963 -5.7593500000000004 8.4440e-09
     2 2|3 0.86086304364935784 0.13360642334434200 6.4432799999999997 1.1692e-10
    
     1 1|2 ***
     1 2|3 ***
     2 1|2 ***
     2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
    
     Error Structure:
     Estimate Std. Error z value Pr(>|z|)
     corr 1 2 0.905795171446422409 0.038855435949034602 23.31193 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.75479706538110091785, 0.86086304364935783973, -0.75479706538110091785, 0.86086304364935783973), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(numeric(0)), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.90579517144642240911), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.13105603927553965971, 0.13360642334434202905, 0.13105603927553965971, 0.13360642334434202905), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(numeric(0)), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.038855435949034601573), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -153.66397119528727444, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 313.51350940088383368, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 321.57073678022845797, tolerance = tolerance))
     >
     > #polychor
     > res <- mvord:::mvord(formula = MMO(Y) ~ 1,
     + data = df,
     + link = mvprobit(),
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + control= mvord.control(solver="BFGS",se=T))
     Note: First threshold for each response is fixed to 0 in order to ensure identifiability!
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 1, data = df, error.structure = cor_general(~1),
     link = mvprobit(), threshold.constraints = c(1, 1), control = mvord.control(solver = "BFGS",
     se = T))
    
     Formula: MMO(Y) ~ 1
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit fix1first 100 2 -153.56 315.46 326.32 54
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 0.00000000000000000 0.00000000000000000 NA NA
     1 2|3 1.61576075463128843 0.15607113685033858 10.35272 < 2.22e-16 ***
     2 1|2 0.00000000000000000 0.00000000000000000 NA NA
     2 2|3 1.61576075463128843 0.15607113685033858 10.35272 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     (Intercept) 1 0.73786225415091911 0.13500906210231617 5.4652799999999999
     (Intercept) 2 0.77255634796640460 0.14063945841907260 5.4931700000000001
     Pr(>|z|)
     (Intercept) 1 4.6218e-08 ***
     (Intercept) 2 3.9478e-08 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.906373932680408223 0.038916052671567343 23.290489999999998
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(0.0000000000000000000, 1.6157607546978449697, 0.0000000000000000000, 1.6157607546978449697), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.73786225358709867095, 0.77255634837212872057), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.90637393266264265623), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.00000000000000000000, 0.15607113684416776267, 0.00000000000000000000, 0.15607113684416776267), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.13500906208906049422, 0.14063945843173727979), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.038916052675185497378), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -153.5640565887688922, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 315.46144651087109878, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 326.31632228582151356, tolerance = tolerance))
     > #######################################################################
     > ## cor_general(~factor)
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + error.structure = cor_general(~f2),
     + threshold.constraints = c(1,1),
     + coef.constraints = c(1,1),
     + contrasts = list(f2 = "contr.sum"))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cor_general(~f2),
     link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), contrasts = list(f2 = "contr.sum"), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -134.17 283.39 303 38
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -0.90552506790519915 0.18265647646915181 -4.9575300000000002 7.1395e-07
     1 2|3 1.00420547521745429 0.17467164336594765 5.7491000000000003 8.9718e-09
     2 1|2 -0.90552506790519915 0.18265647646915181 -4.9575300000000002 7.1395e-07
     2 2|3 1.00420547521745429 0.17467164336594765 5.7491000000000003 8.9718e-09
    
     1 1|2 ***
     1 2|3 ***
     2 1|2 ***
     2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
     X1 1 0.64793126338493945 0.14153585035881874 4.5778600000000003 4.6976e-06
     X2 1 -0.42893128334048874 0.13989943737309904 -3.0659999999999998 0.0021695
    
     X1 1 ***
     X2 1 **
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr f21 1 2 0.736678597234154964 0.178238394308490727 4.1331100000000003
     corr f22 1 2 0.918298263052344188 0.068858656199954238 13.3359900000000007
     corr f23 1 2 0.837106119712607066 0.131036082998728903 6.3883599999999996
     Pr(>|z|)
     corr f21 1 2 3.5789e-05 ***
     corr f22 1 2 < 2.22e-16 ***
     corr f23 1 2 1.6767e-10 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate,
     + c(-0.90552506790519915469, 1.00420547521745429087, -0.90552506790519915469, 1.00420547521745429087), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate,
     + c(0.64793126338493944871, -0.42893128334048874484), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate,
     + c(0.73667859723415496376, 0.91829826305234418804, 0.83710611971260706632), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`,
     + c(0.18265647646915181279, 0.17467164336594759311, 0.18265647646915181279, 0.17467164336594759311), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`,
     + c(0.14153585035881874332, 0.13989943737309898375), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`,
     + c(0.17823839430849069965, 0.06885865619995423792, 0.13103608299872890330), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -134.17046176709035876, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 283.39468697504094052, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 303.00349482656417877, tolerance = tolerance))
     >
     > ##################################################################################
     > res <- mvord:::mvord(formula = MMO(Y, i, j) ~ 0 + X1 + X2,
     + data = df,
     + link = mvlogit(df = 8L),
     + error.structure = cor_general(~1))
     > #threshold.constraints = c(1,1),
     > #coef.constraints = c(1,1))
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y, i, j) ~ 0 + X1 + X2, data = df,
     error.structure = cor_general(~1), link = mvlogit(df = 8L))
    
     Formula: MMO(Y, i, j) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvlogit flexible 100 2 -132.87 285.51 311.28 549
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -1.61992766449521630 0.33738755297078843 -4.8013899999999996 1.5757e-06
     1 2|3 1.83247567901219277 0.36742755252813691 4.9873099999999999 6.1225e-07
     2 1|2 -1.61590173471720733 0.33963204977553102 -4.7577999999999996 1.9571e-06
     2 2|3 1.74613527775605459 0.30112278370317763 5.7987500000000001 6.6812e-09
    
     1 1|2 ***
     1 2|3 ***
     2 1|2 ***
     2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
     X1 1 1.28924211459887883 0.30395055981633762 4.2416200000000002 2.2191e-05
     X1 2 0.88420408803769901 0.29740836295750378 2.9730300000000001 0.0029488
     X2 1 -0.84291046430350114 0.30080363775252439 -2.8022000000000000 0.0050756
     X2 2 -0.71568808646347004 0.25697672370480928 -2.7850299999999999 0.0053523
    
     X1 1 ***
     X1 2 **
     X2 1 **
     X2 2 **
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value Pr(>|z|)
     corr 1 2 0.883893409407905395 0.071174682857726018 12.41865 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > # mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-1.6170817306633420, 1.7855897338762188, -1.6170817306633420, 1.7855897338762188), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$coefficients$Estimate, c(1.07242200717115987, 1.07242200717115987, -0.76715925377701732, -0.76715925377701732), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.85317690560688531), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.28997405649500174, 0.27427389826231802, 0.28997405649500174, 0.27427389826231802), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.24111402270822993, 0.24111402270822993, 0.24156664773225886, 0.24156664773225886), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.063316529381183581), tolerance = tolerance))
     > # mvord:::check(all.equal(logLik(res), -135.41665313840898, tolerance = tolerance))
     > # mvord:::check(all.equal(AIC(res), 281.35962206629165, tolerance = tolerance))
     > # mvord:::check(all.equal(BIC(res), 295.07104409780789, tolerance = tolerance))
     >
     > ##################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cov_general(~1),
     + threshold.constraints = c(1,1),
     + threshold.values = list(c(-1,NA),
     + c(-1,NA)),
     + coef.constraints = c(1,1))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cov_general(~1),
     link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), threshold.values = list(c(-1, NA), c(-1, NA)), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit fix1first 100 2 -134.64 282.04 298.67 37
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -1.00000000000000000 0.00000000000000000 NA NA
     1 2|3 1.08266798008337450 0.24926043418151714 4.3435199999999998 1.4022e-05
     2 1|2 -1.00000000000000000 0.00000000000000000 NA NA
     2 2|3 1.08266798008337450 0.24926043418151714 4.3435199999999998 1.4022e-05
    
     1 1|2
     1 2|3 ***
     2 1|2
     2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
     X1 1 0.68834946142093179 0.15128036741853496 4.5501600000000000 5.3606e-06
     X2 1 -0.45406517655051332 0.15704134440910930 -2.8913700000000002 0.0038356
    
     X1 1 ***
     X2 1 **
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.856153637803174328 0.061822221469819631 13.8486399999999996
     sigma 1 1.001934583868918871 0.185187960683921543 5.4103700000000003
     sigma 2 1.089828199315766355 0.197593289723851151 5.5155099999999999
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     sigma 1 6.2896e-08 ***
     sigma 2 3.4777e-08 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-1.000000000000000000000, 1.082667980083374503764, -1.000000000000000000000, 1.082667980083374503764), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.6883494614209317852271, -0.4540651765505133163892), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.8561536378031743277361, 1.0019345838689188710191, 1.0898281993157663549709), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.0000000000000000000000, 0.2492604341815170820862, 0.0000000000000000000000, 0.2492604341815170820862), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1512803674185349001036, 0.1570413444091093568833), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.06182222146981963123435, 0.18518796068392143205905, 0.19759328972385112321852), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -134.6391112878561102661, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 282.0441800225207202857, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 298.6729258905298252103, tolerance = tolerance))
     > ########################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_equi(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = cbind(c(1,1),c(1,2)))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cor_equi(~1),
     link = mvprobit(), coef.constraints = cbind(c(1, 1), c(1,
     2)), threshold.constraints = c(1, 1), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -134.84 282.45 299.08 35
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -0.96275665235833874 0.16703865743957083 -5.7636799999999999 8.2301e-09
     1 2|3 1.03388029135981419 0.15203091606548905 6.8004600000000002 1.0429e-11
     2 1|2 -0.96275665235833885 0.16703865743957083 -5.7636799999999999 8.2301e-09
     2 2|3 1.03388029135981419 0.15203091606548905 6.8004600000000002 1.0429e-11
    
     1 1|2 ***
     1 2|3 ***
     2 1|2 ***
     2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
     X1 1 0.63820467410973414 0.13670492809262952 4.6684799999999997 3.0343e-06
     X2 1 -0.44676561119532299 0.15900879662665299 -2.8096899999999998 0.0049589
     X2 2 -0.40750155519121800 0.13850979189699267 -2.9420400000000000 0.0032606
    
     X1 1 ***
     X2 1 **
     X2 2 **
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     (Intercept) 1.27290057125464440 0.23221611847930684 5.4815300000000002
     Pr(>|z|)
     (Intercept) 4.2165e-08 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.9627566523584676350112, 1.0338802913596250032668, -0.9627566523584676350112, 1.0338802913596247812222), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.6382046741095841468905, -0.4467656111955374820255, -0.4075015551913261924177), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$error.structure$Estimate, c(1.272900571254629964457), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.1670386574395647250046, 0.1520309160654948199554, 0.1670386574395647250046, 0.1520309160654948199554), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1367049280926301846328, 0.1590087966266460206555, 0.1385097918969951935608), tolerance = tolerance))
     > # mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.2322161184793045674013), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -134.8432321319771745038, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 282.4524217107628487611, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 299.0811675787719536856, tolerance = tolerance))
     >
     > res2 <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_equi(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = list(X2 = cbind(c(1,1,0,0), c(0,0,1,1)),
     + X1 = matrix(rep(1,4), ncol = 1)))
     >
     > mvord:::check(all.equal(res$beta, res2$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res2$sebeta, tolerance = tolerance))
     >
     > ########################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_ar1(~ 1 + X1),
     + threshold.constraints = c(1,1),
     + coef.constraints = c(1,1))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cor_ar1(~1 +
     X1), link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -133.94 280.64 297.27 35
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -0.95722951152632973 0.16558288168117641 -5.7809699999999999 7.4272e-09
     1 2|3 1.03746793968359041 0.15306952616694816 6.7777599999999998 1.2206e-11
     2 1|2 -0.95722951152632973 0.16558288168117641 -5.7809699999999999 7.4272e-09
     2 2|3 1.03746793968359041 0.15306952616694816 6.7777599999999998 1.2206e-11
    
     1 1|2 ***
     1 2|3 ***
     2 1|2 ***
     2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
     X1 1 0.65304208674305497 0.13584016202001462 4.8074300000000001 1.5288e-06
     X2 1 -0.42273136793561233 0.13820217716268313 -3.0587900000000001 0.0022223
    
     X1 1 ***
     X2 1 **
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     (Intercept) 1.29891554197068548 0.24517466007549649 5.29792000000000041
     X1 0.29350892362908221 0.29416988646861303 0.99775000000000003
     Pr(>|z|)
     (Intercept) 1.1713e-07 ***
     X1 0.3184
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.9572295115261755249492, 1.0374679396833161870717, -0.9572295115261755249492, 1.0374679396833161870717), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.6530420867428219366957, -0.4227313679354475217664), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(1.298915541971320308789, 0.293508923629326790028), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.1655828816811988069002, 0.1530695261669812456962, 0.1655828816811988069002, 0.1530695261669812456962), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1358401620200463122412, 0.1382021771627104100855), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.2451746600755327076815, 0.2941698864686516090572), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -133.938827675498004055, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 280.6436127978045078635, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 297.272358665813612788, tolerance = tolerance))
     >
     > res2 <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_ar1(~1 + X1),
     + threshold.constraints = c(1,1),
     + coef.constraints = list(matrix(rep(1,4), ncol = 1), matrix(rep(1,4),
     + ncol = 1)))
     >
     > mvord:::check(all.equal(res$beta, res2$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res2$sebeta, tolerance = tolerance))
     >
     > ########################################################################################
     > res <- mvord(formula = MMO2(Y1,Y2) ~ 0 + X1 + X2,
     + data = data_toy_example,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = cbind(c(1,2),c(NA,1)))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord(formula = MMO2(Y1, Y2) ~ 0 + X1 + X2, data = data_toy_example,
     error.structure = cor_general(~1), link = mvprobit(), coef.constraints = cbind(c(1,
     2), c(NA, 1)), threshold.constraints = c(1, 1), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO2(Y1, Y2) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -137.65 288.06 304.69 63
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     Y1 1|2 -0.89892942527341468 0.16246077762271194 -5.5332100000000004 3.1442e-08
     Y1 2|3 0.98413365940589015 0.15981657213960945 6.1578900000000001 7.3718e-10
     Y2 1|2 -0.89892942527341468 0.16246077762271194 -5.5332100000000004 3.1442e-08
     Y2 2|3 0.98413365940589004 0.15981657213960945 6.1578900000000001 7.3718e-10
    
     Y1 1|2 ***
     Y1 2|3 ***
     Y2 1|2 ***
     Y2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
     X1 1 0.681944859288810656 0.174961892922594314 3.89767999999999981 9.712e-05
     X1 2 0.468376583108646438 0.157106006693618488 2.98127999999999993 0.0028705
     X2 1 -0.052179801987755611 0.089798215416072275 -0.58108000000000004 0.5611876
    
     X1 1 ***
     X1 2 **
     X2 1
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr Y1 Y2 0.890390657144677355 0.063175533024017544 14.093920000000001
     Pr(>|z|)
     corr Y1 Y2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.8989294252736504953205, 0.9841336594057152886705, -0.8989294252736503842982, 0.9841336594057153996928), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.68194485928859438494953, 0.46837658310856356003171, -0.05217980198794509860694), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.8903906571446290607597), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.1624607776227008359982, 0.1598165721396140603883, 0.1624607776227008359982, 0.1598165721396140603883), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.17496189292260000414103, 0.15710600669362786985239, 0.08979821541610534529898), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.0631755330240382911855), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -137.6494615446075613363, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 288.064880536023622426, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 304.6936264040327273506, tolerance = tolerance))
     >
     >
     > res2 <- mvord:::mvord(formula = MMO2(Y1,Y2) ~ 0 + X1 + X2,
     + data = data_toy_example,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = list(X1 = cbind(c(1,1,0,0), c(0,0,1,1)),
     + X2 = matrix(c(rep(0,2),rep(1,2)), ncol = 1)))
     >
     > mvord:::check(all.equal(res$beta, res2$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res2$sebeta, tolerance = tolerance))
     >
     > ########################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + offset(X2),
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,2),
     + coef.constraints = list(X1 = cbind(c(1,1,0,0), c(0,0,1,1))))
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + offset(X2), data = df,
     error.structure = cor_general(~1), link = mvprobit(), coef.constraints = list(X1 = cbind(c(1,
     1, 0, 0), c(0, 0, 1, 1))), threshold.constraints = c(1,
     2), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + offset(X2)
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -194.33 403.71 423.31 78
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -0.76726526230273051 0.12023341892343348 -6.3814599999999997 1.7540e-10
     1 2|3 0.89190483010004851 0.13058942975202986 6.8298399999999999 8.5010e-12
     2 1|2 -0.79675344626870870 0.11172347359918323 -7.1314799999999998 9.9297e-13
     2 2|3 0.87294080476149372 0.13662238116424227 6.3894399999999996 1.6649e-10
    
     1 1|2 ***
     1 2|3 ***
     2 1|2 ***
     2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
     X1 1 0.53370081959015414 0.13636430505433184 3.9137900000000001 9.086e-05 ***
     X1 2 0.34385576067174201 0.10634167794017362 3.2334999999999998 0.0012228 **
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value Pr(>|z|)
     corr 1 2 0.940808527820466201 0.033697206769499194 27.91948 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-0.7672652623027305107684, 0.8919048301000485068357, -0.7967534462687086982413, 0.8729408047614937160574), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.5337008195901541407480, 0.3438557606717420056519), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.9408085278204662005308), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.1202334189234334493879, 0.1305894297520299174309, 0.1117234735991832761393, 0.1366223811642423013879), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1363643050543318635537, 0.1063416779401736012023), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.03369720676949919391241), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -194.3258468555947047207, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 403.705457152049632441, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 423.314265003572927526, tolerance = tolerance))
     >
     > res2 <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1,
     + offset = list(df$X2[1:100], df$X2[101:200]),
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,2),
     + coef.constraints = list(X1 = cbind(c(1,1,0,0), c(0,0,1,1))))
     >
     > mvord:::check(all.equal(res$beta, res2$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res2$sebeta, tolerance = tolerance))
     >
     > res3 <- mvord:::mvord(formula = MMO2(Y1,Y2) ~ 0 + X1 + offset(X2),
     + data = data_toy_example,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,2),
     + coef.constraints = list(X1 = cbind(c(1,1,0,0), c(0,0,1,1))))
     > mvord:::check(all.equal(res$beta, res3$beta, tolerance = tolerance))
     > mvord:::check(all.equal(res$sebeta, res3$sebeta, tolerance = tolerance))
     >
     > ########################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvlogit(),
     + #solver = "newuoa",
     + #se = TRUE,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = cbind(c(NA, NA), c(1,2)),
     + coef.values = cbind(c(1, 1), c(NA,NA)))
     > #rho <- res$rho
     >
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df, error.structure = cor_general(~1),
     link = mvlogit(), coef.constraints = cbind(c(NA, NA), c(1,
     2)), coef.values = cbind(c(1, 1), c(NA, NA)), threshold.constraints = c(1,
     1))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvlogit flexible 100 2 -135.52 281.57 295.28 155
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -1.58328863163709155 0.25851084740245189 -6.1246499999999999 9.0883e-10
     1 2|3 1.75864285951527011 0.24251593328227383 7.2516600000000002 4.1170e-13
     2 1|2 -1.58328863163709177 0.25851084740245189 -6.1246499999999999 9.0883e-10
     2 2|3 1.75864285951526989 0.24251593328227383 7.2516600000000002 4.1170e-13
    
     1 1|2 ***
     1 2|3 ***
     2 1|2 ***
     2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
     X2 1 -0.77306116071759878 0.26732807943954601 -2.8918100000000000 0.0038303 **
     X2 2 -0.72430321202027748 0.24403984758754560 -2.9679700000000002 0.0029977 **
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.855894431283549140 0.060599200129253024 14.123860000000001
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     >
     > options(digits = 22)
     >
     > tolerance2 <- 1e-4
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-1.583288915548024533564, 1.758642501421955106622, -1.583288915548024533564, 1.758642501421955106622), tolerance = tolerance2))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(-0.7730608423417825170176, -0.7243032675736262859800), tolerance = tolerance2))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.8558944465934688050623), tolerance = tolerance2))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.2585108482987530376107, 0.2425159370850978601819, 0.2585108482987530376107, 0.2425159370850978601819), tolerance = tolerance2))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.2673280928737744588375, 0.2440398419906193439033), tolerance = tolerance2))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.06059919526643173237623), tolerance = tolerance2))
     > mvord:::check(all.equal(logLik(res)[[1]], -135.5210996495467554723, tolerance = tolerance2))
     > mvord:::check(all.equal(AIC(res), 281.568515088567210114, tolerance = tolerance2))
     > mvord:::check(all.equal(BIC(res), 295.2799371200834457341, tolerance = tolerance2))
     >
     > # res1 <- mvord:::mvord(formula = MMO(Y) ~ 1 + X1 + X2,
     > # data = df,
     > # index = c("i", "j"),
     > # link = mvprobit(),
     > # solver = "BFGS",
     > # se = TRUE,
     > # error.structure = cor_equi(~1),
     > # threshold.constraints = c(1,1),
     > # coef.constraints = list(Intercept = matrix(rep(1,6), ncol = 1),
     > # X2 = cbind(c(1,1,1,0,0,0), c(0,0,0,1,1,1)), X1 = matrix(rep(1,6), ncol = 1)))
     >
     >
     > ########################################################################################
     > df$X3 <- cut(df$X2, c(-Inf, -0.2, 0.2, Inf))
     >
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X3,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = T,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + coef.constraints = c(1,1))
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X3, data = df, error.structure = cor_general(~1),
     link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X3
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 100 2 -136.05 284.87 301.5 34
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -1.21270698725061821 0.21131407873814662 -5.7388800000000000 9.5302e-09
     1 2|3 0.75845153901163942 0.18860790313556303 4.0213099999999997 5.7874e-05
     2 1|2 -1.21270698725061821 0.21131407873814662 -5.7388800000000000 9.5302e-09
     2 2|3 0.75845153901163942 0.18860790313556303 4.0213099999999997 5.7874e-05
    
     1 1|2 ***
     1 2|3 ***
     2 1|2 ***
     2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     X1 1 0.582346304092615452 0.133874274950941285 4.34994999999999976
     X3(-0.2,0.2] 1 0.065949564659511767 0.394250784463482329 0.16728000000000001
     X3(0.2, Inf] 1 -0.639867558058318964 0.239596299509496208 -2.67060999999999993
     Pr(>|z|)
     X1 1 1.3617e-05 ***
     X3(-0.2,0.2] 1 0.8671512
     X3(0.2, Inf] 1 0.0075714 **
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value Pr(>|z|)
     corr 1 2 0.857975765661140977 0.062391829886069061 13.75141 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     >
     > options(digits = 22)
     >
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-1.212706987250618206886, 0.758451539011639419563, -1.212706987250618206886, 0.758451539011639419563), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.58234630409261545214150, 0.06594956465951176682871, -0.63986755805831896370961), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.8579757656611409766256), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.2113140787381464491546, 0.1886079031355630819533, 0.2113140787381464491546, 0.1886079031355630819533), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1338742749509412854891, 0.3942507844634823288565, 0.2395962995094961522913), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.06239182988606906121731), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -136.0526219854373835005, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 284.8712014176832667545, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 301.499947285692371679, tolerance = tolerance))
     > ########################################################################################
     > res <- mvord:::mvord(formula = MMO(Y) ~ 1 + X1 * X2,
     + data = df,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = T,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,1),
     + threshold.values = list(c(-1,NA),
     + c(-1,NA)),
     + coef.constraints = c(1,1))
     > res.summary <- summary(res, short = TRUE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 1 + X1 * X2, data = df, error.structure = cor_general(~1),
     link = mvprobit(), coef.constraints = c(1, 1), threshold.constraints = c(1,
     1), threshold.values = list(c(-1, NA), c(-1, NA)), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 1 + X1 * X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit fix1first 100 2 -134.63 282.02 298.65 34
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -1.00000000000000000 0.00000000000000000 NA NA
     1 2|3 1.00269513578732417 0.22781659168150073 4.4013299999999997 1.0759e-05
    
     1 1|2
     1 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value
     (Intercept) 1 -0.021913582738025385 0.167100624889530508 -0.13114000000000001
     X1 1 0.621621046245967146 0.138540865251019801 4.48690999999999995
     X2 1 -0.427624354284268537 0.140667605271754703 -3.03996000000000022
     X1:X2 1 -0.102567255999721288 0.150381573286269055 -0.68205000000000005
     Pr(>|z|)
     (Intercept) 1 0.8956645
     X1 1 7.2262e-06 ***
     X2 1 0.0023661 **
     X1:X2 1 0.4952094
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value Pr(>|z|)
     corr 1 2 0.853574805224538435 0.063090989221389948 13.52927 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     >
     > options(digits = 22)
     >
     > mvord:::check(all.equal(res.summary$thresholds$Estimate, c(-1.000000000000000000000, 1.002695135787324165477), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(-0.02191358273802538475516, 0.62162104624596714597118, -0.42762435428426853745165, -0.10256725599972128792903), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.8535748052245384354109), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.0000000000000000000000, 0.2278165916815008962271), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1671006248895306467439, 0.1385408652510198845853, 0.1406676052717547309445, 0.1503815732862690546234), tolerance = tolerance))
     > #mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.063097735777156036), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.06309098922138993426056), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -134.6255603757719541136, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 282.0170781983524079806, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 298.6458240663615129051, tolerance = tolerance))
     > ########################################################################################
     > df_NA <- df[-c(1,90:110),]
     >
     >
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df_NA,
     + #index = c("i", "j"),
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + #se = T,
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,2)
     + #coef.constraints = list(X1 = cbind(c(1,1,1,0,0,0), c(0,0,0,1,1,1)))
     + )
     >
     > res.summary <- summary(res, short = FALSE)
    
     Call: mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2, data = df_NA, error.structure = cor_general(~1),
     link = mvprobit(), threshold.constraints = c(1, 2), control = mvord.control(solver = "BFGS"))
    
     Formula: MMO(Y) ~ 0 + X1 + X2
    
     link threshold nsubjects ndim logPL CLAIC CLBIC fevals
     mvprobit flexible 99 2 -119.45 258.71 284.4 51
    
     Thresholds:
     Estimate Std. Error z value Pr(>|z|)
     1 1|2 -1.02012637737458167 0.19996641229612080 -5.1014900000000001 3.3699e-07
     1 2|3 1.14180364925258937 0.22379090693936640 5.1021000000000001 3.3590e-07
     2 1|2 -0.90662822300392831 0.19059651595702587 -4.7567899999999996 1.9669e-06
     2 2|3 0.99945117269143702 0.17403069740199306 5.7429600000000001 9.3036e-09
    
     1 1|2 ***
     1 2|3 ***
     2 1|2 ***
     2 2|3 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Coefficients:
     Estimate Std. Error z value Pr(>|z|)
     X1 1 0.83733616646652753 0.18898470050441929 4.4307100000000004 9.3924e-06
     X1 2 0.49821075302610573 0.16630384260854328 2.9957900000000000 0.0027374
     X2 1 -0.44740279449461673 0.18763575165138335 -2.3844200000000000 0.0171060
     X2 2 -0.35390389963359220 0.14489624593603051 -2.4424600000000001 0.0145874
    
     X1 1 ***
     X1 2 **
     X2 1 *
     X2 2 *
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
    
     Error Structure:
     Estimate Std. Error z value
     corr 1 2 0.910650307821035288 0.070553566078537591 12.907220000000001
     Pr(>|z|)
     corr 1 2 < 2.22e-16 ***
     ---
     Signif. codes:
     0 '***' 0.001 '**' 0.01 '*' 0.050000000000000003 '.' 0.10000000000000001 ' ' 1
     >
     > options(digits = 22)
     >
     > # paste(format(res.summary$thresholds$Estimate), collapse = ",")
     > # paste(format(res.summary$coefficients$Estimate), collapse = ",")
     > # paste(format(res.summary$error.structure$Estimate), collapse = ",")
     > mvord:::check(all.equal(res.summary$thresholds$Estimate,
     + c(-1.0201263773742383911269, 1.1418036492524186176212, -0.9066282230038859024646, 0.9994511726914282467860), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$Estimate, c(0.8373361664658567349306, 0.4982107530258053085248, -0.4474027944949834356692, -0.3539038996339812781500),
     + tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$Estimate, c(0.9106503078209118307029), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$thresholds$`Std. Error`, c(0.1999664122961394840949, 0.2237909069393826111405, 0.1905965159570501277209, 0.1740306974019930341679), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$coefficients$`Std. Error`, c(0.1889847005044749950198, 0.1663038426085674814647, 0.1876357516513920686840, 0.1448962459360401411335), tolerance = tolerance))
     > mvord:::check(all.equal(res.summary$error.structure$`Std. Error`, c(0.07055356607861483497768), tolerance = tolerance))
     > mvord:::check(all.equal(logLik(res)[[1]], -119.4526280916837492896, tolerance = tolerance))
     > mvord:::check(all.equal(AIC(res), 258.7052561833675099479, tolerance = tolerance))
     > mvord:::check(all.equal(BIC(res), 284.3969426996999345647, tolerance = tolerance))
     > ##########################################################################################
     > df_23 <- df
     >
     > df_23$Y[which(df_23$Y[1:100] == 3)] <- 2
     >
     >
     > res <- mvord:::mvord(formula = MMO(Y) ~ 0 + X1 + X2,
     + data = df_23,
     + link = mvprobit(),
     + control = mvord.control(solver = "BFGS"),
     + error.structure = cor_general(~1),
     + threshold.constraints = c(1,2))
     Error in check_response_missings(rho) :
     Model is not identifiable. For at least one dimension, not all response categories are observed.
     Either remove the non-observed category or set constraints on the threshold parameters to ensure identifiability.
     Calls: <Anonymous> -> check_response_missings
     Execution halted
Flavor: r-devel-windows-ix86+x86_64

Version: 1.0.0
Check: re-building of vignette outputs
Result: NOTE
    Error(s) in re-building vignettes:
    --- re-building ‘vignette_mvord2.Rmd’ using rmarkdown
    Warning in engine$weave(file, quiet = quiet, encoding = enc) :
     Pandoc (>= 1.12.3) and/or pandoc-citeproc not available. Falling back to R Markdown v1.
    Loading required package: minqa
    Loading required package: BB
    Loading required package: ucminf
    Loading required package: dfoptim
    --- finished re-building ‘vignette_mvord2.Rmd’
    
    --- re-building ‘vignette_mvord.Rnw’ using Sweave
    Loading required package: minqa
    Loading required package: BB
    Loading required package: ucminf
    Loading required package: dfoptim
    Error: processing vignette 'vignette_mvord.Rnw' failed with diagnostics:
    Running 'texi2dvi' on 'vignette_mvord.tex' failed.
    LaTeX errors:
    !pdfTeX error: pdflatex (file bbm10): Font bbm10 at 657 not found
     ==> Fatal error occurred, no output PDF file produced!
    --- failed re-building 'vignette_mvord.Rnw'
    
    SUMMARY: processing the following file failed:
     ‘vignette_mvord.Rnw’
    
    Error: Vignette re-building failed.
    Execution halted
Flavor: r-oldrel-osx-x86_64