TeX 55.1%
Python 44.2%
1# The HMTS method (Ishaak, Ouwehand & Remøy 2024) — step-by-step description for benchmark G23**Source (primary).** Ishaak, F., Ouwehand, P. & Remøy, H. (2024). "Constructing Limited-Revisable and4Stable CPPIs for Small Domains." *Journal of Official Statistics* 40(3), 380–408.5DOI 10.1177/0282423X241246617. CC BY-NC. Full text obtained on 2026-10-21 from the TU Delft repository6(final published version, `https://repository.tudelft.nl/file/File_f298b447-720b-44eb-84ab-8437ee06d026`);7SAGE blocks automated access. Page numbers below are the journal page numbers printed in that PDF.8Received 9 March 2023, accepted 18 December 2023 (p. 407).910**Source (implementation).** CRAN package `REPS` 1.1.1 (Statistics Netherlands; authors include Ishaak and11Ouwehand; maintainer V. Gajadhar), files `R/calculate_hmts.R`, `R/hmts_helpers.R`, `R/ts_helpers.R`,12obtained from `https://cran.r-project.org/src/contrib/REPS_1.1.1.tar.gz` on 2026-10-21. The package13cites the JOS 2024 article as the source of HMTS and (in roxygen comments) an "Ishaak, Ouwehand, Remoy14& De Haan (2023)" version that could not be located bibliographically.1516Everything below is paraphrased or quoted from these two sources. Where a detail is not stated, this is17said explicitly.1819---2021## 0. Purpose and context (pp. 380–386)2223- Target: a constant-quality commercial property price index (CPPI) for *small domains* with24 heterogeneous objects, for official statistics. The authors state that hedonic time-dummy pooling fixes25 small samples but "involves backward revision", and that the rolling time dummy method (RTD) "partly26 resolves this issue by pooling the data in a specified time window" but applied to commercial property27 "still often results in a volatile index" (p. 381).28- Three *practical* properties sought (§2.2, pp. 384–386): (1) a stable index ("not to remove market29 volatility, but only to remove volatility that is intrinsic to the index method or non-existing30 volatility", p. 384); (2) limited revision — they note EU member states "are by legislation bound to31 only one preliminary period" for the HPI (p. 385); (3) early detection of turning points, defined as a32 "structural" change in price development from positive to negative or vice versa (p. 385).33- Three *methodological* properties examined (§2.1, pp. 383–384): identity test, time-reversal test,34 circularity test.35- Acronym: **HMTS = "Hedonic Multilateral Time series Splice"** (p. 395). The package expands it as36 "Hedonic Multilateral Time series re-estimation Splicing".3738## 1. Overview: a four-step procedure (§3, p. 386)3940"(1) hedonic imputations, (2) multilateral imputations, (3) time series re-estimated imputations, and41(4) window splicing." The authors write "This is the first study that combines above techniques to42construct a price index." Figure 1 (p. 386) is a schematic of the four steps.4344## 2. Step 1 — Hedonic imputations (§3.1, pp. 386–388)4546**(a) Observation/hedonic equation.** A log-linear hedonic regression estimated **separately for each47period** t by OLS (Equation (1), p. 387, "based on Eurostat (2013b, 50)"):4849 ln p_{t,i} = α_t + Σ_{k=1}^{K} β_t^k x_{t,i}^k + ε_{t,i}, ε_{t,i} ~ N(0, σ²)5051with p_{t,i} the price of property i sold in period t, α_t the period-t intercept, β_t^k the period-t52coefficient of characteristic k, x_{t,i}^k the characteristic, ε_{t,i} the error.5354The paper then recalls Laspeyres and Paasche hedonic double-imputation (HDI) indexes, Equations (2)–(3)55(p. 387), and defines the *average imputed price* of period t with base period b, Equation (4) (p. 389):5657 p̂_{t(b)} = exp(α̂_t) · exp( Σ_{k=1}^{K} β̂_t^k · x̄_b^k )5859i.e. the period-t coefficients applied to the **mean characteristics of base period b**60(x̄_b^k = Σ_i x_{b,i}^k / N_b, p. 388). In words: "the mix of sold real estate in period [b] is tracked61through time. Its prices are imputed in the next periods as if they were sold in another period" (p. 388).62Only the imputed prices are used; "Equations (2) (Laspeyres) and (3) (Paasche) are not performed to63construct the HMTS index" (p. 388).6465**Implementation detail (REPS `calculate_hedonic_imputationmatrix`, hmts_helpers.R):** the model is66`log(price) ~ numerical_variables + categorical_variables`, fitted with `stats::lm.fit` on the period-t67rows (`X_dyn`, `y_dyn`); coefficients that are `NA` are set to 0; the imputation is `X_base %*% fitcoef`68on the *individual* rows of the base period b and the series value is `exp(mean(preds))` — i.e. the69geometric mean of the individual imputed prices of base-period properties, which equals the exp of the70fitted value at the base-period means for a linear-in-parameters model. Rows with missing values in any71model variable are dropped (`na.omit`). Factor levels are fixed on the whole dataset so that the design72matrix is comparable across periods. Numerical variables are *not* log-transformed by the package; the73vignette recommends the user log floor area beforehand. In the paper's application the dependent variable74and the numerical regressors (floor area, distance to train station) are in natural logs (Appendix 175note, p. 407).7677**What the paper does *not* say:** no minimum number of observations per period is given for Step 1; no78treatment of periods where OLS is rank-deficient is specified beyond the package's `NA → 0` rule. The79authors state as a limitation that "in cases where there are (almost) no observations, the HMTS index in80its current form does not work. Each period should have a few observations with explanatory variables81that show variation" and suggest, untested, a "multiperiod regression" fallback (§5.1, p. 404).8283## 3. Step 2 — Multilateral imputations (§3.2, pp. 388–390)8485Every reporting period t is imputed with **every possible base period b**, filling the full T×T matrix86of p̂_{t(b)} (Table 3, p. 390: rows = base period, columns = reporting period). "This results in multiple87series and all of these cover the complete time span, hence the term 'multilateral'" (p. 389). The idea88is "inspired by the multilateral GEKS ... method as described by Willenborg (2017, 2018)" (p. 389).8990A single period-t imputation would be the geometric mean over base periods, Equation (5), p. 389:9192 P̄_t = ( Π_{b=0}^{T} p̂_{t(b)} )^{1/T}, t = 0, 1, …, T9394and a circular index Equation (6): I_{0,t} = P̄_t / P̄_0. **But** "Equations (5) and (6) are not performed95at this stage for the final index. The imputations resulting from Equation (4) are used in the next step"96(p. 389).9798**Implementation detail:** in REPS the outer loop is over base periods (`current_period`), the inner loop99over reporting periods; for a base period b the series p̂_{·(b)} spans reporting periods 1…T (or, in100production mode, only up to the periods allowed by the update scheme — see §5 below).101102## 4. Step 3 — Time-series re-estimation of each base-period series (§3.3, pp. 390–392)103104**(b) The state-space model.** Each base-period series p̂_{t(b)}, t = 0…T, is treated as a noisy105observed series and smoothed with a univariate structural time-series model. The paper first presents the106random-walk-plus-noise (local level) model, Equations (7)–(8), p. 391 (citing Commandeur & Koopman 2007,107p. 9), then states (p. 391): "many models were tested and the best model was chosen by looking at the most108suitable fit. The model showed, among which, the lowest AIC (presented in Appendix 2). **The model of109choice is a local linear trend model.** ... Furthermore, the model sets a hyperparameter of the level110component equal to zero, while the slope component remains stochastic. This model is also known as a111**smooth trend model**, and since the slope disturbances affect the level component indirectly, generally112yield a smoother estimate of the level component."113114Equations (9)–(11), p. 392 (per base period b):115116 p̂_{t(b)} = μ_{t(b)} + ε_{t(b)}, ε_{t(b)} ~ NID(0, σ_ε²) (9) observation117 μ_{t+1(b)} = μ_{t(b)} + ν_{t(b)} (10) level (no disturbance)118 ν_{t+1(b)} = ν_{t(b)} + ζ_{t(b)}, ζ_{t(b)} ~ NID(0, σ_ζ²) (11) slope119120with μ_{t(b)} the unobserved level and ν_{t(b)} the slope. The paper writes "μ_{t(b)} = the unobserved121level in period t for the series with base period b; ν_{t(b)} = the slope component in period t122(deterministic)". So there are **two hyperparameters per series: σ_ε² (observation noise) and σ_ζ²123(slope noise); the level disturbance variance is fixed at 0.**124125**Hyperparameters:** "The model can be estimated by applying a two-step procedure. First, the126hyperparameters, that is, the variances of the error terms, are estimated via a Maximum Likelihood127procedure. Then the Kalman filter (Kalman 1960) is applied to estimate the state variables (level and128slope)" (p. 392). Filtering is forward, then "a backward pass is applied ... This is called smoothing. In129our case we are interested in these smoothed estimates" (p. 392). The re-estimated imputation is the130smoothed level, Equation (12), p. 392:131132 p̃_{t(b)} = p̂_{t(b)} − ε̂_{t(b)} ( = smoothed μ̂_{t(b)} )133134**Pooling across domains / across base-period series:** **none.** Each base-period series of each domain135is re-estimated independently with its own ML hyperparameters (one model fit per base period; Appendix 2,136p. 408, reports log-likelihood and AIC *per base period*). There is no cross-domain or cross-series137hierarchy; the "borrowing" in HMTS is purely temporal (adjacent periods of the same series) plus the138averaging across base periods in Step 4.139140**Implementation detail (REPS `ts_helpers.R`, exact):**141- The series is **log-transformed** before modelling (`original_series <- log(original_series)`), so the142 state-space model is fitted to ln p̂_{t(b)}, and the trend is exponentiated back (`exp(signalsubconf[,1])`).143 (The paper's Equations (9)–(12) are written in price levels; the package works in logs.)144- Model: `KFAS::SSModel(series ~ SSMtrend(2, Q = list(0, matrix(NA))), H = matrix(NA))` — a degree-2 trend145 (local linear trend) with **level variance fixed at 0** and **slope variance NA (estimated)**, observation146 variance `H` estimated. This matches the paper's smooth-trend model.147- Starting values: all hyperparameters start at `log(0.1)`; a one-pass Kalman run computes a "scale148 factor" `(1/(n−d)) Σ v_t' F_t^{-1} v_t` (d = number of diffuse initial states) and the log-start values149 are shifted by `0.5·log(scale)` (citing Commandeur & Koopman eq. 11.7 / SsfPack manual 9.4).150- Estimation: `KFAS::fitSSM(..., method = "BFGS", hessian = TRUE)` (maximum likelihood); then151 `KFAS::KFS(filtering = "state", smoothing = "state")`; the returned trend is the **smoothed level152 signal** (`KFAS::signal(states = c("level","regression"))`), with 95% bands `± 1.96·sqrt(var)` computed153 but not used for the index.154- Diagnostics stored per base period: log-likelihood, log-likelihood/n, AIC, AIC/n, number of diffuse155 states q, number of estimated variances w, σ and σ² for "meas" and "slope" with t-values, and first/last156 period level and slope estimates with standard errors. (Appendix 2 of the paper, p. 408, prints exactly157 these LLM and AIC values per base period for 2008Q1–Q4 and 2022Q1–Q4; e.g. 2022Q4: LL 84.3, AIC −160.6,158 n = 60.)159- The roxygen comments say the time-series step "can be turned off with a parameter" (and mention160 `state_space_model = NULL`), but in REPS 1.1.1 as shipped `calculate_hedonic_imputationmatrix()` calls161 `calculate_trend_line_kfas()` unconditionally and `calculate_hmts()` exposes no such switch; the162 `imputation` argument of `calculate_hedonic_index()` applies to Laspeyres/Paasche only. The HMTS-relevant163 arguments of `calculate_hedonic_index(method = "hmts")` are `number_preliminary_periods` (package default164 3 if omitted), `production_since` (default `NULL`) and `resting_points`; `window_length` is the165 rolling-time-dummy argument (package default 5), not an HMTS argument.166167**What is not stated:** no seasonal component is used (start-value slots for "seas" and "scaling" exist168in the code but the selected model has none); no diffuse-initialisation details beyond KFAS defaults; no169statement on how many periods are needed for ML estimation to be reliable (the empirical application170starts the state-space re-estimation "as if 2015 Q1 was the first moment of index construction" with17128 quarters available, Appendix 2 note, p. 408, and Table 5 footnote says the scheme "is implemented in172practice with around twenty starting periods", p. 394).173174## 5. Step 4 — Window splice (§3.4, pp. 392–395)175176**(c) Why:** both the multilateral step (adding a reporting period adds a base-period series, changing177the geometric mean for *every* period) and the smoothing step (smoothed values use preceding *and*178succeeding periods) revise the whole series each period (p. 393). "The main idea of window splicing is179to keep a limited number of periods provisional instead of the entire time series" (p. 393).180181**Two-part rule (pp. 393–395):**1821831. **Update scheme for the re-estimated series (Table 5, p. 394).** Fix a window of w periods. When184 period t is first computed it is provisional; "In this period and the next two periods [for w = 3],185 the figure for period t is subject to revision. In period t+1, the series with base period t is186 extended with reporting period t+1. In period t+2, the series with base period t is extended again187 with reporting period t+2. In both calculations, the figure for period t is revised. This will be the188 last time, since the series with base period t is fixed from that point on" (p. 393). So **each189 base-period series is re-estimated only until it is w periods long, then frozen**; the number of190 provisional periods equals w − 1.1912. **Index transformation, Equation (17), p. 394:** within each base-period series, I_{0,t(b)} =192 p̃_{t(b)} / p̃_{0(b)}, so that series with different base periods become comparable.1933. **Window splice, Equation (18), p. 394:** the final index for reporting period t is the **geometric194 average of the w index values with base periods t−(w−1), …, t**:195196 I^{HMTS,w}_{0,t} = ( Π_{b = t−(w−1)}^{t} I_{0,t(b)} )^{1/w}, t = w−1, w, …, T197198 e.g. for w = 3, I_{0,t} is the geometric mean of I_{0,t(t−2)}, I_{0,t(t−1)}, I_{0,t(t)} (p. 394).199 For the first periods (t ≤ w−2) the window is truncated at the start, Equation (19), p. 395:200 I^{w}_{0,t} = ( Π_{b=0}^{w−1} I_{0,t(b)} )^{1/w}. Table 6 (p. 395) shows the three-period example.201202**Limited-revisable — how many periods revise:** with window length w, a period is revised in the w−1203publications following its first release; "The HMTS index with three preliminary periods (four period204window splice)" is the configuration used in the empirical section (§4.2, p. 396), i.e. **w = 4, three205preliminary periods, every period published four times** (first estimate, then after +1, +2, +3 periods;206§4.3, p. 397). In the Dutch office-market simulation around the 2015 turning point, first estimates were207on average 0.82 index points from the final, second estimates 0.29, third 0.03; "two or three preliminary208periods seem to be appropriate" (p. 397). Note: the paper's "three preliminary periods" is the package209parameter `number_preliminary_periods` (default 3 in the package doc; the vignette example uses 1 with210`window_length = 4`).211212**Choice of w:** not fixed by the paper. "The choice of window length is rather important ... greater213transitivity goes along with a longer window. On the other hand, a shorter window minimizes the number214of revisions" (p. 395); they refer to Hill et al. (2022) for a window-length criterion and announce "Another215way to look for an optimal window length—and specifically applied to the HMTS—is presented in section2164.2" (p. 395) — **in the published §4.2 the only window-related content is the statement above about217revisions by estimate number; no formal selection criterion for w is given.**218219**Implementation detail (REPS `calculate_hmts_index`):** `number_preliminary_periods` controls the220update scheme; when `production_since` is `NULL` all periods are treated as final (no splice — the221whole matrix is averaged: `window <- imputations[, -1]`). With `production_since` set, the package222rebuilds the series period by period ("To simulate a series, where 1 period a time expires (as in223production)"), keeping for each reporting period the index values of the base-period series that are224still within the window, then takes `calculate_geometric_average(na.omit(...))` per period and rescales to225`reference_period` = 100. The window splice is "a specified (parameter) window ... of index figures that226continues in the calculation" and "can be turned off".227228## 6. (d) What "small domains" means operationally229230- The paper gives **no numerical threshold** defining a small domain. The empirical domain is Dutch231 office-building transactions 2008–2022, quarterly, "approximately 150 to 1,100 per quarter" (§4.1,232 p. 396); a robustness run uses "a randomly selected subset of the data with fifty observations per233 quarter" (§4.3, p. 398), which still produced a usable, "somewhat less stable" index. Appendix 1234 (p. 407) shows n = 461 (2008Q1) and n = 785 (2022Q4) in the period regressions.235- The only stated lower bound is qualitative (§5.1, p. 404): each period needs "a few observations with236 explanatory variables that show variation", otherwise the method "does not work".237- Hedonic variables in the application (Appendix 1, p. 407): log floor area, log distance to train238 station, neighbourhood segment dummies (B–E), "In Amsterdam?" dummy, building-age class dummies;239 adjusted R² 0.85 (2008Q1) and 0.76 (2022Q4). The authors report that "the HMTS method performs well even240 with the most basic hedonic variables: floor area, location, and building age" (p. 396).241242## 7. Properties and evaluation reported (§4.3–4.5, pp. 397–403)243244- Bootstrap 95% intervals (500 resamples with replacement) show "improved reliability compared to a245 Fisher HDI" (p. 398). Leave-one-out influence: individual-observation shares range −0.4 to +0.2 index246 points for HMTS vs −0.67 to +0.27 for Fisher HDI; mean absolute share 0.03 vs 0.06 (pp. 398–399).247- Identity test: **not satisfied** (the index returns to 100 only slowly after prices are reset; "The248 risk is that [it] responds with a lag to extreme price changes", p. 400). Time reversal: satisfied to249 seven decimals *without* the splice, violated with it (p. 401). Circularity: lost after Steps 3–4, but250 indices with different start years "differ the most at the start of the series" and are otherwise251 similar (pp. 401–402).252- Versus RTD (window 3, movement splice, as in de Haan 2015): the RTD "still shows volatility similar to253 the Fisher HDI index" and "shows drift, caused by 'weak' periods in 2013"; longer RTD windows move it254 toward Fisher (pp. 402–403).255- Lead/lag: pairwise Granger tests between HMTS and (i) HDI with first base, (ii) multilateral HDI, (iii)256 Fisher HDI — "there is no clear lag of the HMTS. If any, the HMTS leads the other indices" (p. 403).257- Turning-point detection in the 2015 simulation: the first two estimates of the turning-point quarter258 pointed down, the third was corrected up — "the desired accuracy lags two periods" (p. 397).259260## 8. Checklist for implementing benchmark G faithfully (derived strictly from the above)2612621. Per period t: OLS of ln price on the hedonic variables (period-specific coefficients).2632. For every base period b and reporting period t: p̂_{t(b)} = exp(mean of X_b β̂_t) (log scale mean of264 base-period rows under period-t coefficients).2653. For every b: fit a smooth-trend model (level variance = 0, slope variance and observation variance by266 ML; KFAS `SSMtrend(2, Q=list(0,NA))`, `H=NA`, BFGS) to **ln p̂_{·(b)}**, take the smoothed level,267 exponentiate → p̃_{t(b)}.2684. Index each series on its own first period: I_{0,t(b)} = p̃_{t(b)} / p̃_{0(b)}.2695. Production/vintage mode: a base-period series b is extended only while it has < w periods; afterwards270 frozen. Final index for t = geometric mean of I_{0,t(b)} for b = t−(w−1)…t (truncated at the start).271 Paper configuration: w = 4 (three preliminary periods). Rescale to the reference period.2726. Report: revisions by estimate number (1st…wth), bootstrap CI (500 draws), leave-one-out influence,273 identity/time-reversal/circularity checks, Granger lead/lag vs Fisher HDI — as in §4.274275**Open points that the sources do not settle (must be chosen and documented by us):** minimum276observations per period for Step 1; how to handle a period with no valid regression (paper suggests but277does not test a multi-period regression); the window length w (paper: 4 in the application, no278criterion); whether to run the state-space step on logs (package) or levels (paper equations) — the279package behaviour (logs) is what Statistics Netherlands ships; the number of "starting periods" before280the splice regime begins (paper: "around twenty", p. 394; 28 quarters in the application).281