test_that("basic model fits", {
  set.seed(1); d <- data.frame(y = rnorm(20), x = rnorm(20))
  m <- fit_biomix(y ~ x, d)
  expect_s3_class(m, "biomix")
  expect_true(is.numeric(haic(m)))
  expect_true(is.data.frame(variance_partition(m)))
  expect_true(is.numeric(icc_biomix(m)))
  expect_true(is.data.frame(model_summary(m)$variance))
  expect_true(is.data.frame(diagnose_biomix(m)))
})

test_that("random model fits when nlme is available", {
  
  skip_if_not_installed("nlme")
  
  set.seed(123)
  
  dat <- data.frame(
    id = factor(rep(1:10, each = 5)),
    x = rep(1:5, times = 10)
  )
  
  dat$y <- 10 +
    2 * dat$x +
    rnorm(10)[dat$id] +
    rnorm(50)
  
  m <- fit_biomix(
    y ~ x,
    data = dat,
    random = ~ 1 | id
  )
  
  expect_s3_class(m, "biomix")
  expect_true(is.finite(icc_biomix(m)))
  expect_true(is.finite(haic(m)))
})

library(testthat)
library(BioMixModel)


# ---------------------------------------------------------
# Test data
# ---------------------------------------------------------

set.seed(123)

dat <- data.frame(
  id = factor(rep(1:10, each = 5)),
  x = rep(1:5, times = 10)
)

dat$y <- 10 +
  2 * dat$x +
  rnorm(10)[dat$id] +
  rnorm(50, 0, 1)


# ---------------------------------------------------------
# fit_biomix()
# ---------------------------------------------------------

test_that("fit_biomix fits ordinary linear model", {
  
  m <- fit_biomix(
    y ~ x,
    data = dat
  )
  
  expect_s3_class(m, "biomix")
  expect_s3_class(m$fixed, "lm")
  expect_null(m$random)
  expect_equal(m$method, "REML")
  expect_equal(m$formula, y ~ x)
})


test_that("fit_biomix fits mixed model with REML", {
  
  skip_if_not_installed("nlme")
  
  m <- fit_biomix(
    y ~ x,
    data = dat,
    random = ~ 1 | id
  )
  
  expect_s3_class(m, "biomix")
  expect_s3_class(m$fixed, "lme")
  expect_equal(m$method, "REML")
  expect_equal(m$random_formula, ~1 | id)
})


test_that("fit_biomix fits mixed model with ML", {
  
  skip_if_not_installed("nlme")
  
  m <- fit_biomix(
    y ~ x,
    data = dat,
    random = ~ 1 | id,
    method = "ML"
  )
  
  expect_s3_class(m, "biomix")
  expect_s3_class(m$fixed, "lme")
  expect_equal(m$method, "ML")
})


test_that("fit_biomix rejects invalid method", {
  
  expect_error(
    fit_biomix(
      y ~ x,
      data = dat,
      method = "INVALID"
    )
  )
})


# ---------------------------------------------------------
# variance_partition()
# ---------------------------------------------------------

test_that("variance_partition works for ordinary model", {
  
  m <- fit_biomix(
    y ~ x,
    data = dat
  )
  
  vp <- variance_partition(m)
  
  expect_true(is.data.frame(vp))
  expect_true(all(c(
    "component",
    "variance",
    "proportion"
  ) %in% names(vp)))
  
  expect_true(all(is.finite(vp$variance)))
  expect_true(all(is.finite(vp$proportion)))
})


test_that("variance_partition works for mixed model", {
  
  skip_if_not_installed("nlme")
  
  m <- fit_biomix(
    y ~ x,
    data = dat,
    random = ~ 1 | id
  )
  
  vp <- variance_partition(m)
  
  expect_true(is.data.frame(vp))
  expect_true("Residual" %in% vp$component)
  expect_true(all(vp$variance >= 0))
  expect_equal(
    sum(vp$proportion),
    1,
    tolerance = 1e-8
  )
})


# ---------------------------------------------------------
# icc_biomix()
# ---------------------------------------------------------

test_that("icc_biomix returns finite ICC for mixed model", {
  
  skip_if_not_installed("nlme")
  
  m <- fit_biomix(
    y ~ x,
    data = dat,
    random = ~ 1 | id
  )
  
  icc <- icc_biomix(m)
  
  expect_length(icc, 1)
  expect_true(is.finite(icc))
  expect_gte(icc, 0)
  expect_lte(icc, 1)
})


test_that("icc_biomix returns zero for ordinary model", {
  
  m <- fit_biomix(
    y ~ x,
    data = dat
  )
  
  icc <- icc_biomix(m)
  
  expect_length(icc, 1)
  expect_equal(icc, 0)
})


# ---------------------------------------------------------
# haic()
# ---------------------------------------------------------

test_that("haic returns finite value", {
  
  m <- fit_biomix(
    y ~ x,
    data = dat
  )
  
  h <- haic(m)
  
  expect_length(h, 1)
  expect_true(is.finite(h))
})


test_that("haic works for mixed model", {
  
  skip_if_not_installed("nlme")
  
  m <- fit_biomix(
    y ~ x,
    data = dat,
    random = ~ 1 | id
  )
  
  h <- haic(m)
  
  expect_length(h, 1)
  expect_true(is.finite(h))
})


# ---------------------------------------------------------
# compare_biomix()
# ---------------------------------------------------------

test_that("compare_biomix compares multiple models", {
  
  m1 <- fit_biomix(
    y ~ x,
    data = dat
  )
  
  m2 <- fit_biomix(
    y ~ 1,
    data = dat
  )
  
  result <- compare_biomix(m1, m2)
  
  expect_true(is.data.frame(result))
  expect_equal(nrow(result), 2)
  
  expect_true(all(c(
    "model",
    "AIC",
    "HAIC",
    "logLik"
  ) %in% names(result)))
  
  expect_true(all(is.finite(result$AIC)))
  expect_true(all(is.finite(result$HAIC)))
  expect_true(all(is.finite(result$logLik)))
})


test_that("compare_biomix requires a model", {
  
  expect_error(
    compare_biomix(),
    "Supply at least one model"
  )
})


# ---------------------------------------------------------
# model_summary()
# ---------------------------------------------------------

test_that("model_summary returns expected components", {
  
  m <- fit_biomix(
    y ~ x,
    data = dat
  )
  
  result <- model_summary(m)
  
  expect_type(result, "list")
  
  expect_true(all(c(
    "n",
    "AIC",
    "HAIC",
    "ICC",
    "variance"
  ) %in% names(result)))
  
  expect_equal(result$n, nrow(dat))
  expect_true(is.finite(result$AIC))
  expect_true(is.finite(result$HAIC))
  expect_true(is.finite(result$ICC))
  expect_true(is.data.frame(result$variance))
})


# ---------------------------------------------------------
# diagnose_biomix()
# ---------------------------------------------------------

test_that("diagnose_biomix returns diagnostics", {
  
  m <- fit_biomix(
    y ~ x,
    data = dat
  )
  
  result <- diagnose_biomix(m)
  
  expect_true(is.data.frame(result))
  
  expect_true(all(c(
    "mean",
    "sd",
    "min",
    "max",
    "n"
  ) %in% names(result)))
  
  expect_true(is.finite(result$mean))
  expect_true(is.finite(result$sd))
  expect_true(is.finite(result$min))
  expect_true(is.finite(result$max))
  expect_equal(result$n, nrow(dat))
})


# ---------------------------------------------------------
# plot_biomix()
# ---------------------------------------------------------

test_that("plot_biomix produces residual plot", {
  
  m <- fit_biomix(
    y ~ x,
    data = dat
  )
  
  expect_invisible(
    plot_biomix(m, type = "residuals")
  )
})


test_that("plot_biomix produces QQ plot", {
  
  m <- fit_biomix(
    y ~ x,
    data = dat
  )
  
  expect_invisible(
    plot_biomix(m, type = "qq")
  )
})


test_that("plot_biomix rejects invalid plot type", {
  
  m <- fit_biomix(
    y ~ x,
    data = dat
  )
  
  expect_error(
    plot_biomix(m, type = "invalid")
  )
})

test_that("fit_biomix rejects invalid formula", {
  
  expect_error(
    fit_biomix(
      "y ~ x",
      data.frame(
        y = rnorm(10),
        x = rnorm(10)
      )
    )
  )
})

test_that("fit_biomix rejects invalid data", {
  
  expect_error(
    fit_biomix(
      y ~ x,
      list(
        y = rnorm(10),
        x = rnorm(10)
      )
    )
  )
})

test_that("fit_biomix rejects invalid random formula", {
  
  skip_if_not_installed("nlme")
  
  dat <- data.frame(
    y = rnorm(20),
    x = rnorm(20),
    id = factor(rep(1:10, each = 2))
  )
  
  expect_error(
    fit_biomix(
      y ~ x,
      data = dat,
      random = "1 | id"
    )
  )
})