Merge branch 'develop'

DESCRIPTION

@@ -1,21 +1,23 @@
 Package: rjd3bench
 Type: Package
 Title: Interface to 'JDemetra+ 3.x' time series analysis software
-Version: 2.0.1
+Version: 2.1.0
 Authors@R: c(
- person("Jean", "Palate", role = c("aut", "cre"),
+ person("Jean", "Palate", role = "aut",
  email = "jean.palate@nbb.be"),
+ person("Corentin", "Lemasson", role = c("cre", "ctb"),
+ email ="corentin.lemasson@nbb.be"),
  person("Tanguy", "Barthelemy", role = c("ctb", "art"),
  email ="tanguy.barthelemy@insee.fr"))
 Description: R Interface to 'JDemetra+ 3.x' (<https://github.com/jdemetra>) time series analysis software.
 Depends:
  R (>= 4.1.0)
 Imports:
  rJava (>= 1.0-6),
- rjd3toolkit (>= 3.2.2),
+ rjd3toolkit (>= 3.2.4),
  RProtoBuf (>= 0.4.20)
 Remotes: 
- github::rjdverse/rjd3toolkit@*release
+ github::rjdverse/rjd3toolkit
 SystemRequirements: Java (>= 17)
 License: EUPL
 URL: https://github.com/rjdverse/rjd3bench,

NEWS.md

@@ -7,13 +7,21 @@ to [Semantic Versioning](https://semver.org/spec/v2.0.0.html).
 
 ## [Unreleased]
 
-## [2.0.1] - 2024-07-12
 
+## [2.1.0] - 2024-07-18
+
+### Added
+
+* Add output on residuals in temporaldisaggregation() function
+
+
+## [2.0.1] - 2024-07-12
 
 ### Changed 
 
 *new jars related to version [1.2.1](https://github.com/jdemetra/jdplus-benchmarking/releases/tag/v1.2.1)*
 
+
 ## [2.0.0] - 2023-12-12
 
 ### Added
@@ -28,7 +36,8 @@ to [Semantic Versioning](https://semver.org/spec/v2.0.0.html).
 * Initial commit
 
 
-[Unreleased]: https://github.com/rjdverse/rjd3bench/compare/v2.0.1...HEAD
+[Unreleased]: https://github.com/rjdverse/rjd3bench/compare/v2.1.0...HEAD
+[2.1.0]: https://github.com/rjdverse/rjd3bench/compare/v2.0.1...v2.1.0
 [2.0.1]: https://github.com/rjdverse/rjd3bench/compare/v2.0.0...v2.0.1
 [2.0.0]: https://github.com/rjdverse/rjd3bench/releases/tag/v1.0.0...2.0.0
 [1.0.0]: https://github.com/rjdverse/rjd3bench/releases/tag/v1.0.0

R/tempdisagg.R

@@ -60,11 +60,11 @@ temporaldisaggregation<-function(series, constant = TRUE, trend = FALSE, indic
  }
  } else if (is.ts(indicators)){
  jlist[[1]]<-rjd3toolkit::.r2jd_tsdata(indicators)
- } else {
+ } else{
  stop("Invalid indicators")
  }
  jindicators<-.jarray(jlist, contents.class = "jdplus/toolkit/base/api/timeseries/TsData")
- } else {
+ } else{
  jindicators<-.jnull("[Ljdplus/toolkit/base/api/timeseries/TsData;")
  }
  jrslt<-.jcall("jdplus/benchmarking/base/r/TemporalDisaggregation", "Ljdplus/benchmarking/base/core/univariate/TemporalDisaggregationResults;",
@@ -92,8 +92,8 @@ temporaldisaggregation<-function(series, constant = TRUE, trend = FALSE, indic
  regeffect=rjd3toolkit::.proc_ts(jrslt, "regeffect"),
  smoothingpart=rjd3toolkit::.proc_numeric(jrslt, "smoothingpart"),
  parameter=rjd3toolkit::.proc_numeric(jrslt, "parameter"),
- eparameter=rjd3toolkit::.proc_numeric(jrslt, "eparameter")
- # res= TODO
+ eparameter=rjd3toolkit::.proc_numeric(jrslt, "eparameter"),
+ residuals= .proc_residuals(jrslt) # temporary solution (see function below)
  )
  likelihood<-rjd3toolkit::.proc_likelihood(jrslt, "likelihood.")
 
@@ -180,7 +180,7 @@ temporaldisaggregationI<-function(series, indicator,
 print.JD3TempDisagg<-function(x, ...){
  if (is.null(x$regression$model)){
  cat("Invalid estimation")
- } else {
+ } else{
  cat("Model:", x$regression$type, "\n")
  print(x$regression$model)
 
@@ -205,7 +205,7 @@ print.JD3TempDisagg<-function(x, ...){
 print.JD3TempDisaggI<-function(x, ...){
  if (is.null(x$estimation$parameter)){
  cat("Invalid estimation")
- } else {
+ } else{
  model<-data.frame(coef = c(round(x$regression$a, 4), round(x$regression$b, 4)))
  row.names(model)<-c("a", "b")
  print(model)
@@ -254,7 +254,7 @@ summary_disagg<-function(object){
  if (is.null(object)){
  cat("Invalid estimation")
 
- } else {
+ } else{
  cat("\n")
  cat("Likelihood statistics", "\n")
  cat("\n")
@@ -298,7 +298,7 @@ summary.JD3TempDisaggI<-function(object, ...){
  if (is.null(object)){
  cat("Invalid estimation")
 
- } else {
+ } else{
  cat("\n")
  cat("Likelihood statistics", "\n")
  cat("\n")
@@ -336,7 +336,7 @@ plot.JD3TempDisagg<-function(x, ...){
  if (is.null(x)){
  cat("Invalid estimation")
 
- } else {
+ } else{
  td_series <- x$estimation$disagg
  reg_effect <- x$estimation$regeffect
  smoothing_effect <- td_series - reg_effect
@@ -364,9 +364,51 @@ plot.JD3TempDisaggI<-function(x, ...){
  if (is.null(x)){
  cat("Invalid estimation")
 
- } else {
+ } else{
  td_series <- x$estimation$disagg
  ts.plot(td_series, gpars=list(xlab="", ylab="disaggragated series", xaxt="n"))
  axis(side=1, at=start(td_series)[1]:end(td_series)[1])
  }
 }
+
+# TEMPORARY SOLUTION
+# For the next release, we should use proto and move the functions to rjd3toolkit
+.proc_residuals <- function (jrslt){
+
+ z<-rjd3toolkit::.jd3_object(jrslt, "TD", TRUE)
+
+ full_residuals <- get_result_item(z,"residuals.fullresiduals")
+
+ extr_normality <- list(get_result_item(z,"residuals.mean"),
+ get_result_item(z,"residuals.skewness"),
+ get_result_item(z,"residuals.kurtosis"),
+ get_result_item(z,"residuals.doornikhansen"))
+ extr_independence <- get_result_item(z,"residuals.ljungbox")
+ extr_randomness <- list(get_result_item(z,"residuals.nruns"),
+ get_result_item(z,"residuals.lruns"),
+ get_result_item(z,"residuals.nudruns"),
+ get_result_item(z,"residuals.ludruns"))
+ linearity_test <- tryCatch(rjd3toolkit::ljungbox(full_residuals^2, k = 8, lag = 1, mean = TRUE), error=function(err) NaN)
+
+
+ normality <- matrix(unlist(extr_normality), nrow = 4, ncol = 2, byrow = TRUE,
+ dimnames = list(c("mean", "skewness", "kurtosis", "test(doornikhansen)"), c("value", "p-value")))
+ independence <- matrix(unlist(extr_independence), nrow = 1, ncol = 2, byrow = TRUE,
+ dimnames = list(c("ljung_box"), c("value", "p-value")))
+ randomness <- matrix(unlist(extr_randomness), nrow = 4, ncol = 2, byrow = TRUE,
+ dimnames = list(c("Runs around the mean: number", "Runs around the mean: length", "Up and Down runs: number", "Up and Down runs: length"), c("value", "p-value")))
+ linearity <- matrix(unlist(linearity_test), nrow = 1, ncol = 2, byrow = TRUE,
+ dimnames = list(c("ljung_box on squared residuals"), c("value", "p-value")))
+
+ return(list(full_residuals=full_residuals,
+ tests=list(normality=round(normality,4),
+ independence=round(independence,4),
+ randomness=round(randomness,4),
+ linearity=round(linearity,4))))
+
+}
+
+get_result_item <- function(jd3_obj, item){
+ return(tryCatch(rjd3toolkit::result(jd3_obj, item),
+ error=function(err) NaN))
+}

vignettes/rjd3bench.Rmd

@@ -1,5 +1,5 @@
 ---
-title: "Temporal disaggregation and Benchmarking methods based on JDemetra+"
+title: "Temporal disaggregation and Benchmarking methods based on JDemetra+ v3.x"
 output: 
  html_vignette:
  toc: true
@@ -8,29 +8,27 @@ output:
  toc: true
  toc_depth: 3
 vignette: >
- %\VignetteIndexEntry{Part of the interface to 'JDemetra+ 3.0' Seasonal adjustment software related to Temporal disaggregation and Benchmarking methods}
+ %\VignetteIndexEntry{Temporal disaggregation and Benchmarking methods based on JDemetra+ v3.x}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
 abstract: The package **rjd3bench** provides a variety of methods for temporal disaggregation, benchmarking, reconciliation and calendarization. It is part of the interface to 'JDemetra+ 3.0' Seasonal adjustement software. Methods of temporal disaggregation and benchmarking are used to derive high frequency time series from low frequency time series with or without the help of high frequency information. Consistency of the high frequency series with the low frequency series can be achieved by either the sum, the average, the first or last value or any other user-defined conversion mode. In addition to temporal constraints, reconciliation methods deals with contemporaneous consistency while adjusting multiple time series. Finally, calendarization method can be used when time series data do not coincide with calendar periods. 
 ---
 
-```{r, include = FALSE}
+```{r setup vignette, include = FALSE}
 knitr::opts_chunk$set(
  collapse = TRUE,
+ eval = FALSE,
  comment = "#>"
 )
 ```
-\
- 
+
 
 # Introduction
 
 The methods implemented in the package **rjd3bench** intend to bridge the gap when there is a lack of high frequency time series or when there are temporal and/or contemporaneous inconsistencies between the high frequency series and the corresponding low frequency series, the latter being the benchmark. Although this can be an issue in any fields of research dealing with time series, methods of temporal disaggregation, benchmarking, reconciliation and calendarization are mostly used in official statistics, especially in the production of macroeconomic statistics. For example, National Accounts are often compiled according to two frequencies of production: annual series, the low frequency data, based on precise and detailed sources and quarterly series, the high frequency data, which usually rely on less accurate sources but give information on a timelier basis. The use of temporal disaggregation, benchmarking, reconciliation or calendarization methods allow to achieve consistency between annual and quarterly national accounts over time. 
 
 The package is an R interface to the highly efficient algorithms and modeling developed in the official 'JDemetra+ 3.0' Seasonal adjustement software. It provides a wide variety of methods and intend to include those suggested in the *ESS guidelines on temporal disaggregation, benchmarking and reconciliation (Eurostat, 2018)*. 
 
-\
- 
 
 # Data
 
@@ -39,17 +37,14 @@ We illustrate the various methods using two sets of time series:
 * The **retail** data set contains monthly figures over retail activity of various categories of goods and services from 1992 to 2010.
 * The **qna_data** contains a list of two data sets. The first data set 'B1G_Y_data' includes three annual benchmark series which are the Belgian annual value added on the period 2009-2020 in chemical industry (CE), construction (FF) and transport services (HH). The second data set 'TURN_Q_data' includes the corresponding quarterly indicators which are (modified) production indicators derived from VAT statistics and covering the period 2009Q1-2021Q4.
 
-```{r, data_used}
+```{r data, echo = TRUE, eval = FALSE}
 library("rjd3toolkit") # for **retail** data
 data("qna_data")
 ```
 
-\
- 
-
 # Application of methods
 
-```{r setup}
+```{r package loading, echo = TRUE, eval = FALSE}
 library("rjd3bench")
 ```
 
@@ -60,15 +55,15 @@ library("rjd3bench")
 Eurostat (2018) recommends the use of regression-based models for the purpose of temporal disaggregation. Among them, Chow-Lin method and its variants proposed by Fernandez and Litterman are largely used in practice. Those can be called through the **temporaldisaggregation()** function in the **rjd3bench** package.
 
 Here is an example of code to use Fernandez variant as method to disaggregate annual value added in construction sector using a quarterly production index as indicator.
-```{r, chow_lin_example}
+```{r chow_lin example, echo = TRUE, eval = FALSE}
 Y<-ts(qna_data$B1G_Y_data[, "B1G_FF"], frequency=1, start=c(2009, 1))
 x<-ts(qna_data$TURN_Q_data[, "TURN_INDEX_FF"], frequency=4, start=c(2009, 1))
 td_cl<-rjd3bench::temporaldisaggregation(Y, indicators=x, model = "Rw")
 y_cl<-td_cl$estimation$disagg # the disaggregated series
 ```
 
 The variable *td_cl* is of class **'JD3TempDisagg'** and contains the most important information about the regression including the estimates of model coefficients and their covariance matrix, the decomposition of the disaggregated series and information about the likelihood. There exist a specific print(), summary() and plot() functions for objects of class **'JD3TempDisagg'**. print() shows the fitted model, summary() also shows extra information about the model and the likelihood. plot() displays the decomposition of the disaggregated series between regression and the smoothing effect. 
-```{r, chow_lin_example_cont, fig.width=6}
+```{r chow_lin example cont, echo = TRUE, eval = FALSE}
 summary(td_cl) # example of call to summary() 
 plot(td_cl)
 ```
@@ -97,7 +92,7 @@ where the annual constraint are taken care of by considering a cumulated series
 This alternative representation of Denton PFD method is interesting as it allows more flexibility in the model such as the inclusion of outliers (level shift in the Benchmark to Indicator ratio) - that could otherwise induce unintended wave effects - and the possibility to fix the BI ratio (and therefore also the disaggregated) at some specific periods. Outliers and their intensity are defined by changing the value of the 'innovation variances'. Following the principle of movement preservation inherent to Denton, the model-based Denton PFD method constitutes an interesting alternative for temporal disaggregation, giving a thorough analysis of the data (and the Benchmark-to-Indicator (BI) ratio in particular) and a clear strategy for extrapolation.
 
 Model-based denton can be called through the **denton.modelbased()** function. Here is an example of code to apply it to the same data as before.
-```{r, mb_denton_example, fig.height=6, fig.width=6}
+```{r mb_denton example, echo = TRUE, eval = FALSE}
 td_mbd<-rjd3bench::denton_modelbased(Y, x, outliers = list("2020-01-01"=100, "2020-04-01"=100))
 y_mbd<-td_mbd$estimation$disagg
 plot(td_mbd)
@@ -113,7 +108,7 @@ Denton method relies on the principle of movement preservation. There exist a fe
 preservation: additive first difference (AFD), proportional first difference (PFD), additive second difference (ASD), proportional second difference (PSD), etc. The different variants can be called through the **denton()** function in the **rjd3bench** package. 
 
 Here are a few examples using the table 'B1G_Y_data' from **qna_data**:
-```{r, denton_example}
+```{r denton example, echo = TRUE, eval = FALSE}
 Y<-ts(qna_data$B1G_Y_data[, "B1G_HH"], frequency=1, start=c(2009, 1))
 
 y_d1<-rjd3bench::denton(t=Y, nfreq=4) # example of denton PFD without high frequency series
@@ -129,7 +124,7 @@ The **denton()** function returns the high frequency series benchmarked with den
 This method corresponds to the method of Cauley and Trager, using the solution proposed by Di Fonzo and Marini. It can be called through the **grp()** function.
 
 Here is an example on how to use it:
-```{r, grp_example}
+```{r grp example, echo = TRUE, eval = FALSE}
 y_grp<-rjd3bench::grp(s=x, t=Y)
 ```
 
@@ -140,7 +135,7 @@ The **grp()** function returns the high frequency series benchmarked with grp me
 Cubic splines are piecewise cubic functions that are linked together in a way to guarantee smoothness at data points. Additivity constraints are added for benchmarking purpose and sub-period estimates are derived from each spline. When a sub-period indicator (or disaggregated series) is used, cubic splines are no longer drawn based on the low frequency data but the Benchmark-to-Indicator (BI ratio) is the one being smoothed. Sub-period estimates are then simply the product between the smoothed high frequency BI ratio and the indicator.
 
 The method can be called through the **cubicspline()** function. Here are a few examples on how to use it:
-```{r, cubicspline_example}
+```{r cubicspline example, echo = TRUE, eval = FALSE}
 y_cs1<-rjd3bench::cubicspline(t=Y, nfreq=4) # example of cubic spline without high frequency series (smoothing)
 
 x<-y_cs1+rnorm(n=length(y_cs1), mean=0, sd=10)