Fixed tolerances in the tests, added tests on TDcond, mixCond and BUI…

…S to cover cases not covered.

tests/testthat/dataForTests/Monthly-Count_ts.csv

@@ -0,0 +1,120 @@
+4
+4
+1
+3
+2
+2
+3
+0
+2
+3
+2
+3
+4
+1
+2
+4
+5
+0
+2
+2
+4
+1
+6
+4
+0
+2
+1
+4
+2
+3
+3
+3
+1
+3
+0
+3
+0
+1
+4
+2
+1
+2
+0
+5
+2
+5
+4
+2
+5
+2
+1
+1
+1
+3
+0
+3
+3
+1
+1
+2
+2
+5
+3
+2
+3
+1
+1
+3
+3
+1
+0
+1
+1
+2
+1
+3
+0
+1
+2
+0
+2
+1
+1
+2
+3
+2
+1
+0
+0
+1
+2
+0
+1
+4
+4
+3
+1
+2
+3
+2
+2
+1
+1
+1
+4
+5
+3
+3
+2
+0
+2
+3
+3
+2
+2
+2
+0
+1
+2
+3

tests/testthat/dataForTests/generate_monthlyCountData.R

@@ -0,0 +1,12 @@
+# Generate the monthly count time series for the mixCond and TDcond tests
+# CHANGE THE WORKING DIRECTORY BEFORE RUNNING
+rm(list=ls())
+library(bayesRecon)
+
+set.seed(42)
+vals <- stats::rpois(12*10,lambda = 2)
+
+
+write.table(vals,file="./Monthly-Count_ts.csv",row.names = FALSE,sep=',',
+            col.names = FALSE,quote = FALSE)
+

tests/testthat/test-reconc_BUIS_gaussian.R

@@ -106,4 +106,48 @@ test_that("Monthly, in_type=='samples', distr='discrete'",{
   expect_equal(abs(m) < 1.5e-2, TRUE)
 })
 
+test_that("Monthly simple, in_type=='params', distr='nbinom'",{
+
+  # Read samples from dataForTests (reproducibility)
+  vals <- read.csv(file = "dataForTests/Monthly-Count_ts.csv", header = FALSE)
+
+  # Create a count time series with monthly observations for 10 years
+  y <- ts(data=vals,frequency = 12)
+
+  # Create the aggregated yearly time series
+  y_agg <- temporal_aggregation(y,agg_levels = c(1,12))
+
+  # We use a marginal forecast that computes for each month 
+  # the empirical mean and variance
+  # the forecast is a negative binomial with those params
+  fc_bottom <- list()
+  for(i in seq(12)){
+    mm <- mean(y_agg$`f=12`[seq(i,120,12)])
+    vv <- max(var(y_agg$`f=12`[seq(i,120,12)]), mm+0.5)
+    #cat("i: ",i, "mean: ",mm, "var: ",vv, "size: ",mm^2/(vv-mm), "prob: ",mm/vv, "\n")
+
+    fc_bottom[[i]] <- list(size=mm^2/(vv-mm),mu=mm)
+  }
+
+  # We compute the empirical mean and variance of the yearly ts
+  # we forecast with a negative binomial with those parameters
+  mm <- mean(y_agg$`f=1`)
+  vv <- var(y_agg$`f=1`)
+  fc_upper <- list(size=mm^2/(vv-mm), prob= mm/vv)
+
+  # Obtain the aggregation matrix for this hierarchy
+  rec_mat <- get_reconc_matrices(c(1,12),12)
+
+  base_forecasts = append(list(fc_upper),fc_bottom)
+  res.buis_params = reconc_BUIS(rec_mat$A, base_forecasts, in_type = "params", distr = "nbinom", seed=42)
+
+
+  fc_upper_gauss <- list(mu=mm, Sigma = matrix(vv))
+  res.mixCond <- reconc_MixCond(rec_mat$A, fc_bottom, fc_upper_gauss, bottom_in_type = "params", distr = 'nbinom')
+  upp_pmf <- PMF.from_samples(as.integer(res.buis_params$upper_reconciled_samples))
+
+  expect_equal(res.mixCond$upper_reconciled$pmf[[1]],upp_pmf,tolerance = 0.1)
+
+})
+
 ##############################################################################

tests/testthat/test-reconc_MixCond.R

@@ -23,7 +23,7 @@ test_that("reconc_MixCond simple example", {
   ## upper
   fc_upper <- list(mu = means[1:10],
                    Sigma = diag(vars[1:10]))
-
+  
   ## bottom
   fc_bottom <- list()
   for(i in seq(ncol(A))){
@@ -50,7 +50,65 @@ test_that("reconc_MixCond simple example", {
   bott_rec_means_pmf <- unlist(lapply(res.MixCond_pmf$bottom_reconciled$pmf,PMF.get_mean))
   bott_rec_vars_pmf <- unlist(lapply(res.MixCond_pmf$bottom_reconciled$pmf,PMF.get_var))
 
-  expect_equal(bott_rec_means,bott_rec_means_pmf, tolerance = 0.1)
+  expect_equal(bott_rec_means,bott_rec_means_pmf, tolerance = 0.01)
   expect_equal(bott_rec_vars,bott_rec_vars_pmf, tolerance = 0.1)
+
+})
+
+test_that("reconc_MixCond and reconc_TDcond with temporal hier and params", {
+
+
+  # Read samples from dataForTests (reproducibility)
+  vals <- read.csv(file = "dataForTests/Monthly-Count_ts.csv", header = FALSE)
+
+  # Create a count time series with monthly observations for 10 years
+  y <- ts(data=vals,frequency = 12)
+
+  # Create the aggregated yearly time series
+  y_agg <- temporal_aggregation(y,agg_levels = c(1,12))
+
+  # We use a marginal forecast that computes for each month 
+  # the empirical mean and forecasts a Poisson with that value
+  fc_bottom <- list()
+  for(i in seq(12)){
+    fc_bottom[[i]] <- list(lambda=mean(y_agg$`f=12`[seq(i,120,12)]))
+  }
+
+  # We compute the empirical mean and variance of the yearly ts
+  # we forecast with a Gaussian with those parameters
+  fc_upper <- list(mu=mean(y_agg$`f=1`), Sigma=matrix(var(y_agg$`f=1`)))
+
+  # Obtain the aggregation matrix for this hierarchy
+  rec_mat <- get_reconc_matrices(c(1,12),12)
+
+  # Do a couple of checks on S and A
+  expect_no_error(.check_S(rec_mat$S))
+  expect_error(.check_S(rec_mat$A))
+  expect_true(.check_BU_matr(rec_mat$A))
+  expect_false(.check_BU_matr(rec_mat$S))
+
+  # We can reconcile with reconc_MixCond
+  res.mixCond <- reconc_MixCond(rec_mat$A, fc_bottom, fc_upper, bottom_in_type = "params", distr = 'poisson')
+
+  # We can reconcile with reconc_TDcond
+  res.TDcond <- reconc_TDcond(rec_mat$A, fc_bottom, fc_upper, bottom_in_type = "params", distr = 'poisson')
+
+  # Summary of the upper reconciled with TDcond 
+  pmfSum <- PMF.summary(res.TDcond$upper_reconciled$pmf[[1]])
+  # We expect that the reconciled mean is very similar to the initial mean (should be equal)
+  expect_equal(pmfSum$Mean,fc_upper$mu,tolerance = 0.01)
+
+  # Check that all bottom and upper reconciled PMF sum to 1
+  check_pmf_bott_mixCond <- sum(unlist(lapply(res.mixCond$bottom_reconciled$pmf, function(x){sum(x)}))) 
+  check_pmf_upp_mixCond <- sum(unlist(lapply(res.mixCond$upper_reconciled$pmf, function(x){sum(x)}))) 
+  expect_equal(check_pmf_bott_mixCond,12)
+  expect_equal(check_pmf_upp_mixCond,1)
+
+  # Check that all bottom and upper reconciled PMF sum to 1
+  check_pmf_bott_TDcond <- sum(unlist(lapply(res.TDcond$bottom_reconciled$pmf, function(x){sum(x)}))) 
+  check_pmf_upp_TDcond <- sum(unlist(lapply(res.TDcond$upper_reconciled$pmf, function(x){sum(x)}))) 
+  expect_equal(check_pmf_bott_TDcond,12)
+  expect_equal(check_pmf_upp_TDcond,1)
 
+
 })

tests/testthat/test-sample_funs.R

@@ -114,11 +114,11 @@ test_that("MVN density works", {
                  2,3,4,
                  0.5,0.5,0.5,
                  0,1,-1), ncol=3,byrow=TRUE)
-
+  
   res <- .MVN_density(x=xx,mu=mu,Sigma=Sigma)
 
   true_val <- c(8.742644e-04, 1.375497e-11, 3.739985e-03, 1.306453e-01)
-  expect_equal(res,true_val,tolerance = 0.1)
+  expect_equal(res,true_val, tolerance = 1e-6)
 
   # Check if block-evaluation works
   xx <- matrix(runif(3*1e4),ncol=3,byrow=TRUE)