tests/testthat/test-smk-ds.gamlss.R

In datashield/dsBaseClient: DataSHIELD Client Functions

#-------------------------------------------------------------------------------
# Copyright (c) 2019-2022 University of Newcastle upon Tyne. All rights reserved.
#
# This program and the accompanying materials
# are made available under the terms of the GNU Public License v3.0.
#
# You should have received a copy of the GNU General Public License
# along with this program.  If not, see <http://www.gnu.org/licenses/>.
#-------------------------------------------------------------------------------

#
# Set up
#

context("ds.gamlss::smk::setup")

connect.studies.dataset.gamlss(list("e3_bw", "e3_gac_None", "hs_zbmi_who", "hs_child_age_None",
                                    "h_mbmi_None", "hs_correct_raven", "hs_wgtgain_None"))

test_that("setup", {
    ds_expect_variables(c("D"))
})

#
# Tests
#

context("ds.gamlss::smk::mod1")
test_that("gamlss, birth weight", {
    model_e3_bw.DS <- ds.gamlss(formula = 'e3_bw ~ e3_gac_None',
                              sigma.formula = 'e3_bw ~ e3_gac_None',
                              data = 'D', family = 'NO()', centiles = TRUE,
                              xvar = 'D$e3_gac_None', newobj = 'z_scores_e3_bw')

    expect_length(model_e3_bw.DS, 3)
    expect_length(model_e3_bw.DS$study1, 2)
    expect_length(model_e3_bw.DS$study2, 2)
    expect_length(model_e3_bw.DS$study3, 2)
    expect_equal(model_e3_bw.DS$study1$results$family, c("NO","Normal"))
    expect_true("character" %in% class(model_e3_bw.DS$study1$results$family))
    expect_equal(model_e3_bw.DS$study1$results$parameters, c("mu","sigma"))
    expect_true("character" %in% class(model_e3_bw.DS$study1$results$parameters))
    expect_true("call" %in% class(model_e3_bw.DS$study1$results$call))
    expect_equal(model_e3_bw.DS$study1$results$y, "The response variable is not disclosed!")
    expect_true("character" %in% class(model_e3_bw.DS$study1$results$y))
    expect_true("list" %in% class(model_e3_bw.DS$study1$results$control))
    expect_equal(model_e3_bw.DS$study1$results$control$c.crit, 0.001)
    expect_equal(model_e3_bw.DS$study1$results$control$n.cyc, 20)
    expect_equal(model_e3_bw.DS$study1$results$control$mu.step, 1)
    expect_equal(model_e3_bw.DS$study1$results$control$sigma.step, 1)
    expect_equal(model_e3_bw.DS$study1$results$control$nu.step, 1)
    expect_equal(model_e3_bw.DS$study1$results$control$tau.step, 1)
    expect_equal(model_e3_bw.DS$study1$results$control$gd.tol, Inf)
    expect_equal(model_e3_bw.DS$study1$results$control$iter, 0)
    expect_equal(model_e3_bw.DS$study1$results$control$trace, TRUE)
    expect_equal(model_e3_bw.DS$study1$results$control$autostep, TRUE)
    expect_equal(model_e3_bw.DS$study1$results$control$save, TRUE)
    expect_true("numeric" %in% class(model_e3_bw.DS$study1$results$weights))
    expect_equal(model_e3_bw.DS$study1$results$G.deviance, 9071.036, tolerance=1e-07)
    expect_equal(model_e3_bw.DS$study1$results$N, 607)
    expect_true("expression" %in% class(model_e3_bw.DS$study1$results$rqres))
    expect_equal(model_e3_bw.DS$study1$results$iter, 4)
    expect_equal(model_e3_bw.DS$study1$results$type, "Continuous")
    expect_true("character" %in% class(model_e3_bw.DS$study1$results$type))
    expect_true("call" %in% class(model_e3_bw.DS$study1$results$method))
    expect_equal(model_e3_bw.DS$study1$results$contrasts, NULL)
    expect_true("NULL" %in% class(model_e3_bw.DS$study1$results$contrasts))
    expect_equal(model_e3_bw.DS$study1$results$converged, TRUE)
    expect_true("logical" %in% class(model_e3_bw.DS$study1$results$converged))
    expect_equal(model_e3_bw.DS$study1$results$residuals, "The residuals of the model are not disclosed!")
    expect_true("character" %in% class(model_e3_bw.DS$study1$results$residuals))
    expect_equal(model_e3_bw.DS$study1$results$noObs, 607)
    expect_equal(model_e3_bw.DS$study1$results$mu.fv, "The fitted values of the mu model are not disclosed!")
    expect_equal(model_e3_bw.DS$study1$results$mu.link, "identity")
    expect_true("character" %in% class(model_e3_bw.DS$study1$results$mu.link))
    expect_true("terms" %in% class(model_e3_bw.DS$study1$results$mu.terms))
    expect_true("formula" %in% class(model_e3_bw.DS$study1$results$mu.terms))
    expect_equal(as.numeric(model_e3_bw.DS$study1$results$mu.coefficients["e3_gac_None"]), 178.6317, tolerance=1e-5)
    expect_equal(as.numeric(model_e3_bw.DS$study1$results$mu.coefficients["(Intercept)"]), -3641.618, tolerance=1e-5)
    expect_true("numeric" %in% class(model_e3_bw.DS$study1$results$mu.offset))
    expect_true("formula" %in% class(model_e3_bw.DS$study1$results$mu.formula))
    expect_true("list" %in% class(model_e3_bw.DS$study1$results$mu.xlevels))
    expect_equal(model_e3_bw.DS$study1$results$mu.df, 2)
    expect_equal(model_e3_bw.DS$study1$results$mu.nl.df, 0)
    expect_equal(model_e3_bw.DS$study1$results$mu.pen, 0)
    expect_equal(model_e3_bw.DS$study1$results$P.deviance, 9071.036, tolerance=1e-5)
    expect_equal(model_e3_bw.DS$study1$results$df.fit, 4)
    expect_equal(model_e3_bw.DS$study1$results$pen, 0)
    expect_equal(model_e3_bw.DS$study1$results$df.residual, 603)
    expect_equal(model_e3_bw.DS$study1$results$aic, 9079.036, tolerance=1e-5)
    expect_equal(model_e3_bw.DS$study1$results$sbc, 9096.671, tolerance=1e-5)
    expect_true("matrix" %in% class(model_e3_bw.DS$study1$centiles))
    expect_true("array" %in% class(model_e3_bw.DS$study1$centiles))
    expect_equal(as.numeric(model_e3_bw.DS$study1$centiles[1,'cent']), 0.4)
    expect_equal(as.numeric(model_e3_bw.DS$study1$centiles[2,'cent']), 2.0)
    expect_equal(as.numeric(model_e3_bw.DS$study1$centiles[3,'cent']), 10.0)
    expect_equal(as.numeric(model_e3_bw.DS$study1$centiles[4,'cent']), 25.0)
    expect_equal(as.numeric(model_e3_bw.DS$study1$centiles[5,'cent']), 50.0)
    expect_equal(as.numeric(model_e3_bw.DS$study1$centiles[6,'cent']), 75.0)
    expect_equal(as.numeric(model_e3_bw.DS$study1$centiles[7,'cent']), 90.0)
    expect_equal(as.numeric(model_e3_bw.DS$study1$centiles[8,'cent']), 98.0)
    expect_equal(as.numeric(model_e3_bw.DS$study1$centiles[9,'cent']), 99.6)
    expect_equal(as.numeric(model_e3_bw.DS$study1$centiles[1,'per']), 0.4942339, tolerance=1e-7)
    expect_equal(as.numeric(model_e3_bw.DS$study1$centiles[2,'per']), 2.1416804, tolerance=1e-7)
    expect_equal(as.numeric(model_e3_bw.DS$study1$centiles[3,'per']), 8.4019769, tolerance=1e-7)
    expect_equal(as.numeric(model_e3_bw.DS$study1$centiles[4,'per']), 24.3822076, tolerance=1e-7)
    expect_equal(as.numeric(model_e3_bw.DS$study1$centiles[5,'per']), 50.9060956, tolerance=1e-7)
    expect_equal(as.numeric(model_e3_bw.DS$study1$centiles[6,'per']), 75.4530478, tolerance=1e-7)
    expect_equal(as.numeric(model_e3_bw.DS$study1$centiles[7,'per']), 90.6095552, tolerance=1e-7)
    expect_equal(as.numeric(model_e3_bw.DS$study1$centiles[8,'per']), 97.3640857, tolerance=1e-7)
    expect_equal(as.numeric(model_e3_bw.DS$study1$centiles[9,'per']), 99.3410214, tolerance=1e-7)
})

context("ds.gamlss::smk::mod2")
test_that("gamlss, Z BMI", {
  model_zbmi_who.DS <- ds.gamlss(formula = 'hs_zbmi_who ~ hs_child_age_None+h_mbmi_None',
                                 sigma.formula = 'hs_zbmi_who ~ hs_child_age_None',
                                 data = 'D', family = 'NO()', centiles = TRUE,
                                 xvar = 'D$hs_child_age_None', newobj = 'z_scores_hs_zbmi_who') 
  
  expect_length(model_zbmi_who.DS, 3)
  expect_length(model_zbmi_who.DS$study1, 2)
  expect_length(model_zbmi_who.DS$study2, 2)
  expect_length(model_zbmi_who.DS$study3, 2)
  expect_equal(model_zbmi_who.DS$study1$results$family, c("NO","Normal"))
  expect_true("character" %in% class(model_zbmi_who.DS$study1$results$family))
  expect_equal(model_zbmi_who.DS$study1$results$parameters, c("mu","sigma"))
  expect_true("character" %in% class(model_zbmi_who.DS$study1$results$parameters))
  expect_true("call" %in% class(model_zbmi_who.DS$study1$results$call))
  expect_equal(model_zbmi_who.DS$study1$results$y, "The response variable is not disclosed!")
  expect_true("character" %in% class(model_zbmi_who.DS$study1$results$y))
  expect_true("list" %in% class(model_zbmi_who.DS$study1$results$control))
  expect_equal(model_zbmi_who.DS$study1$results$control$c.crit, 0.001)
  expect_equal(model_zbmi_who.DS$study1$results$control$n.cyc, 20)
  expect_equal(model_zbmi_who.DS$study1$results$control$mu.step, 1)
  expect_equal(model_zbmi_who.DS$study1$results$control$sigma.step, 1)
  expect_equal(model_zbmi_who.DS$study1$results$control$nu.step, 1)
  expect_equal(model_zbmi_who.DS$study1$results$control$tau.step, 1)
  expect_equal(model_zbmi_who.DS$study1$results$control$gd.tol, Inf)
  expect_equal(model_zbmi_who.DS$study1$results$control$iter, 0)
  expect_equal(model_zbmi_who.DS$study1$results$control$trace, TRUE)
  expect_equal(model_zbmi_who.DS$study1$results$control$autostep, TRUE)
  expect_equal(model_zbmi_who.DS$study1$results$control$save, TRUE)
  expect_true("numeric" %in% class(model_zbmi_who.DS$study1$results$weights))
  expect_equal(model_zbmi_who.DS$study1$results$G.deviance, 1979.397, tolerance=1e-05)
  expect_equal(model_zbmi_who.DS$study1$results$N, 607)
  expect_true("expression" %in% class(model_zbmi_who.DS$study1$results$rqres))
  expect_equal(model_zbmi_who.DS$study1$results$iter, 3)
  expect_equal(model_zbmi_who.DS$study1$results$type, "Continuous")
  expect_true("character" %in% class(model_zbmi_who.DS$study1$results$type))
  expect_true("call" %in% class(model_zbmi_who.DS$study1$results$method))
  expect_equal(model_zbmi_who.DS$study1$results$contrasts, NULL)
  expect_true("NULL" %in% class(model_zbmi_who.DS$study1$results$contrasts))
  expect_equal(model_zbmi_who.DS$study1$results$converged, TRUE)
  expect_true("logical" %in% class(model_zbmi_who.DS$study1$results$converged))
  expect_equal(model_zbmi_who.DS$study1$results$residuals, "The residuals of the model are not disclosed!")
  expect_true("character" %in% class(model_zbmi_who.DS$study1$results$residuals))
  expect_equal(model_zbmi_who.DS$study1$results$noObs, 607)
  expect_equal(model_zbmi_who.DS$study1$results$mu.fv, "The fitted values of the mu model are not disclosed!")
  expect_equal(model_zbmi_who.DS$study1$results$mu.link, "identity")
  expect_true("character" %in% class(model_zbmi_who.DS$study1$results$mu.link))
  expect_true("terms" %in% class(model_zbmi_who.DS$study1$results$mu.terms))
  expect_true("formula" %in% class(model_zbmi_who.DS$study1$results$mu.terms))
  expect_equal(as.numeric(model_zbmi_who.DS$study1$results$mu.coefficients["h_mbmi_None"]), 0.04635936, tolerance=1e-5)
  expect_equal(as.numeric(model_zbmi_who.DS$study1$results$mu.coefficients["hs_child_age_None"]), -0.10898234, tolerance=1e-5)
  expect_equal(as.numeric(model_zbmi_who.DS$study1$results$mu.coefficients["(Intercept)"]), -0.06404270, tolerance=1e-4)
  expect_true("numeric" %in% class(model_zbmi_who.DS$study1$results$mu.offset))
  expect_true("formula" %in% class(model_zbmi_who.DS$study1$results$mu.formula))
  expect_true("list" %in% class(model_zbmi_who.DS$study1$results$mu.xlevels))
  expect_equal(model_zbmi_who.DS$study1$results$mu.df, 3)
  expect_equal(model_zbmi_who.DS$study1$results$mu.nl.df, 0)
  expect_equal(model_zbmi_who.DS$study1$results$mu.pen, 0)
  expect_equal(model_zbmi_who.DS$study1$results$P.deviance, 1979.397, tolerance=1e-5)
  expect_equal(model_zbmi_who.DS$study1$results$df.fit, 5)
  expect_equal(model_zbmi_who.DS$study1$results$pen, 0)
  expect_equal(model_zbmi_who.DS$study1$results$df.residual, 602)
  expect_equal(model_zbmi_who.DS$study1$results$aic, 1989.397, tolerance=1e-5)
  expect_equal(model_zbmi_who.DS$study1$results$sbc, 2011.439, tolerance=1e-5)
  expect_true("matrix" %in% class(model_zbmi_who.DS$study1$centiles))
  expect_true("array" %in% class(model_zbmi_who.DS$study1$centiles))
  expect_equal(as.numeric(model_zbmi_who.DS$study1$centiles[1,'cent']), 0.4)
  expect_equal(as.numeric(model_zbmi_who.DS$study1$centiles[2,'cent']), 2.0)
  expect_equal(as.numeric(model_zbmi_who.DS$study1$centiles[3,'cent']), 10.0)
  expect_equal(as.numeric(model_zbmi_who.DS$study1$centiles[4,'cent']), 25.0)
  expect_equal(as.numeric(model_zbmi_who.DS$study1$centiles[5,'cent']), 50.0)
  expect_equal(as.numeric(model_zbmi_who.DS$study1$centiles[6,'cent']), 75.0)
  expect_equal(as.numeric(model_zbmi_who.DS$study1$centiles[7,'cent']), 90.0)
  expect_equal(as.numeric(model_zbmi_who.DS$study1$centiles[8,'cent']), 98.0)
  expect_equal(as.numeric(model_zbmi_who.DS$study1$centiles[9,'cent']), 99.6)
  expect_equal(as.numeric(model_zbmi_who.DS$study1$centiles[1,'per']), 0.1647446, tolerance=1e-7)
  expect_equal(as.numeric(model_zbmi_who.DS$study1$centiles[2,'per']), 0.9884679, tolerance=1e-7)
  expect_equal(as.numeric(model_zbmi_who.DS$study1$centiles[3,'per']), 8.2372323, tolerance=1e-7)
  expect_equal(as.numeric(model_zbmi_who.DS$study1$centiles[4,'per']), 24.5469522, tolerance=1e-7)
  expect_equal(as.numeric(model_zbmi_who.DS$study1$centiles[5,'per']), 54.5304778, tolerance=1e-7)
  expect_equal(as.numeric(model_zbmi_who.DS$study1$centiles[6,'per']), 78.7479407, tolerance=1e-7)
  expect_equal(as.numeric(model_zbmi_who.DS$study1$centiles[7,'per']), 89.6210873, tolerance=1e-7)
  expect_equal(as.numeric(model_zbmi_who.DS$study1$centiles[8,'per']), 96.0461285, tolerance=1e-7)
  expect_equal(as.numeric(model_zbmi_who.DS$study1$centiles[9,'per']), 99.1762768, tolerance=1e-7)
})


context("ds.gamlss::smk::mod3")
test_that("gamlss, RAVEN TEST", {

  ds.dataFrameSort(df.name = 'D', sort.key.name = 'D$hs_correct_raven', newobj = 'D2')

  model_correct_raven.DS <- ds.gamlss(formula = 'hs_correct_raven ~ cs(hs_child_age_None)',
                                      sigma.formula = 'hs_correct_raven ~ (hs_child_age_None)',
                                      nu.formula = 'hs_correct_raven ~ hs_child_age_None',
                                      tau.formula = 'hs_correct_raven ~ hs_child_age_None',
                                      data = 'D2', family = 'BCT()', centiles = TRUE,
                                      xvar = 'D$hs_child_age_None', newobj = 'z_scores_hs_correct_raven')
  
  expect_length(model_correct_raven.DS, 3)
  expect_length(model_correct_raven.DS$study1, 2)
  expect_length(model_correct_raven.DS$study2, 2)
  expect_length(model_correct_raven.DS$study3, 2)
  expect_equal(model_correct_raven.DS$study1$results$family, c("BCT","Box-Cox t"))
  expect_true("character" %in% class(model_correct_raven.DS$study1$results$family))
  expect_equal(model_correct_raven.DS$study1$results$parameters, c("mu","sigma","nu","tau"))
  expect_true("character" %in% class(model_correct_raven.DS$study1$results$parameters))
  expect_true("call" %in% class(model_correct_raven.DS$study1$results$call))
  expect_equal(model_correct_raven.DS$study1$results$y, "The response variable is not disclosed!")
  expect_true("character" %in% class(model_correct_raven.DS$study1$results$y))
  expect_true("list" %in% class(model_correct_raven.DS$study1$results$control))
  expect_equal(model_correct_raven.DS$study1$results$control$c.crit, 0.001)
  expect_equal(model_correct_raven.DS$study1$results$control$n.cyc, 20)
  expect_equal(model_correct_raven.DS$study1$results$control$mu.step, 1)
  expect_equal(model_correct_raven.DS$study1$results$control$sigma.step, 1)
  expect_equal(model_correct_raven.DS$study1$results$control$nu.step, 1)
  expect_equal(model_correct_raven.DS$study1$results$control$tau.step, 1)
  expect_equal(model_correct_raven.DS$study1$results$control$gd.tol, Inf)
  expect_equal(model_correct_raven.DS$study1$results$control$iter, 0)
  expect_equal(model_correct_raven.DS$study1$results$control$trace, TRUE)
  expect_equal(model_correct_raven.DS$study1$results$control$autostep, TRUE)
  expect_equal(model_correct_raven.DS$study1$results$control$save, TRUE)
  expect_true("numeric" %in% class(model_correct_raven.DS$study1$results$weights))
  expect_equal(model_correct_raven.DS$study1$results$G.deviance, 3723.97, tolerance=1e-01)
  expect_equal(model_correct_raven.DS$study1$results$N, 607)
  expect_true("expression" %in% class(model_correct_raven.DS$study1$results$rqres))
  expect_equal(model_correct_raven.DS$study1$results$iter, 7)
  expect_equal(model_correct_raven.DS$study1$results$type, "Continuous")
  expect_true("character" %in% class(model_correct_raven.DS$study1$results$type))
  expect_true("call" %in% class(model_correct_raven.DS$study1$results$method))
  expect_equal(model_correct_raven.DS$study1$results$contrasts, NULL)
  expect_true("NULL" %in% class(model_correct_raven.DS$study1$results$contrasts))
  expect_equal(model_correct_raven.DS$study1$results$converged, TRUE)
  expect_true("logical" %in% class(model_correct_raven.DS$study1$results$converged))
  expect_equal(model_correct_raven.DS$study1$results$residuals, "The residuals of the model are not disclosed!")
  expect_true("character" %in% class(model_correct_raven.DS$study1$results$residuals))
  expect_equal(model_correct_raven.DS$study1$results$noObs, 607)
  expect_equal(model_correct_raven.DS$study1$results$mu.fv, "The fitted values of the mu model are not disclosed!")
  expect_equal(model_correct_raven.DS$study1$results$mu.link, "identity")
  expect_true("character" %in% class(model_correct_raven.DS$study1$results$mu.link))
  expect_true("terms" %in% class(model_correct_raven.DS$study1$results$mu.terms))
  expect_true("formula" %in% class(model_correct_raven.DS$study1$results$mu.terms))
  expect_equal(as.numeric(model_correct_raven.DS$study1$results$mu.coefficients["cs(hs_child_age_None)"]), 2.631276, tolerance=1e-1)
  expect_equal(as.numeric(model_correct_raven.DS$study1$results$mu.coefficients["(Intercept)"]), 4.180341, tolerance=1e-1)
  expect_true("numeric" %in% class(model_correct_raven.DS$study1$results$mu.offset))
  expect_true("formula" %in% class(model_correct_raven.DS$study1$results$mu.formula))
  expect_true("list" %in% class(model_correct_raven.DS$study1$results$mu.xlevels))
  expect_equal(model_correct_raven.DS$study1$results$mu.df, 5.000742, tolerance=1e-7)
  expect_equal(model_correct_raven.DS$study1$results$mu.nl.df, 3.000742, tolerance=1e-6)
  expect_equal(model_correct_raven.DS$study1$results$mu.pen, 0.8507004, tolerance=1e-2)
  expect_equal(model_correct_raven.DS$study1$results$P.deviance, 3724.82, tolerance=1e-5)
  expect_equal(model_correct_raven.DS$study1$results$df.fit, 11.00074, tolerance=1e-6)
  expect_equal(model_correct_raven.DS$study1$results$pen, 0.8506994, tolerance=1e-2
)
  expect_equal(model_correct_raven.DS$study1$results$df.residual, 595.9993, tolerance=1e-7)
  expect_equal(model_correct_raven.DS$study1$results$aic, 3745.971, tolerance=1e-5)
  expect_equal(model_correct_raven.DS$study1$results$sbc, 3794.468, tolerance=1e-5)
  expect_true("matrix" %in% class(model_correct_raven.DS$study1$centiles))
  expect_true("array" %in% class(model_correct_raven.DS$study1$centiles))
  expect_equal(as.numeric(model_correct_raven.DS$study1$centiles[1,'cent']), 0.4)
  expect_equal(as.numeric(model_correct_raven.DS$study1$centiles[2,'cent']), 2.0)
  expect_equal(as.numeric(model_correct_raven.DS$study1$centiles[3,'cent']), 10.0)
  expect_equal(as.numeric(model_correct_raven.DS$study1$centiles[4,'cent']), 25.0)
  expect_equal(as.numeric(model_correct_raven.DS$study1$centiles[5,'cent']), 50.0)
  expect_equal(as.numeric(model_correct_raven.DS$study1$centiles[6,'cent']), 75.0)
  expect_equal(as.numeric(model_correct_raven.DS$study1$centiles[7,'cent']), 90.0)
  expect_equal(as.numeric(model_correct_raven.DS$study1$centiles[8,'cent']), 98.0)
  expect_equal(as.numeric(model_correct_raven.DS$study1$centiles[9,'cent']), 99.6)
  expect_equal(as.numeric(model_correct_raven.DS$study1$centiles[1,'per']), 0.000000, tolerance=1e-7)
  expect_equal(as.numeric(model_correct_raven.DS$study1$centiles[2,'per']), 1.153213, tolerance=1e-6)
  expect_equal(as.numeric(model_correct_raven.DS$study1$centiles[3,'per']), 13.014827, tolerance=1e-7)
  expect_equal(as.numeric(model_correct_raven.DS$study1$centiles[4,'per']), 26.359143, tolerance=1e-7)
  expect_equal(as.numeric(model_correct_raven.DS$study1$centiles[5,'per']), 47.775947, tolerance=1e-7)
  expect_equal(as.numeric(model_correct_raven.DS$study1$centiles[6,'per']), 74.135091, tolerance=1e-7)
  expect_equal(as.numeric(model_correct_raven.DS$study1$centiles[7,'per']), 90.115321, tolerance=1e-7)
  expect_equal(as.numeric(model_correct_raven.DS$study1$centiles[8,'per']), 98.352554, tolerance=1e-7)
  expect_equal(as.numeric(model_correct_raven.DS$study1$centiles[9,'per']), 100.000000, tolerance=1e-7)
})

context("ds.gamlss::smk::mod4")
test_that("gamlss, WEIGHT GAIN", {
  model_wgtgain.DS <- ds.gamlss(formula = 'hs_wgtgain_None ~ (e3_bw)',
                                sigma.formula = 'hs_wgtgain_None ~ (e3_bw)',
                                data = 'D', family = 'NO()', centiles = TRUE,
                                xvar = 'D$hs_wgtgain_None', newobj = 'z_scores_hs_wgtgain_None') 
  
  expect_length(model_wgtgain.DS, 3)
  expect_length(model_wgtgain.DS$study1, 2)
  expect_length(model_wgtgain.DS$study2, 2)
  expect_length(model_wgtgain.DS$study3, 2)
  expect_equal(model_wgtgain.DS$study1$results$family, c("NO","Normal"))
  expect_true("character" %in% class(model_wgtgain.DS$study1$results$family))
  expect_equal(model_wgtgain.DS$study1$results$parameters, c("mu","sigma"))
  expect_true("character" %in% class(model_wgtgain.DS$study1$results$parameters))
  expect_true("call" %in% class(model_wgtgain.DS$study1$results$call))
  expect_equal(model_wgtgain.DS$study1$results$y, "The response variable is not disclosed!")
  expect_true("character" %in% class(model_wgtgain.DS$study1$results$y))
  expect_true("list" %in% class(model_wgtgain.DS$study1$results$control))
  expect_equal(model_wgtgain.DS$study1$results$control$c.crit, 0.001)
  expect_equal(model_wgtgain.DS$study1$results$control$n.cyc, 20)
  expect_equal(model_wgtgain.DS$study1$results$control$mu.step, 1)
  expect_equal(model_wgtgain.DS$study1$results$control$sigma.step, 1)
  expect_equal(model_wgtgain.DS$study1$results$control$nu.step, 1)
  expect_equal(model_wgtgain.DS$study1$results$control$tau.step, 1)
  expect_equal(model_wgtgain.DS$study1$results$control$gd.tol, Inf)
  expect_equal(model_wgtgain.DS$study1$results$control$iter, 0)
  expect_equal(model_wgtgain.DS$study1$results$control$trace, TRUE)
  expect_equal(model_wgtgain.DS$study1$results$control$autostep, TRUE)
  expect_equal(model_wgtgain.DS$study1$results$control$save, TRUE)
  expect_true("numeric" %in% class(model_wgtgain.DS$study1$results$weights))
  expect_equal(model_wgtgain.DS$study1$results$G.deviance, 3973.189, tolerance=1e-07)
  expect_equal(model_wgtgain.DS$study1$results$N, 607)
  expect_true("expression" %in% class(model_wgtgain.DS$study1$results$rqres))
  expect_equal(model_wgtgain.DS$study1$results$iter, 3)
  expect_equal(model_wgtgain.DS$study1$results$type, "Continuous")
  expect_true("character" %in% class(model_wgtgain.DS$study1$results$type))
  expect_true("call" %in% class(model_wgtgain.DS$study1$results$method))
  expect_equal(model_wgtgain.DS$study1$results$contrasts, NULL)
  expect_true("NULL" %in% class(model_wgtgain.DS$study1$results$contrasts))
  expect_equal(model_wgtgain.DS$study1$results$converged, TRUE)
  expect_true("logical" %in% class(model_wgtgain.DS$study1$results$converged))
  expect_equal(model_wgtgain.DS$study1$results$residuals, "The residuals of the model are not disclosed!")
  expect_true("character" %in% class(model_wgtgain.DS$study1$results$residuals))
  expect_equal(model_wgtgain.DS$study1$results$noObs, 607)
  expect_equal(model_wgtgain.DS$study1$results$mu.fv, "The fitted values of the mu model are not disclosed!")
  expect_equal(model_wgtgain.DS$study1$results$mu.link, "identity")
  expect_true("character" %in% class(model_wgtgain.DS$study1$results$mu.link))
  expect_true("terms" %in% class(model_wgtgain.DS$study1$results$mu.terms))
  expect_true("formula" %in% class(model_wgtgain.DS$study1$results$mu.terms))
  expect_equal(as.numeric(model_wgtgain.DS$study1$results$mu.coefficients["e3_bw"]), 0.002408255, tolerance=1e-5)
  expect_equal(as.numeric(model_wgtgain.DS$study1$results$mu.coefficients["(Intercept)"]), 5.230886349, tolerance=1e-5)
  expect_true("numeric" %in% class(model_wgtgain.DS$study1$results$mu.offset))
  expect_true("formula" %in% class(model_wgtgain.DS$study1$results$mu.formula))
  expect_true("list" %in% class(model_wgtgain.DS$study1$results$mu.xlevels))
  expect_equal(model_wgtgain.DS$study1$results$mu.df, 2)
  expect_equal(model_wgtgain.DS$study1$results$mu.nl.df, 0)
  expect_equal(model_wgtgain.DS$study1$results$mu.pen, 0)
  expect_equal(model_wgtgain.DS$study1$results$P.deviance, 3973.189, tolerance=1e-5)
  expect_equal(model_wgtgain.DS$study1$results$df.fit, 4)
  expect_equal(model_wgtgain.DS$study1$results$pen, 0)
  expect_equal(model_wgtgain.DS$study1$results$df.residual, 603)
  expect_equal(model_wgtgain.DS$study1$results$aic, 3981.189, tolerance=1e-5)
  expect_equal(model_wgtgain.DS$study1$results$sbc, 3998.823, tolerance=1e-5)
  expect_true("matrix" %in% class(model_wgtgain.DS$study1$centiles))
  expect_true("array" %in% class(model_wgtgain.DS$study1$centiles))
  expect_equal(as.numeric(model_wgtgain.DS$study1$centiles[1,'cent']), 0.4)
  expect_equal(as.numeric(model_wgtgain.DS$study1$centiles[2,'cent']), 2.0)
  expect_equal(as.numeric(model_wgtgain.DS$study1$centiles[3,'cent']), 10.0)
  expect_equal(as.numeric(model_wgtgain.DS$study1$centiles[4,'cent']), 25.0)
  expect_equal(as.numeric(model_wgtgain.DS$study1$centiles[5,'cent']), 50.0)
  expect_equal(as.numeric(model_wgtgain.DS$study1$centiles[6,'cent']), 75.0)
  expect_equal(as.numeric(model_wgtgain.DS$study1$centiles[7,'cent']), 90.0)
  expect_equal(as.numeric(model_wgtgain.DS$study1$centiles[8,'cent']), 98.0)
  expect_equal(as.numeric(model_wgtgain.DS$study1$centiles[9,'cent']), 99.6)
  expect_equal(as.numeric(model_wgtgain.DS$study1$centiles[1,'per']), 0.000000, tolerance=1e-7)
  expect_equal(as.numeric(model_wgtgain.DS$study1$centiles[2,'per']), 1.317957, tolerance=1e-6)
  expect_equal(as.numeric(model_wgtgain.DS$study1$centiles[3,'per']), 6.919275, tolerance=1e-7)
  expect_equal(as.numeric(model_wgtgain.DS$study1$centiles[4,'per']), 25.864909, tolerance=1e-7)
  expect_equal(as.numeric(model_wgtgain.DS$study1$centiles[5,'per']), 56.342669, tolerance=1e-7)
  expect_equal(as.numeric(model_wgtgain.DS$study1$centiles[6,'per']), 74.958814, tolerance=1e-7)
  expect_equal(as.numeric(model_wgtgain.DS$study1$centiles[7,'per']), 90.280066, tolerance=1e-7)
  expect_equal(as.numeric(model_wgtgain.DS$study1$centiles[8,'per']), 94.892916, tolerance=1e-7)
  expect_equal(as.numeric(model_wgtgain.DS$study1$centiles[9,'per']), 98.682043, tolerance=1e-7)
})

#
# Done
#

context("ds.gamlss::smk::shutdown")

test_that("shutdown", {
    ds_expect_variables(c("D", "D2", "z_scores_e3_bw", "z_scores_hs_correct_raven", "z_scores_hs_wgtgain_None",
                          "z_scores_hs_zbmi_who"))
})

disconnect.studies.dataset.gamlss()

context("ds.gamlss::smk::done")