# gfpop Vignette ### Vincent Runge #### LaMME, Evry University ### March 2, 2022 > [Quick Start](#qs) > [Some examples](#se) > [Graph construction](#gc) > [Supplementary R functions](#suppl) ## Quick Start we present a basic use of the main functions of the `gfpop` package. More details about optional arguments are given in later sections. We install the package from Github: ```{r} #devtools::install_github("vrunge/gfpop") library(gfpop) ``` We simulate some univariate gaussian data (`n = 1000` points) with relative change-point positions `0.1, 0.3, 0.5, 0.8, 1` and means `1, 2, 1, 3, 1` with a variance equal to `1`. ```{r} n <- 1000 myData <- dataGenerator(n, c(0.1,0.3,0.5,0.8,1), c(1,2,1,3,1), sigma = 1) ``` We define the graph of constraints to use for the dynamic programming algorithm. A simple case is the up-down constraint with a penalty here equal to a classic `2 log(n)`. ```{r} myGraph <- graph(penalty = 2*log(n), type = "updown") ``` The gfpop function gives the result of the segmentation using `myData` and `myGraph` as parameters. We choose a gaussian cost. ```{r} gfpop(data = myData, mygraph = myGraph, type = "mean") ``` The vector `changepoints` gives the last index of each segment. It always ends with the length of the vector `vectData`. The vector `states` contains the states in which lies each mean. The length of this vector is the same as the length of `changepoint`. The vector `forced` is a boolean vector. A forced element means that two consecutive means have been forced to satisfy the constraint. For example, the "up" edge with parameter c is forced if m(i+1) - m(i) = c. The vector `parameters` contains the inferred means/parameters of the successive segments. The number `globalCost` is equal to the non-penalized cost, that is the value of the fit to the data ignoring the penalties for adding changes. ## Some examples ### Isotonic regression The isotonic regression infers a sequence of nondecreasing means. ```{r} n <- 1000 mydata <- dataGenerator(n, c(0.1, 0.2, 0.3, 0.4, 0.6, 0.8, 1), c(0, 0.5, 1, 1.5, 2, 2.5, 3), sigma = 1) myGraphIso <- graph(penalty = 2*log(n), type = "isotonic") gfpop(data = mydata, mygraph = myGraphIso, type = "mean") ``` In this example, we use in `gfpop` function a robust biweight gaussian cost with `K = 1` and the `min` parameter in order to infer means greater than `0.5`. ### Fixed number of change-points This algorithm is called segment neighborhood in the change-point litterature. In this example, we fixed the number of segments at 3 with an isotonic constraint. The graph contains two "up" edges with no cycling. ```{r} n <- 1000 mydata <- dataGenerator(n, c(0.1, 0.2, 0.3, 0.4, 0.6, 0.8, 1), c(0, 0.5, 1, 1.5, 2, 2.5, 3), sigma = 1) beta <- 0 myGraph <- graph( Edge(0, 1,"up", beta), Edge(1, 2, "up", beta), Edge(0, 0, "null"), Edge(1, 1, "null"), Edge(2, 2, "null"), StartEnd(start = 0, end = 2)) gfpop(data = mydata, mygraph = myGraph, type = "mean") ``` ### Robust up-down with constrained starting and ending states In presence of outliers we need a robust loss (biweight). We can also force the starting and ending state and a minimal gap between the means (here equal to `1`) ```{r} n <- 1000 chgtpt <- c(0.1, 0.3, 0.5, 0.8, 1) myData <- dataGenerator(n, chgtpt, c(0, 1, 0, 1, 0), sigma = 1) myData <- myData + 5 * rbinom(n, 1, 0.05) - 5 * rbinom(n, 1, 0.05) beta <- 2 * log(n) myGraph <- graph( Edge("Dw", "Up", type = "up", penalty = beta, gap = 1, K = 3), Edge("Up", "Dw", type = "down", penalty = beta, gap = 1, K = 3), Edge("Dw", "Dw", type = "null", K = 3), Edge("Up", "Up", type = "null", K = 3), StartEnd(start = "Dw", end = "Dw")) gfpop(data = myData, mygraph = myGraph, type = "mean") ``` If we skip all these constraints and use a standard fpop algorithm, the result is the following ```{r} myGraphStd <- graph(penalty = 2*log(n), type = "std") gfpop(data = myData, mygraph = myGraphStd, type = "mean") ``` ### abs edge With a unique `"abs"` edge, we impose a difference between the means of size at least 1. ```{r} n <- 10000 myData <- dataGenerator(n, c(0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1), c(0, 1, 0, 2, 1, 2, 0, 1, 0, 1), sigma = 0.5) beta <- 2*log(n) myGraph <- graph( Edge(0, 0,"abs", penalty = beta, gap = 1), Edge(0, 0,"null")) gfpop(data = myData, mygraph = myGraph, type = "mean") ``` Notice that some of the edges are forced, the vector `forced` contains non-zero values. ### Exponential decay The null edge corresponds to an exponential decay state if its parameter is not equal to 1. ```{r} n <- 1000 mydata <- dataGenerator(n, c(0.2, 0.5, 0.8, 1), c(5, 10, 15, 20), sigma = 1, gamma = 0.966) beta <- 2*log(n) myGraphDecay <- graph( Edge(0, 0, "up", penalty = beta), Edge(0, 0, "null", 0, decay = 0.966) ) g <- gfpop(data = mydata, mygraph = myGraphDecay, type = "mean") g ``` and we plot the result ```{r, expdecay} gamma <- 0.966 len <- diff(c(0, g$changepoints)) signal <- NULL for(i in length(len):1) {signal <- c(signal, g$parameters[i]*c(1, cumprod(rep(1/gamma,len[i]-1))))} signal <- rev(signal) ylimits <- c(min(mydata), max(mydata)) plot(mydata, type ='p', pch ='+', ylim = ylimits) par(new = TRUE) plot(signal, type ='l', col = 4, ylim = ylimits, lwd = 3) ``` ## Graph construction In the `gfpop` package, graphs are represented by a dataframe with 9 features and build with the R functions `Edge`, `Node`, `StartEnd` and `graph`. ```{r} emptyGraph <- graph() emptyGraph ``` `state1` is the starting node of an edge, `state2` its ending node. `type` is one of the available edge type (`"null"`, `"std"`, `"up"`, `"down"`, `"abs"`). `penalty` is a nonnegative parameter: the additional cost $\beta_i$ to consider when we move within the graph using a edge (or stay on the same node). `parameter` is annother nonnegative parameter, a characteristics of the edge, depending of its type (it is a decay if type is "null" and a gap otherwise). `K` and `a` are robust parameters. `min` and `max` are used to constrain the rang of value for the node parameter. We add edges into a graph as follows ```{r} myGraph <- graph( Edge("E1", "E1", "null"), Edge("E1", "E2", "down", 3.1415, gap = 1.5) ) myGraph ``` we can only add edges to this dataframe using the object `Edge`. The graph can contain information on the starting and/or ending edge to use with the `StartEnd` function. ```{r} beta <- 2 * log(1000) myGraph <- graph( Edge("Dw", "Dw", "null"), Edge("Up", "Up", "null"), Edge("Dw", "Up", "up", penalty = beta, gap = 1), Edge("Dw", "Dw", "down", penalty = beta), Edge("Up", "Dw", "down", penalty = beta), StartEnd(start = "Dw", end = "Dw")) myGraph ``` Some graphs are often used: they are defined by default in the `graph` function. To use these graphs, we specify a string `type` equal to `"std"`, `"isotonic"`, `"updown"` or `"relevant"`. For example, ```{r} myGraphIso <- graph(penalty = 12, type = "isotonic") myGraphIso ``` The function `Node` can be used to restrict the range of value for parameter associated to a node (called also a vertex). For example the following graph is an isotonic graph with inferred parameters between 0 et 1 only. ```{r} myGraph <- graph( Edge("Up", "Up", "up", penalty = 3.1415), Edge("Up", "Up"), Node("Up", min = 0, max = 1) ) myGraph ``` ## Supplementary R functions ### Data generator function the `dataGenerator` function is used to simulate `n` data-points from a distribution of `type` equal to `"mean"`, `"poisson"`, `"exp"`, `"variance"` or `"negbin"`. Standard deviation parameter `sigma` and decay `gamma` are specific to the Gaussian mean model. `size` is linked to the R `rnbinom` function from R stats package. ### Standard deviation estimation We often need to estimate the standard deviation from the observed data to normalize the data or choose the edge penalties. The `sdDiff` returns such an estimation with the default HALL method [Hall et al., 1990] well suited for time series with change-points. [Back to Top](#top)