L0TFinv
is a toolkit for L0-regularized sparse
approximation, which is used for change point detection. The package
excels in both efficiency and accuracy of trend estimation and
changepoint detection in segmented functions.
Trend filtering is a typical method for nonparametric regression, including change point detection and time series segmentation. The commonly used trend filtering models is the L1 trend filtering model \((a)\) based on the difference matrix \(\boldsymbol{D}^{(q+1)}\), as illustrated below.
\[ \min _{\boldsymbol{\alpha} \in \mathbb{R}^n} \frac{1}{2}\|\boldsymbol{y}-\boldsymbol{\alpha}\|_2^2 + \lambda\|\boldsymbol{D}^{(q+1)} \boldsymbol{\alpha}\|_{\ell_1}, \quad q=0,1,2, \ldots. \quad (a) \] In previous studies, algorithms solving trend filtering problems \((a)\) necessitate the computation of \(((\boldsymbol{D}^{(q+1)})^T \boldsymbol{D}^{(q+1)})^{-1}\). When \(n\) is large, just fitting the matrix into memory becomes an issue. L0 trend filtering \((b)\) has a advantage over other trend filtering methods, especially in the detection of change points. The expression for L0 trend filtering is as follows:
\[ \min _{\boldsymbol{\alpha} \in \mathbb{R}^n} \frac{1}{2}\|\boldsymbol{y}-\boldsymbol{\alpha}\|_2^2 + \lambda\|\boldsymbol{D}^{(q+1)} \boldsymbol{\alpha}\|_{\ell_0}. \quad (b) \] In L0 trend filtering \((b)\), the positions of non-zero elements in the L0 norm correspond with the locations of change points. We consider two subsets: the active set \(A\) for non-zero elements and the inactive set \(I\) for zero elements. Despite this, computing \(((\boldsymbol{D}^{(q+1)}_I)^T \boldsymbol{D}^{(q+1)}_I)^{-1}\) remains a task involving a substantial matrix. We explore transforming the problem \((b)\) into a L0-regularized sparse format \((c)\) by introducing an artificial design matrix \(\boldsymbol{X}^{(q+1)}\) that corresponds to the difference matrix, thereby reformulating the L0 trend filtering problem into the following format.
\[ \min _{\boldsymbol{\beta} \in \mathbb{R}^n} \frac{1}{2}\|\boldsymbol{y}-\boldsymbol{X}^{(q+1)}\boldsymbol{\beta}\|_2^2 + \lambda \sum_{i=q+2}^n |\boldsymbol{\beta}_i|_{\ell_0}. \quad (c) \] Due to the connection between L0 constraint problems and L0 penalty problems, and considering that the sparsity of \(\boldsymbol{\beta}\) is is more meaningful in practical applications than the selection of the hyperparameter \(\lambda\). We focus on the constraint that reflects our aim to achieve an estimated trend with a given number of change points. So we transform the L0 penalty problem \((c)\) into the L0 constraint problem \((d)\).
In our practical approach, we consider the maximum number of change points \(k_{\text{max}}\) as a constraint, transforming the aforementioned L0 penalty problem \((c)\) into the following L0 constraint problem.
\[ \text{ minimize }\frac{1}{2}\|\boldsymbol{y}-\boldsymbol{X}^{(q+1)}\boldsymbol{\beta}\|_2^2,\quad \text{ subject to } \sum_{i=q+2}^n |\boldsymbol{\beta}_i|_{\ell_0} \leq k_{\text{max}}. \quad (d) \]
As previously discussed, the computation of \(((\boldsymbol{D}^{(q+1)}_I)^T \boldsymbol{D}^{(q+1)}_I)^{-1}\) remains a substantial challenge, particularly under constraints of sparsity where \(|A| \ll n\), resulting in a very large matrix. However, when transforming the problem \((b)\) into the problem \((d)\), we only need to calculate a relatively small matrix \(((\boldsymbol{X}^{(q+1)}_A)^T \boldsymbol{X}^{(q+1)}_A)^{-1}\), which involves computational simplifications. For such L0 constraint problems \((d)\), we employ a splicing-based approach to design algorithms for processing. Many computational tricks are used to speed up the algorithms and improve the solution quality including warm starts, the active set search, splicing methods and the simplification of \(((\boldsymbol{X}^{(q+1)}_A)^T \boldsymbol{X}^{(q+1)}_A)^{-1}\).
In this vignette, we provide a tutorial on using the R interface. Particularly, we will demonstrate how to use L0TFinv’s main functions for fitting models, change point detection, trend estimation and visualization.
L0TFinv can be installed by executing:
Please proceed to install and load the subsequent R packages.
L0TFinv
package is well suited to analyzing time series
with underlying piecewise constant or piecewise linear trends. First we
load L0TFinv:
We will start by fitting a simple L0-regularized sparse model and then move on to present the results of change point detection and trend estimation.
Before demonstrating how L0TFinv works, we will first describe the structure matrix \(\boldsymbol{D}^{(q+1)}\) for the original problem \((b)\) and the matrix \(\boldsymbol{X}^{(q+1)}\) for the transformed problem \((d)\).
The elements of \(\boldsymbol{D}^{(q+1)}\) are related to
combination numbers, while the elements of \(\boldsymbol{X}^{(q+1)}\) are associated
with the cumulative sum function. Additionally, \(\boldsymbol{D}^{(q+1)}\boldsymbol{X}^{(q+1)}\)
has a special structure, details of which can be found in
XMat
function. Another special structure matrix is denoted
as \((\boldsymbol{X}^{(q+1)}_A)^T
\boldsymbol{X}^{(q+1)}_A\) and the solMat function simplifies the
computation of this matrix inverse ( \(A\) is the active set ), commonly employed
in splicing algorithms. For more details, see the solMat
function.
#> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
#> [1,] 1 -2 1 0 0 0 0 0 0 0
#> [2,] 0 1 -2 1 0 0 0 0 0 0
#> [3,] 0 0 1 -2 1 0 0 0 0 0
#> [4,] 0 0 0 1 -2 1 0 0 0 0
#> [5,] 0 0 0 0 1 -2 1 0 0 0
#> [6,] 0 0 0 0 0 1 -2 1 0 0
#> [7,] 0 0 0 0 0 0 1 -2 1 0
#> [8,] 0 0 0 0 0 0 0 1 -2 1
#> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
#> [1,] 1 0 0 0 0 0 0 0 0 0
#> [2,] 2 1 0 0 0 0 0 0 0 0
#> [3,] 3 2 1 0 0 0 0 0 0 0
#> [4,] 4 3 2 1 0 0 0 0 0 0
#> [5,] 5 4 3 2 1 0 0 0 0 0
#> [6,] 6 5 4 3 2 1 0 0 0 0
#> [7,] 7 6 5 4 3 2 1 0 0 0
#> [8,] 8 7 6 5 4 3 2 1 0 0
#> [9,] 9 8 7 6 5 4 3 2 1 0
#> [10,] 10 9 8 7 6 5 4 3 2 1
#> [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
#> [1,] 0 0 1 0 0 0 0 0 0 0
#> [2,] 0 0 0 1 0 0 0 0 0 0
#> [3,] 0 0 0 0 1 0 0 0 0 0
#> [4,] 0 0 0 0 0 1 0 0 0 0
#> [5,] 0 0 0 0 0 0 1 0 0 0
#> [6,] 0 0 0 0 0 0 0 1 0 0
#> [7,] 0 0 0 0 0 0 0 0 1 0
#> [8,] 0 0 0 0 0 0 0 0 0 1
In this part, we generate two types of underlying trends for time series segmentation, piecewise constant or piecewise linear. The observation point \(\boldsymbol{y}\) is generated as follows:
\[ \boldsymbol{y}=\boldsymbol{X}^{(q+1)}\boldsymbol{\beta}+\boldsymbol{\epsilon},\quad q=0,1. \] The matrix \(\boldsymbol{X}^{(q+1)}\) is as previously defined, with \(\boldsymbol{\epsilon}\) representing iid standard normal entries, and \(\boldsymbol{\beta}\) being the coefficients that contain the information of the segmented functions. The constraints on the change points of the trend are transformed into sparsity constraints on the \(\boldsymbol{\beta}\) parameters. Most of the \(\boldsymbol{\beta}\)’s positional components are zero, with non-zero values only at the corresponding change points.
For both types of trends, we require inputs for \(n\) (the number of data points), \(sigma\) (the variance of the normal distribution), \(seed\) (to ensure reproducibility by fixing the seed), and \(tau\) (the positions of the change points, standardized to the interval [0, 1] for ease of processing).
The change points of \(length(tau)\) divides the interval \([0,1]\) into \(length(tau)+1\) segments. Additionally, we require a vector \(h\) of \(length(tau)+1\) to describe the information. For a piecewise constant trend, we need to input \(length(tau)+1\) constant values to represent the piecewise constant function, as illustrated below.
tau = c(0.2, 0.3, 0.5, 0.65, 0.85)
h = c(-1,3,-2,0,4,-3)
BlocksData <- SimuBlocksInv(n = 500, sigma = 0.1, seed = 50, tau = tau ,h = h)
plot(BlocksData$x, BlocksData$y, xlab="", ylab="") ## The piecewise linear simulated data
lines(BlocksData$x, BlocksData$y0, col = "red") ## The underlying trend
#> [1] 100 150 250 325 425
#> [1] 0.20 0.30 0.50 0.65 0.85
For a piecewise linear trend, input a vector \(h\) of \(length(tau)+1\) to represent the slope of each linear segment. It is also necessary to define the initial value, denoted as \(a_0\).
tau1 = c(0.1, 0.3, 0.4, 0.7, 0.9)
h1 = c(-2, 5, -3, 2, -1, 4)
a0 = -10
WaveData <- SimuWaveInv(n = 500, sigma = 0.1, seed = 50, tau = tau1, h = h1, a0 = a0)
plot(WaveData$x, WaveData$y, xlab="", ylab="")
lines(WaveData$x, WaveData$y0, col = "red")
#> [1] 50 150 200 350 450
#> [1] 0.1 0.3 0.4 0.7 0.9
There are two ways to fit models, fixed change points or optimal change points. Both methods require inputting data points, the order \(q\) ( piecewise constant or piecewise linear fitting ) and the maximum possible interval \([first , last]\) for change points to occur ( with the entire range being \([0 , 1]\) ).
For the case of the fixed change points, we need to additionally define the given number \(k\) of change points. Since it involves warm starts, the final results will include scenarios where the number of change points varies from 1 to \(k\).
FitBlocks.fix <- L0TFinv.fix(y=BlocksData$y, k=10, q=0, first=0.01, last=1)
FitWave.fix <- L0TFinv.fix(y=WaveData$y, k=8, q=1, first=0, last=0.99)
For the case of the optimal change points, we need to additionally define the maximum number \(k_{max}\) of change points. The final results will include scenarios where the number of change points varies from 1 to \(k_{max}\). And we need to specify the \(sic\) or \(bic\) penalty term to select the optimal number.
The outcomes of model fitting can be primarily categorized into three
components, the results of change point detection (\(A\)), the estimation of trends ( \(\hat{y}\) ), and the estimation of model
parameters ( \(\hat{\boldsymbol{\beta}}\) ). The results
obtained from the warm start will be denoted as \(A.all\), \(beta.all\), and \(y.all\) respectively. coef
function can extract and output the model’s results for a given number
\(k\) of change points. If the value of
\(k\) is not provided, it will output
\(A.all\), \(beta.all\) and \(y.all\).
#> $A
#> [1] 100 150 250 325 425 493
#>
#> $beta
#> [1] -1.0113527 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [7] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [13] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [19] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [25] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [31] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [37] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [43] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [49] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [55] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [61] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [67] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [73] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [79] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [85] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [91] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [97] 0.0000000 0.0000000 0.0000000 0.0000000 3.9825577 0.0000000
#> [103] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [109] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [115] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [121] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [127] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [133] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [139] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [145] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [151] -4.9625842 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [157] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [163] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [169] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [175] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [181] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [187] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [193] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [199] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [205] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [211] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [217] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [223] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [229] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [235] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [241] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [247] 0.0000000 0.0000000 0.0000000 0.0000000 1.9968180 0.0000000
#> [253] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [259] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [265] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [271] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [277] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [283] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [289] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [295] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [301] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [307] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [313] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [319] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [325] 0.0000000 3.9888904 0.0000000 0.0000000 0.0000000 0.0000000
#> [331] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [337] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [343] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [349] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [355] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [361] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [367] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [373] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [379] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [385] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [391] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [397] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [403] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [409] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [415] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [421] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 -6.9705018
#> [427] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [433] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [439] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [445] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [451] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [457] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [463] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [469] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [475] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [481] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [487] 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
#> [493] 0.0000000 -0.1299847 0.0000000 0.0000000 0.0000000 0.0000000
#> [499] 0.0000000 0.0000000
#>
#> $yhat
#> [1] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [6] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [11] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [16] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [21] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [26] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [31] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [36] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [41] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [46] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [51] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [56] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [61] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [66] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [71] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [76] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [81] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [86] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [91] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [96] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [101] 2.971204980 2.971204980 2.971204980 2.971204980 2.971204980
#> [106] 2.971204980 2.971204980 2.971204980 2.971204980 2.971204980
#> [111] 2.971204980 2.971204980 2.971204980 2.971204980 2.971204980
#> [116] 2.971204980 2.971204980 2.971204980 2.971204980 2.971204980
#> [121] 2.971204980 2.971204980 2.971204980 2.971204980 2.971204980
#> [126] 2.971204980 2.971204980 2.971204980 2.971204980 2.971204980
#> [131] 2.971204980 2.971204980 2.971204980 2.971204980 2.971204980
#> [136] 2.971204980 2.971204980 2.971204980 2.971204980 2.971204980
#> [141] 2.971204980 2.971204980 2.971204980 2.971204980 2.971204980
#> [146] 2.971204980 2.971204980 2.971204980 2.971204980 2.971204980
#> [151] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [156] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [161] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [166] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [171] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [176] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [181] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [186] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [191] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [196] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [201] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [206] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [211] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [216] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [221] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [226] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [231] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [236] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [241] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [246] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [251] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [256] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [261] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [266] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [271] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [276] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [281] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [286] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [291] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [296] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [301] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [306] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [311] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [316] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [321] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [326] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [331] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [336] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [341] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [346] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [351] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [356] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [361] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [366] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [371] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [376] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [381] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [386] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [391] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [396] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [401] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [406] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [411] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [416] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [421] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [426] -2.976172601 -2.976172601 -2.976172601 -2.976172601 -2.976172601
#> [431] -2.976172601 -2.976172601 -2.976172601 -2.976172601 -2.976172601
#> [436] -2.976172601 -2.976172601 -2.976172601 -2.976172601 -2.976172601
#> [441] -2.976172601 -2.976172601 -2.976172601 -2.976172601 -2.976172601
#> [446] -2.976172601 -2.976172601 -2.976172601 -2.976172601 -2.976172601
#> [451] -2.976172601 -2.976172601 -2.976172601 -2.976172601 -2.976172601
#> [456] -2.976172601 -2.976172601 -2.976172601 -2.976172601 -2.976172601
#> [461] -2.976172601 -2.976172601 -2.976172601 -2.976172601 -2.976172601
#> [466] -2.976172601 -2.976172601 -2.976172601 -2.976172601 -2.976172601
#> [471] -2.976172601 -2.976172601 -2.976172601 -2.976172601 -2.976172601
#> [476] -2.976172601 -2.976172601 -2.976172601 -2.976172601 -2.976172601
#> [481] -2.976172601 -2.976172601 -2.976172601 -2.976172601 -2.976172601
#> [486] -2.976172601 -2.976172601 -2.976172601 -2.976172601 -2.976172601
#> [491] -2.976172601 -2.976172601 -2.976172601 -3.106157313 -3.106157313
#> [496] -3.106157313 -3.106157313 -3.106157313 -3.106157313 -3.106157313
#> $A
#> [1] 100 150 250 325 425
#>
#> $beta
#> [1] -1.011353 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [8] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [15] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [22] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [29] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [36] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [43] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [50] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [57] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [64] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [71] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [78] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [85] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [92] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [99] 0.000000 0.000000 3.982558 0.000000 0.000000 0.000000 0.000000
#> [106] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [113] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [120] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [127] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [134] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [141] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [148] 0.000000 0.000000 0.000000 -4.962584 0.000000 0.000000 0.000000
#> [155] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [162] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [169] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [176] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [183] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [190] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [197] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [204] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [211] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [218] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [225] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [232] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [239] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [246] 0.000000 0.000000 0.000000 0.000000 0.000000 1.996818 0.000000
#> [253] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [260] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [267] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [274] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [281] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [288] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [295] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [302] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [309] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [316] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [323] 0.000000 0.000000 0.000000 3.988890 0.000000 0.000000 0.000000
#> [330] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [337] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [344] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [351] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [358] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [365] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [372] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [379] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [386] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [393] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [400] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [407] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [414] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [421] 0.000000 0.000000 0.000000 0.000000 0.000000 -6.982634 0.000000
#> [428] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [435] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [442] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [449] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [456] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [463] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [470] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [477] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [484] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [491] 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000
#> [498] 0.000000 0.000000 0.000000
#>
#> $yhat
#> [1] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [6] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [11] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [16] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [21] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [26] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [31] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [36] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [41] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [46] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [51] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [56] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [61] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [66] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [71] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [76] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [81] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [86] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [91] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [96] -1.011352715 -1.011352715 -1.011352715 -1.011352715 -1.011352715
#> [101] 2.971204980 2.971204980 2.971204980 2.971204980 2.971204980
#> [106] 2.971204980 2.971204980 2.971204980 2.971204980 2.971204980
#> [111] 2.971204980 2.971204980 2.971204980 2.971204980 2.971204980
#> [116] 2.971204980 2.971204980 2.971204980 2.971204980 2.971204980
#> [121] 2.971204980 2.971204980 2.971204980 2.971204980 2.971204980
#> [126] 2.971204980 2.971204980 2.971204980 2.971204980 2.971204980
#> [131] 2.971204980 2.971204980 2.971204980 2.971204980 2.971204980
#> [136] 2.971204980 2.971204980 2.971204980 2.971204980 2.971204980
#> [141] 2.971204980 2.971204980 2.971204980 2.971204980 2.971204980
#> [146] 2.971204980 2.971204980 2.971204980 2.971204980 2.971204980
#> [151] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [156] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [161] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [166] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [171] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [176] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [181] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [186] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [191] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [196] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [201] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [206] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [211] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [216] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [221] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [226] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [231] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [236] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [241] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [246] -1.991379256 -1.991379256 -1.991379256 -1.991379256 -1.991379256
#> [251] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [256] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [261] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [266] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [271] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [276] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [281] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [286] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [291] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [296] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [301] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [306] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [311] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [316] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [321] 0.005438767 0.005438767 0.005438767 0.005438767 0.005438767
#> [326] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [331] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [336] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [341] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [346] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [351] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [356] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [361] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [366] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [371] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [376] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [381] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [386] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [391] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [396] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [401] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [406] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [411] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [416] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [421] 3.994329185 3.994329185 3.994329185 3.994329185 3.994329185
#> [426] -2.988304507 -2.988304507 -2.988304507 -2.988304507 -2.988304507
#> [431] -2.988304507 -2.988304507 -2.988304507 -2.988304507 -2.988304507
#> [436] -2.988304507 -2.988304507 -2.988304507 -2.988304507 -2.988304507
#> [441] -2.988304507 -2.988304507 -2.988304507 -2.988304507 -2.988304507
#> [446] -2.988304507 -2.988304507 -2.988304507 -2.988304507 -2.988304507
#> [451] -2.988304507 -2.988304507 -2.988304507 -2.988304507 -2.988304507
#> [456] -2.988304507 -2.988304507 -2.988304507 -2.988304507 -2.988304507
#> [461] -2.988304507 -2.988304507 -2.988304507 -2.988304507 -2.988304507
#> [466] -2.988304507 -2.988304507 -2.988304507 -2.988304507 -2.988304507
#> [471] -2.988304507 -2.988304507 -2.988304507 -2.988304507 -2.988304507
#> [476] -2.988304507 -2.988304507 -2.988304507 -2.988304507 -2.988304507
#> [481] -2.988304507 -2.988304507 -2.988304507 -2.988304507 -2.988304507
#> [486] -2.988304507 -2.988304507 -2.988304507 -2.988304507 -2.988304507
#> [491] -2.988304507 -2.988304507 -2.988304507 -2.988304507 -2.988304507
#> [496] -2.988304507 -2.988304507 -2.988304507 -2.988304507 -2.988304507
In addition to the above results, print
function can
also output the ‘mse’, ‘bic’ and ‘sic’ results, demonstrating the
process of selecting the optimal change points.
#> [1] 4.791976044 1.749890740 1.409428027 0.352043640 0.010276593 0.010062126
#> [7] 0.009951844 0.009678123 0.009571390 0.009535951 0.009459965 0.009432617
#> [13] 0.009316856 0.009232618 0.009169940 0.009141021 0.008966564 0.008892028
#> [19] 0.008863110 0.008792653
#> [1] 808.3299 317.0643 221.3088 -459.8540 -2214.3679 -2212.4839
#> [7] -2205.5650 -2207.0807 -2200.1962 -2189.6218 -2181.1927 -2170.2110
#> [13] -2163.9560 -2156.0680 -2147.0448 -2136.1949 -2133.4005 -2125.1449
#> [19] -2114.3444 -2105.9059
#> [1] 828.8854 347.8976 262.4199 -408.4652 -2152.7014 -2140.5396
#> [7] -2123.3430 -2114.5809 -2097.4187 -2076.5665 -2057.8597 -2036.6003
#> [13] -2020.0674 -2001.9018 -1982.6008 -1961.4731 -1948.4010 -1929.8676
#> [19] -1908.7894 -1890.0731
There is another function called TFmetrics that is widely used in
simulation experiments. It outputs four trend filtering metrics between
the underlying trend and the data points, denoted as ‘MSE’, ‘MAD’, ‘dH’,
and ‘nknot’. For details, see the definition of the
TFmetrics
function.
metrics <- TFmetrics(BlocksData$y0,BlocksData$tau,FitBlocks.opt$yopt,FitBlocks.opt$Aopt/length(BlocksData$y0))
print(metrics)
#> MSE MAD dH nknot
#> 1 0.0001549419 0.0105785 0 5
The plot
function is used for visualizing the output
results, capable of drawing graphs that show how ‘mse’, ‘bic’ and ‘sic’
change with the number of change points, as well as the fitting trend
charts under a given number of change points.
Kim, S. J., Koh, K., Boyd, S., & Gorinevsky, D. (2009). \(\ell_1\) trend filtering. SIAM review, 51(2), 339-360.
Wen, C., Wang, X., & Zhang, A. (2023). \(\ell_0\) trend filtering. INFORMS Journal on Computing, 35(6), 1491-1510.