User adunaic - Cross Validated most recent 30 from stats.stackexchange.com 2019-08-25T20:43:50Z https://stats.stackexchange.com/feeds/user/13409 http://www.creativecommons.org/licenses/by-sa/3.0/rdf https://stats.stackexchange.com/q/387367 2 How small is too small to fit a reasonable long memory model? adunaic https://stats.stackexchange.com/users/13409 2019-01-15T17:26:21Z 2019-01-15T20:15:55Z <p>When looking at papers about long memory they tend to analyze data sets whose length is in the thousands, see <a href="http://www.math.canterbury.ac.nz/~m.reale/pub/Reaetal2011.pdf" rel="nofollow noreferrer">http://www.math.canterbury.ac.nz/~m.reale/pub/Reaetal2011.pdf</a> for an example.</p> <p>My question is to the long memory researchers and practitioners out there. What rule of thumb do you use to decide whether a data set is too small to be able to appropriately detect/estimate long memory?</p> <p>[Naturally the smaller the long memory parameter the more observations you will require to detect it, but i'm after a general rule of thumb rather than exact notions around the number of observations required to detect a specific effect size.]</p> https://stats.stackexchange.com/questions/362126/-/362690#362690 0 Answer by adunaic for Changepoint/Step Detection in Univariate Time Series adunaic https://stats.stackexchange.com/users/13409 2018-08-17T14:26:10Z 2018-08-17T14:26:10Z <p>Some test data that has some similar properties, code is in R:</p> <pre><code>set.seed(1) a=rep(c(1,5,9,15),each=250) x=1:1000 y=a+-0.02*x+rnorm(1000,sd=0.4) </code></pre> <p><a href="https://i.stack.imgur.com/mO1ea.jpg" rel="nofollow noreferrer"><img src="https://i.stack.imgur.com/mO1ea.jpg" alt="enter image description here"></a></p> <p>To analyze this: library(EnvCpt) out=envcpt(y,models="trendcpt") cpts(out\$trendcpt) # gives changes at 250, 500, 750 as simulated.</p> <p>plot(out\$trendcpt)<a href="https://i.stack.imgur.com/7AXo4.jpg" rel="nofollow noreferrer"><img src="https://i.stack.imgur.com/7AXo4.jpg" alt="enter image description here"></a></p> <p>The <code>envcpt</code> function can fit several models and compare the fits with and without changepoints so this is why we specify <code>models="trendcpt"</code> so it only fits the single model.</p> <p>This can be run from Python using <code>rpy2</code> or alternative packages that can call R from Python. Unfortunately we don't have a Python implementation yet.</p> https://stats.stackexchange.com/questions/357687/-/357877#357877 1 Answer by adunaic for Finding the change point before a significant increase adunaic https://stats.stackexchange.com/users/13409 2018-07-18T22:43:05Z 2018-07-18T22:43:05Z <p>You could use the <code>EnvCpt</code> package in <code>R</code> which fits mean, trend, AR and changepoint models. It then gives you the best fit of all models and you can use AIC (or another metric) to choose the best model.</p> <pre><code>library('EnvCpt') out=envcpt(df) out\$summary # gives the fit and number of parameters for each model plot(out,type='aic') # plots the aic values, trendcpt+AR2 model is clearly the best out\$trendar2cpt # gives the fit for trendcpt+AR2 </code></pre> <p>The above gives the best model as a Trend+AR2+cpt model (assuming Normal errors) with a changepoint after 6 observations. I'm not sure you would want to fit a Trend+AR2 model to 6 observations though - but that is your call.</p> <p>The next lowest AIC value is the Trend+AR2 model. Thus indicating that there is not a clear changepoint in this data. This could be that there is not enough data to be confident there is a changepoint at observation 6, or the change is too small (relative to the noise), or a combination.</p> https://stats.stackexchange.com/questions/338021/-/339856#339856 1 Answer by adunaic for Finding changepoints in the movement of a car with respect to a lead car (using ecp package in R) adunaic https://stats.stackexchange.com/users/13409 2018-04-11T12:00:14Z 2018-04-11T12:00:14Z <p>Due to the length of the data and the gradual change you are expecting (and have seen with the final plot) i'm not sure that <code>ecp</code> is the best approach to use. Instead you may want to fit a changepoint model which allows for slopes and large datasets. I suggest you give the <code>EnvCpt</code> package a try. The primary focus of the package is for model selection and so the main function fits 8 different models to allow the user to choose the most appropriate.</p> <p>Having said this you can use the underlying functions to fit a specific model form, such as the trend you see in your plots. Use the <code>EnvCpt:::cpt.reg</code> function to do this. In contrast to <code>ecp</code> this means making model assumptions but it looks like it might give you the segmentation you want.</p> https://stats.stackexchange.com/questions/312851/-/313036#313036 1 Answer by adunaic for Bayesian Information Criterion - Non-physical model selection adunaic https://stats.stackexchange.com/users/13409 2017-11-10T10:59:21Z 2017-11-10T10:59:21Z <p>What you described sounds more like constrained changepoint detection. This is where you impose some constraints on the resulting segmentation profile. Specifically there have been up-down constraints (where a positive jump in mean must be followed by a negative jump):</p> <p><a href="https://arxiv.org/abs/1703.03352" rel="nofollow noreferrer">https://arxiv.org/abs/1703.03352</a></p> <p><a href="https://cran.r-project.org/web/packages/PeakSegOptimal/index.html" rel="nofollow noreferrer">https://cran.r-project.org/web/packages/PeakSegOptimal/index.html</a></p> <p>I'm also aware of others working on different constraints but these have not been published yet.</p> <p>What you want to do is to restrict the fit rather than impose the rules at the end. This is because there may be a slightly sub-optimal fit for one \$k\$ that satisfies the rules but is a better fit than an alternative \$k'\$ where the optimal does satisfy the rules. Thus imposing the rules at the fitting stage results in all \$K\$ segmentations satisfying the rules and then you can use BIC on those as usual. The crux of the attached paper (and works in progress) is that imposing these rules makes obtaining the best segmentation that satisfies the rules difficult.</p> <p>In summary, it is fine to constrain your segmentations to follow specific rules but this should be imposed in the fit and not in the penalization step.</p> https://stats.stackexchange.com/questions/29910/-/300712#300712 1 Answer by adunaic for multiple change point analysis simultaneously for mean AND variance WITHOUT distribution assumption in R adunaic https://stats.stackexchange.com/users/13409 2017-08-31T09:46:37Z 2017-08-31T09:46:37Z <p>For those looking at this thread there is now also the <code>changepoint.np</code> package for R (available on CRAN) which contains a method for a change in distribution and is currently being revamped with more specific robust changepoint methods.</p> https://stats.stackexchange.com/questions/300377/-/300579#300579 0 Answer by adunaic for Identify two major change points in time series data and summarise for several replicates in R adunaic https://stats.stackexchange.com/users/13409 2017-08-30T16:13:22Z 2017-08-30T16:13:22Z <p>My first suggestion would be to use a package for changepoint analysis rather than outlier detection. An example would be the <code>changepoint</code> package in <code>R</code>. Here is the code for the analysis of the above data:</p> <pre><code>library(changepoint) out=cpt.meanvar(dat.ts) plot(out) cpts(out) </code></pre> <p>The mean and variance of your data is changing so it makes sense to use the <code>cpt.meanvar</code> function.</p> <p>In terms of a multivariate comparison, I would advise using a multivariate changepoint detection approach. If you want the changes in the same place in all series then you can use the <code>InspectChangepoint</code> <code>R</code> package which appears to work well.</p> https://stats.stackexchange.com/questions/299750/-/300307#300307 0 Answer by adunaic for How do I do a change point analysis on a sparse data set in python? adunaic https://stats.stackexchange.com/users/13409 2017-08-29T08:51:18Z 2017-08-29T08:51:18Z <p>You can use standard changepoint analysis techniques on data with missing points but you just have to be careful about what you infer when the changepoints are between missing values.</p> <p>As an example consider p02. If an algorithm returned a changepoint at day 3 then this doesn't mean that the change is at day 3, instead it means that the change is somewhere after the data point on day 3 was collected and before the data point on day 5.</p> <p>You will want to be careful as you only have integer values so you shouldn't really use a Normal distribution for your data - especially as you have so few data points. You might be best using a multinomial distribution as you have 5 outcomes.</p> https://stats.stackexchange.com/questions/271084/-/271283#271283 1 Answer by adunaic for How do I split a vector and minimize absolute error efficiently? adunaic https://stats.stackexchange.com/users/13409 2017-04-01T20:44:37Z 2017-04-01T20:44:37Z <p>The changepoint.np R package does this in O(n) time using a range of quantiles. </p> <pre><code>library(changepoint.np) set.seed(1) n&lt;-20L y&lt;-rnorm(n) out&lt;-cpt.np(y) cpts(out) </code></pre> <p>This uses the defaults and identifies that there are no changes in the data. If you want to know which location is the most likely for a changepoint then you can use the CROPS penalty:</p> <pre><code>out1&lt;-cpt.np(y,penalty="CROPS",pen.value=c(3,7)) cpts.full(out1) </code></pre> <p>The CROPS penalty gives you all segmentations within a range of penalties, using <code>cpts.full</code> we can see all the segmentations with a penalty between 3 and 7 in the above. It gives from no changepoints at a penalty of 7 to eight changes for a penalty of 3 and all in between.</p> <p>You can set nquantiles=1 but i'm not sure what quantile that would calculate, you would need to read the original paper here (open access): <a href="https://link.springer.com/article/10.1007/s11222-016-9687-5" rel="nofollow noreferrer">PAPER</a></p> <p>It isn't quite the answer to the question as it isn't a single changepoint, but it is a more general solution, potentially multiple changes and different quantiles, with a faster O(n) computational time than those listed already.</p> https://stats.stackexchange.com/questions/262043/-/262271#262271 0 Answer by adunaic for Forecast time series with breakpoint using Holt-Winters (R) adunaic https://stats.stackexchange.com/users/13409 2017-02-16T11:53:10Z 2017-02-16T11:53:10Z <p>One of the typical assumptions when doing changepoint detection is that the segments are independent. Thus if this is the case, the data prior to the break is uninformative for the data after the break. Thus you are safe to only use the data post-break for training a forecasting model.</p> https://stats.stackexchange.com/questions/238134/-/238360#238360 2 Answer by adunaic for What is change point analysis for At Most One Change(AMOC)? adunaic https://stats.stackexchange.com/users/13409 2016-10-04T12:02:41Z 2016-10-04T12:02:41Z <p>It appears as though you are using the <code>changepoint</code> package in <code>R</code> from the function names mentioned and the <code>AMOC</code> definition. If this is correct then the help files for the package give references:</p> <pre><code> ?cpt.mean </code></pre> <p>gives the following references:</p> <pre><code> Change in Normal mean: Hinkley, D. V. (1970) Inference About the Change-Point in a Sequence of Random Variables, Biometrika 57, 1–17 CUSUM Test: M. Csorgo, L. Horvath (1997) Limit Theorems in Change-Point Analysis, Wiley </code></pre> <p>I acknowledge that it isn't 100% clear that <code>AMOC</code> should use these references. If you have <code>test.stat="Normal"</code> which is the default then the first Hinkley reference is the one you want. If you use <code>test.stat="CUSUM"</code> then you want the second reference.</p> <p>The majority of changepoint techniques start with <code>AMOC</code> or a single changepoint. The <code>changepoint</code> package was designed more for multiple changes but allows single changes using <code>AMOC</code>.</p> <p>The above covers the reference part of your question. For the "what pen.value option should I use" that is like asking "who is the best guitarist?", the answer depends on personal experience.</p> <p>For <code>AMOC</code> you have the option of using asymptotic penalty values, for example, if you want to be 95% confident that a changepoint has occurred then you would use <code>penalty="Asymptotic", pen.value=0.05</code>. If you are happy with any of the other penalty choices i.e. MBIC (default), BIC, SIC, AIC, Hannan-Quinn then you don't need to specify <code>pen.value</code> as it is taken care of by setting the penalty to be one of the options, e.g. <code>cpt.mean(data,penalty='SIC')</code> would be valid if you wanted to use the SIC penalty.</p> <p>If you don't like any of the options then you can set your own penalty but you have to select a value that works for your problem - sadly it is still an open research question as to which penalty is best.</p> <p>Finally, you will get different answers when using the different <code>cpt.mean, cpt.var, cpt.meanvar</code> functions as they are doing different things. <code>cpt.mean</code> looks for a change in the mean value only, assuming a constant variance. <code>cpt.var</code> looks for a change in the variance only, assuming a constant mean. <code>cpt.meanvar</code> looks for a change in both the mean and variance.</p> <p>If you haven't read the paper associated to the package then I suggest you read it as it gives several examples demonstrating how to change the penalties and also how to use the 3 different functions.</p> <p><a href="https://www.jstatsoft.org/article/view/v058i03" rel="nofollow">URL to paper</a></p> https://stats.stackexchange.com/questions/231078/-/231081#231081 0 Answer by adunaic for Multivariate number generation with pairwise correlation adunaic https://stats.stackexchange.com/users/13409 2016-08-22T12:49:57Z 2016-08-22T13:01:22Z <p>The covariance is the same as correlation when the variances are 1. Thus we can generate a sigma matrix as follows with the pairwise correlations equal to 1:</p> <pre><code>sigma=toeplitz(c(1,rep(0.85,29))) </code></pre> <p>Note i'm using 29 as you asked for <code>p=30</code>. Then we can simulate any length of data from this sigma matrix using <code>rmvnorm</code>:</p> <pre><code>set.seed(1) # to be reproducible data=rmvnorm(n=200,sigma=sigma) </code></pre> <p>The default mean is 0 for each dimension so I haven't included it in the above. You can then check the correlation of individual pairs as follows:</p> <pre><code>cor(data[,1],data[,5]); </code></pre> <p>Or you can get the correlation matrix:</p> <pre><code>cor(data); </code></pre> <p>You can see that most of them are around 0.85 ish. You can make it more readable doing:</p> <pre><code>round(cor(data),2); </code></pre> https://stats.stackexchange.com/questions/230772/-/231077#231077 0 Answer by adunaic for Switchpoint calculation in linear models with probabilistic programming adunaic https://stats.stackexchange.com/users/13409 2016-08-22T12:24:29Z 2016-08-22T12:24:29Z <p>The <code>strucchange</code> package in <code>R</code> is there to do just this. Here is a minimal working example presuming you have already done <code>install.packages('strucchange')</code>:</p> <pre><code>library(`strucchange`) set.seed(1) x=1:250 y=c(0.01*x[1:100],1.5-0.02*(x[101:250]-101)) ynoise=y+rnorm(250,0,0.2) ans=breakpoints(ynoise~x) </code></pre> <p>The <code>ans</code> object contains all the information on the fit. For example:</p> <pre><code>ans\$breakpoints gives the changepoint locations (100 in this example) using the BIC information criterion for the decision on number of changepoints. If you want to use another criterion you can get the fitted values under `ans\$RSS.table`. You can get the fit using the standard `lm` function using lm(ynoise~x+x*breakfactor(ans)) </code></pre> <p>I've put <code>x*breakfactor(ans)</code> here as we have a break in both the intercept and trend, just putting <code>breakfactor(ans)</code> here would only put a break in the intercept.</p> https://stats.stackexchange.com/questions/58598/-/226060#226060 0 Answer by adunaic for Test to identify change in median in a time-series adunaic https://stats.stackexchange.com/users/13409 2016-07-28T08:28:35Z 2016-07-28T08:28:35Z <p>There is a new R package, <code>changepoint.np</code> which detects changes in the empirical distribution function. A minimal working example (demonstrated using a Normal distribution but it is a nonparametric test) is:</p> <pre><code>set.seed(1) x=c(rnorm(100),rnorm(100,2)) out=cpt.np(x) # runs the ecdf changepoint method plot(out) # plots the data with changepoints marked cpts(out) # lists the changepoints identified </code></pre> <p>Alternatively, if you want a change in mean and are willing to make distributional assumptions the <code>changepoint</code> package contains the <code>cpt.mean</code> function or <code>cpt.meanvar</code> if the variance is changing too.</p> <p>If you have dependence in the data then you just need to inflate the default penalty otherwise you get spurious changes that are just due to the dependence structure. For that you might want to use the <code>CROPS</code> (changepoints for a range of penalties) options available in both packages.</p> https://stats.stackexchange.com/questions/224644/-/224746#224746 1 Answer by adunaic for Distribution Fitting with multiple changepoints adunaic https://stats.stackexchange.com/users/13409 2016-07-20T14:36:04Z 2016-07-20T14:36:04Z <p>The <code>changepoint</code> package in the <code>R</code> software provides the option for doing the above. To demonstrate with a toy example:</p> <pre><code>library(changepoint) set.seed(1) x=c(rnorm(50,0,1),rnorm(50,5,3),rnorm(50,10,1),rnorm(50,3,10)) out=cpt.meanvar(x,method="PELT") cpts(out) # gives 50,100,150 as changepoints </code></pre> <p>The above uses the <code>cpt.meanvar</code> function as you want to identify a change in mean and variance. The default is to assume a Normal distribution so you are good there. The default is to identify 1 change so we use <code>method='PELT'</code> to tell the function to identify multiple changes uses the <code>PELT</code> algorithm (exact minimization of the objective function in approximately linear time). </p> <p>There are many things you can do with the output, <code>cpts(out)</code> gives you the changepoint locations but you can also plot:</p> <pre><code>plot(out,cpt.width=3) </code></pre> <p>Here i've changed the width of the changepoint lines (means) to be thicker than the default (1), you can also use all your normal plotting arguments.</p> <p>We have used the default penalty <code>MBIC</code> in calculating the segmentation above but if your data are not independent within each segment you need to change this. There is the <code>CROPS</code> penalty which gives you a range of segmentations between two penalty values:</p> <pre><code>out1=cpt.meanvar(x,method='PELT',penalty='CROPS',pen.value=c(10,500)) plot(out1,diagnostic=TRUE) </code></pre> <p>You need to make sure that your <code>pen.value</code> upper limit gives a segmentation with 0 changepoints for the diagnostic to make sense. In the diagnostic plot you get the number of changepoints against the negative-log-likelihood. The idea is that when true changes are added the negative-log-likelihood goes down a large amount, when you add false changes then the negative-log-likelihood show a small improvement in fit. For this simulation it is clear that the appropriate number of changes is 3 but for real data it can be less clear. Once the number of changes has been decided you can plot that number using:</p> <pre><code>plot(out1,ncpts=3,cpt.width=3) </code></pre> https://stats.stackexchange.com/questions/223779/-/223938#223938 0 Answer by adunaic for How to determine correct changepoints from Posterior Probabilities (bcp R package)? adunaic https://stats.stackexchange.com/users/13409 2016-07-15T14:16:26Z 2016-07-15T14:16:26Z <p>This is essentially the difference between the output from a frequentist and Bayesian procedure. To get frequentist type output i.e. point estimates one typically uses the Bayesian Maximum a posteriori (MAP) estimation. This essentially chooses a threshold and returns the modes above that threshold. You are then back to the frequentist world of "what threshold is appropriate?" which is the same as asking "how long is a piece of string?".</p> <p>Typically people use thresholds based on the Bayes Factor, for example:</p> <p><a href="http://www.lancs.ac.uk/~killick/Pub/EckleyFearnheadKillick2010.pdf" rel="nofollow">http://www.lancs.ac.uk/~killick/Pub/EckleyFearnheadKillick2010.pdf</a></p> <p>shows how this choice can be made and:</p> <p><a href="http://eprints.lancs.ac.uk/8189/1/PScpt2.pdf" rel="nofollow">http://eprints.lancs.ac.uk/8189/1/PScpt2.pdf</a></p> <p>provides a histogram of number of changes against posterior probability too indicating that the mode of that histogram should be chosen.</p> https://stats.stackexchange.com/questions/208710/-/208803#208803 2 Answer by adunaic for Change point identification adunaic https://stats.stackexchange.com/users/13409 2016-04-22T14:27:54Z 2016-04-22T14:27:54Z <p>This is an interesting question. In essence the answer will come down to what assumptions you are making and how large the change you record is. If the change is large enough then you will see it immediately. For example, using the <code>changepoint</code> package in the statistical software 'R':</p> <p><code>library(changepoint) set.seed(10) cpts(cpt.mean(c(rnorm(50),10),method='PELT'))</code></p> <p>Here we have simulated 50 data points from a Normal distribution with mean 0 and variance 1 then added a 10 at the end as our "new" data point, maybe an example of removing an object. The <code>cpt.mean</code> function is for detecting changes in mean and defaults to a Normal distribution assumption. I have then put <code>cpts()</code> around it so only the changes are returned.</p> <p>The above correctly identifies the change as being at time 50, i.e. the last data point is from a new regime. This is identified as the new data point is very different from the existing data.</p> <p>In essence (when using distributional assumptions) if the new data point is within the expected range of the distribution (in this case you expect to see a value of 3 or above roughly once every thousand observations and a value of 10 or above is within machine precision of 0) then a change won't be signalled immediately. It may be that your change is too small to be able to detect with a single data point (or that your assumptions mean this is so).</p> https://stats.stackexchange.com/questions/164661/-/175855#175855 2 Answer by adunaic for Interpretation of Multiple Change point results and graph for offline analysis in R adunaic https://stats.stackexchange.com/users/13409 2015-10-07T11:46:41Z 2015-10-07T11:46:41Z <p>What the output is telling you is the following:</p> <ul> <li>There are no changepoints detected with the penalty used (hence blank on "changepoints" field)</li> <li>The penalty value you entered (0.02) was not used as you specified BIC which gives a penalty value of 2.080237.</li> <li>The plot reiterates this with a single horizontal line for the mean with no breaks in it.</li> </ul> <p>If you want to use a penalty of 0.02 then you should put <code>penalty="Manual"</code> in the function call.</p> <p>However, you don't need to play with the penalty as the "Range of Segmentations" in the output tells you that if there is 1 changepoint it is at 1937, for 2 they would be at 10 and 5 etc. You can change Q to get more segmentations listed.</p> https://stats.stackexchange.com/questions/97946/-/99504#99504 3 Answer by adunaic for Changepoints in R adunaic https://stats.stackexchange.com/users/13409 2014-05-21T08:49:18Z 2015-05-21T07:49:00Z <p>IrishStat is correct in that you are trying to identify a change in mean, not a change in variance. Thus in the <code>changepoint</code> package you should be using <code>mean=cpt.mean(results\$P1, method="PELT")</code> instead. As for the <code>bcp</code> package this gives no changes in mean.</p> <p>The <code>cpt.var</code> function gave 3 changes in variance because the variances of each part, calculated using</p> <pre><code>segvar=param.est(var)\$variance segvar &gt;0.02126190 0.09944762 0.00080000 0.07043333 segvar[-1]/segvar[-length(segvar)] &gt;4.677267637 0.008044436 88.041666667 </code></pre> <p>Typically changes in variance are detected with roughly 80% power or more if the ratio of neighbouring variances is greater than 3 (or less than 1/3). The ratio of these variances clearly fits this paradigm which is why the changes were detected but not necessarily in the places you would have expected to identify a change in mean.</p> <p>Note that this is all based on a penalty that only penalizes the number of changepoints. This is why segment lengths of 2/3 observations are detected. If the application suggests segments of small lengths such as these are implausible then I would use a penalty that penalizes segment length too (or set a minimum segment length).</p> <p>See introductory references at <a href="http://www.changepoint.info" rel="nofollow">www.changepoint.info</a> for more background details on changepoint analysis. There is also a list of various changepoint open source software packages there.</p> https://stats.stackexchange.com/questions/139481/-/139569#139569 2 Answer by adunaic for Detecting if samples belong to a given distribution adunaic https://stats.stackexchange.com/users/13409 2015-02-27T09:40:05Z 2015-02-27T09:40:05Z <p>I believe that Anastasia wants to calculate the departure from the Exponential distribution as new data arrives. Thus in essence, she can perform a changepoint test to decide if the recent data obeys the same distribution as previously seen data.</p> <p>There are several options here, if you are sure about the Exponential assumption then you could use a likelihood ratio test statistic. The downside with this is that you are also assuming an Exponential distribution for the new data but just with a different rate parameter. This may work in practice and give you the answers you want but you would have to try it.</p> <p>What might be more beneficial is to use a nonparametric test statistic such as a Kolmogorov-Smirnov test as described by Glen_b. However, I would adovate using this in a changepoint setting so that as a new data point arrives, you see if a change has occurred recently.</p> <p>Both these methods are available in the cpm R package. Unfortunately due to not passing some CRAN checks, this has been archived <a href="http://cran.r-project.org/src/contrib/Archive/cpm/" rel="nofollow">here</a>. I have used the pacakge for a while and can vouch that the methods are correct, the package has just failed a CRAN check on the Solaris operating system which the package author is finding awkward to correct. You can install it manually by downloading the source from the archive. I don't believe that there are any other changepoint packages on CRAN that are capable of detecting changepoints from a data stream.</p> https://stats.stackexchange.com/questions/121019/-/123550#123550 1 Answer by adunaic for Correlation on ordered subset adunaic https://stats.stackexchange.com/users/13409 2014-11-11T13:46:36Z 2014-11-11T13:46:36Z <p>I haven't come across anyone using correlation as a measure to detect where the changepoint occurs. I think that for the data you describe, the easiest thing to do is to perform a change in regression. This will easily identify both the changepoint and the slope (or not) of the pre change and post change mean.</p> <p>In <code>R</code> the <code>strucchange</code> package can easily do this: <a href="http://cran.r-project.org/web/packages/strucchange/index.html" rel="nofollow">http://cran.r-project.org/web/packages/strucchange/index.html</a></p> <p>The documentation, in my opinion, isn't easy to follow so here is a short example:</p> <pre><code>set.seed(987234) y=c(1:50,rep(0,50))+rnorm(100,sd=c(rep(10,50),rep(20,50))) plot(y,pch=19) # looks similar to your data library(strucchange) breakpoints(y~c(1:100)) </code></pre> <p>The last line gives 1 breakpoint with the first observation of the new segment being 51.</p> https://stats.stackexchange.com/questions/120933/-/123544#123544 1 Answer by adunaic for Maximizing Log-Likelihood Estimation for Changepoint Detection adunaic https://stats.stackexchange.com/users/13409 2014-11-11T13:23:00Z 2014-11-11T13:23:00Z <p>It is not clear from the presentation what distributional assumptions are being made in order to calculate the likelihoods.</p> <p>It might be simpler for you to look at the recently published BreakoutDetection package published by the same authors: <a href="https://blog.twitter.com/2014/breakout-detection-in-the-wild" rel="nofollow">https://blog.twitter.com/2014/breakout-detection-in-the-wild</a>.</p> <p>But if you are more interested in learning about changepoint detection and how to use likelihoods then read on.</p> <p>Normally you have a time series <code>y</code> which we assumes has <code>n</code> observations. In order to use likelihoods we need to make some assumptions about the distribution that <code>y</code> comes from. It is often assumed that the data come from Normal distribution (although this isn't always appropriate). Following a distributional assumption you need to decide which parameters of the distribution are allowed to change, e.g. mean, variance, both. </p> <p>If all the parameters change then you proceed by splitting your data into 2 halves, before change and after change, and use maximum likelihoods to fit the parameters to each half. In this way it is like a normal analysis where you just have data points and you fit a model to that data. If not all parameters can change then you need to estimate those that don't change from the whole data (not always an easy task, especially when parameters are linked).</p> <p>The trick with changepoint analysis is that you don't know where the change is, so you have to calculate the likelihood for each possible changepoint location and take the most likely as the hypothesized changepoint location. It is this location and likelihood that you then test to see if the change is significant by comparing the likelihood ratio to a threshold to see if the change is significant (comparing to c in slide 16).</p> https://stats.stackexchange.com/questions/116363/-/118558#118558 1 Answer by adunaic for How to characterize abrupt change? adunaic https://stats.stackexchange.com/users/13409 2014-10-10T08:43:31Z 2014-10-10T08:43:31Z <p>The area of statistics that you are looking for is changepoint analysis. There is a website <a href="http://www.changepoint.info" rel="nofollow">here</a> that will give you an overview of the area and also have a page for software.</p> <p>If you are an <code>R</code> user then i'd recommend the <code>changepoint</code> package for changes in mean and the <code>strucchange</code> package for changes in regression. If you want to be Bayesian then the <code>bcp</code> package is good too.</p> <p>In general you have to choose a threshold which indicates the strength of the changes you are looking for. There are, of course, threshold choices that people advocate in certain situations and you can use asymptotic confidence levels or bootstrapping to get confidence too.</p> https://stats.stackexchange.com/questions/115733/-/116461#116461 0 Answer by adunaic for Detecting a step change in time ordered data adunaic https://stats.stackexchange.com/users/13409 2014-09-23T14:08:00Z 2014-09-23T14:08:00Z <p>The changepoint methods such as the <code>strucchange</code> package assume that the series occur in an ordered manner. Thus even though you have several observations at a single time point, the order in which they occur in your data table matters. Having said this the <code>strucchange</code> package works well on the series above identifying the correct (in some sense) observations. Strucchange doesn't work on data tables so converting to a matrix we get</p> <pre><code>set.seed(10) dat &lt;- matrix(c(rep(0,30), rep(5,30), rep(10,30), rep(15,30), rep(20,30), rnorm(30, 5, .5), rnorm(30, 4, .5), rnorm(30, 3, .5), rnorm(30, 2, .5), rnorm(30, 1, .5)),ncol=2) library(strucchange) breakpoints(dat[,2]~dat[,1]) </code></pre> <p>No changepoint is identified. Now for dat2:</p> <pre><code>dat2 &lt;- dat change &lt;- which(dat2[,1] %in% 10) dat2[change,2] &lt;- dat2[change,2] + 1 change &lt;- which(dat2[,1] %in% 15) dat2[change,2] &lt;- dat2[change,2] - 1 breakpoints(dat2[,2]~dat2[,1]) </code></pre> <p>Optimal is 2 changepoints at 30 and 90 (0-30 is the slightly higher mean of 5 in t=0, 31-90 have mean 4 and 91-150 have mean 1).</p> <p>Now for dat3:</p> <pre><code>dat3 &lt;- dat change &lt;- which(dat3[,1] %in% 0) dat3[change,2] &lt;- dat3[change,2] - 3 breakpoints(dat3[,2]~dat3[,1]) </code></pre> <p>Here we have 1 changepoint at 59 which suggests that the first two form an upward trend followed by the remaining which form a downward trend.</p> <p>This behaviour of the incorrect changepoint is due to the fact that only have 1 time observation before the change. Modifying this to 2 we estimate correctly:</p> <pre><code>dat4=rbind(matrix(c(rep(-5,30),rnorm(30,2,0.5)),ncol=2),dat3) breakpoints(dat4[,2]~dat4[,1]) </code></pre> <p>1 changepoint at 60.</p> <p>You are trying to fit a change in the regression relationship. In some sense, repeated observations are not an issue, this gives you are handle on the variability. What is most important is that you have enough information to estimate the relationship prior to the change and also after the change as dat3 and dat4 demonstrate (albeit simply).</p> https://stats.stackexchange.com/questions/103273/-/103523#103523 0 Answer by adunaic for find the point at which the curve significantly shoots up adunaic https://stats.stackexchange.com/users/13409 2014-06-16T09:53:57Z 2014-06-16T09:53:57Z <p>This really depends on what the data looks like.</p> <p>Without a plot and from the description it sounds like the mean increases during the rainy season. If it is just a case of a baseline value of rainfall outside the rainy season and then this switches to another (higher) baseline during rainfall season then you are looking at a change in mean model. You can fit this using the <code>cpt.mean</code> function in the R package <code>changepoint</code>.</p> <p>Alternatively, if there aren't really two baselines and during the rainy season you can see an increase in the mean but it isn't really constant and the variability is higher then you might want to transform your data. The easiest way to do this is to take first differences, i.e. \$x_2-x_1\$ (you can use the <code>diff</code> function in R). Then you find the changepoint in the differences you would use the <code>cpt.var</code> function in the <code>changepoint</code> package.</p> <p>Both of these find changes in the mean and variance respectively. Without seeing the data it is hard to know which one (if either) might be appropriate for your data.</p> https://stats.stackexchange.com/questions/96331/-/96940#96940 1 Answer by adunaic for required: good and straightforward method to detect change points in dependent univariate time series using r adunaic https://stats.stackexchange.com/users/13409 2014-05-08T15:41:29Z 2014-05-08T15:41:29Z <p>As there has been no reply i'll throw a suggestion into the pot.</p> <p>As far as i'm aware there is no non-Bayesian implementation of a dependent changepoint test currently in R.</p> <p>Having said this the <code>strucchange</code> package implements changes in regression and as such you can construct AR models manually in this framework. Take an AR(1) model as an example:</p> <pre><code>set.seed(1) # reproducibility x=c(arima.sim(model=list(ar=0.8),n=50),arima.sim(model=list(ar=-0.8),n=50)) # simulate some data with a change library(strucchange) # load the package breakpoints(x[2:100]~x[1:99]) # detect breaks in an AR1 model </code></pre> <p>The final line fits an AR(1) model to the data and assesses if there is a change in the AR(1) parameter. You have to essentially ignore the first p data points in order to fit an AR(p) model.</p> <p>This is probably not exactly what we want as you are required to fully specify the AR model but i'd argue it is better than nothing.</p> <p>I have some code to fit general AR(p) models to changepoints using likelihoods where the code automatically chooses the best <code>p</code> for a segment but I haven't had time to incorporate this into the <code>changepoint</code> package yet. I can post this to a repository if it would be useful to you.</p> <p>For information there is a list of R packages related to changepoint analysis at: <a href="http://www.changepoint.info/software" rel="nofollow">www.changepoint.info/software</a></p> https://stats.stackexchange.com/questions/96012/-/96073#96073 0 Answer by adunaic for Determine if a % Change occurred adunaic https://stats.stackexchange.com/users/13409 2014-05-02T09:42:26Z 2014-05-02T09:42:26Z <p>This is a standard t-test where you have group1 (pre event) and group2 (post event). You can calculate the means of each group and do a test to see if they are the same or not. If they are statistically different then you can infer what %age change you have had.</p> <p>The changepoint tag is not warranted on this post as the change is at a known time point. If you wanted to test if a change has occurred or not (without specifying the time point) then you could use changepoint analysis to do this.</p> https://stats.stackexchange.com/questions/87421/-/87635#87635 0 Answer by adunaic for Checking that values are piecewise uniform adunaic https://stats.stackexchange.com/users/13409 2014-02-24T09:49:43Z 2014-02-24T09:49:43Z <p>In changepoints if you want to assume a distribution for each segment then you can test if there is a changepoint (or multiple changepoints) using a likelihood ratio test.</p> <p>So in your example you would fit a Uniform distribution to the pre-change and post-change observations separately, calculate the likelihood and repeat this for each possible changepoint location. Then the most likely changepoint location is the one with the highest likelihood. You then test this by comparing the ratio between the null (no change) likelihood and alternative (one change at the maximized location) likelihood to a threshold.</p> <p>The above is the general likelihood setup. I've not come across this done for a Uniform distribution before, but give it a google and see what you find.</p> https://stats.stackexchange.com/questions/84111/-/84192#84192 1 Answer by adunaic for Change detection for beginners adunaic https://stats.stackexchange.com/users/13409 2014-02-02T20:12:50Z 2014-02-02T20:12:50Z <p>The most powerful tests are based on parametric assumptions. If you can assume a pre-change and post-change distribution then you can construct a standard likelihood ratio test. Whether the test can identify a change as small as the one suggested would require testing - even if it can you may get many false alarms. The only way to know is to conduct a simulation study for 1 change (maybe with 1,000 observations pre and post change) to check the sensitivity.</p> <p>Then if you do have a test that is sensitive enough you can generalize the single change to multiple changes using an exact search approach such as PELT or Segment Neighbourhood from the <code>changepoint</code> or <code>strucchange</code> R package. If you require help with this part then give me a shout.</p> https://stats.stackexchange.com/questions/77131/-/77233#77233 1 Answer by adunaic for How to analyze data which might come from a few normal distribution concatenate together in order? adunaic https://stats.stackexchange.com/users/13409 2013-11-21T09:32:00Z 2013-11-21T09:32:00Z <p>Please be careful when using the Binary Segmentation flag in the changepoint package. If you are not familiar with changepoints then you may not be aware that this gives an <em>approximate</em> solution to the changepoint search - meaning that a more accurate segmentation may exist. I would advise using <code>method='PELT'</code> as this provides an exact search approach and thus the most accurate segmentation.</p> <p>You should also be aware that for changes in mean it is hard to detect a change if (mu_1-mu_2)/sigma &lt; 0.7. By hard I mean that the likelihood ratio test has power less than 0.8.</p> https://stats.stackexchange.com/questions/362126/changepoint-step-detection-in-univariate-time-series/362690?cid=681338#362690 Comment by adunaic on Changepoint/Step Detection in Univariate Time Series adunaic https://stats.stackexchange.com/users/13409 2018-08-17T14:27:16Z 2018-08-17T14:27:16Z The downside is that although the slope is constant across the segments, the code fits a new slope after each changepoint. https://stats.stackexchange.com/questions/357687/finding-the-change-point-before-a-significant-increase/357877?cid=673421#357877 Comment by adunaic on Finding the change point before a significant increase adunaic https://stats.stackexchange.com/users/13409 2018-07-20T11:48:36Z 2018-07-20T11:48:36Z The approach above considers the entire data with no choice of window size and the ability to incorporate a minimum distance between the changes through the minseglen argument. The idea is to make it as data driven as possible. https://stats.stackexchange.com/questions/357687/finding-the-change-point-before-a-significant-increase/357877?cid=673420#357877 Comment by adunaic on Finding the change point before a significant increase adunaic https://stats.stackexchange.com/users/13409 2018-07-20T11:47:27Z 2018-07-20T11:47:27Z You can do a moving window test but then you have to choose the window size that you use. Potentially you could get different changes with different window sizes and then how would you make inference? https://stats.stackexchange.com/questions/348777/finding-the-likelihood-of-a-normal-distribution-in-a-change-point-model?cid=657988 Comment by adunaic on Finding the likelihood of a normal distribution in a change point model adunaic https://stats.stackexchange.com/users/13409 2018-05-30T20:27:01Z 2018-05-30T20:27:01Z This is a classic homework question and so the full answer should not be given. https://stats.stackexchange.com/questions/300377/identify-two-major-change-points-in-time-series-data-and-summarise-for-several-r/300579?cid=572373#300579 Comment by adunaic on Identify two major change points in time series data and summarise for several replicates in R adunaic https://stats.stackexchange.com/users/13409 2017-09-02T18:33:35Z 2017-09-02T18:33:35Z Oh and I suppose you can average the first and second changes if you want to but i&#39;m not sure what that tells you as the average probably isn&#39;t optimal in any sense. https://stats.stackexchange.com/questions/300377/identify-two-major-change-points-in-time-series-data-and-summarise-for-several-r/300579?cid=572372#300579 Comment by adunaic on Identify two major change points in time series data and summarise for several replicates in R adunaic https://stats.stackexchange.com/users/13409 2017-09-02T18:32:43Z 2017-09-02T18:32:43Z It is unfortunate that the InspectChangepoint package didn&#39;t work out of the box, have you tried changing the lambda value to be larger? As lambda increases you get less changes. You also might want to set M=100 or something similar as M=0 is the default which uses the Binary Segmentation algorithm whereas M&gt;0 uses the Wild Binary Segmentation algorithm which is preferred. https://stats.stackexchange.com/questions/297239/linear-regression-getting-worse-results-on-training-set-with-additional-paramete?cid=565077 Comment by adunaic on Linear regression getting worse results on training set with additional parameters adunaic https://stats.stackexchange.com/users/13409 2017-08-10T12:33:30Z 2017-08-10T12:33:30Z What is the length of your feature vectors? i.e. what is the size of the problem? https://stats.stackexchange.com/questions/59755/change-point-analysis/59843?cid=565045#59843 Comment by adunaic on Change point analysis adunaic https://stats.stackexchange.com/users/13409 2017-08-10T10:19:20Z 2017-08-10T10:19:20Z You can always write the code yourself by implementing the methods in any of the publications. The question is why you would want to do this when someone has already done it for you? Unless it isn&#39;t available in the language of your choosing that is. https://stats.stackexchange.com/questions/283681/changepoint-detection-in-a-single-unevenly-spaced-time-series?cid=543112 Comment by adunaic on Changepoint detection in a single unevenly spaced time series adunaic https://stats.stackexchange.com/users/13409 2017-06-06T08:57:19Z 2017-06-06T08:57:19Z What are you looking for changes in? If it is the trend then unevenly spaced data matters more than for changes in a flat mean. https://stats.stackexchange.com/questions/271084/how-do-i-split-a-vector-and-minimize-absolute-error-efficiently/271283?cid=520805#271283 Comment by adunaic on How do I split a vector and minimize absolute error efficiently? adunaic https://stats.stackexchange.com/users/13409 2017-04-06T09:13:25Z 2017-04-06T09:13:25Z Apologies, I misread the MAE part. I think you should still try the package as the use of the EDF with a small number of quantiles reduces the effect of outliers. https://stats.stackexchange.com/questions/726/famous-statistical-quotations/1330?cid=457151#1330 Comment by adunaic on Famous statistical quotations adunaic https://stats.stackexchange.com/users/13409 2016-10-14T19:58:58Z 2016-10-14T19:58:58Z This is precisely why I got into the field, i&#39;m nosy! https://stats.stackexchange.com/questions/238134/what-is-change-point-analysis-for-at-most-one-changeamoc/238360?cid=456568#238360 Comment by adunaic on What is change point analysis for At Most One Change(AMOC)? adunaic https://stats.stackexchange.com/users/13409 2016-10-13T15:13:53Z 2016-10-13T15:13:53Z In contrast the PELT algorithm solve the changepoint placement optimization problem exactly. In the SegNeigh option (not mentioned above) all the placement combinations are considered but as you can imagine that is computationally expensive. The PELT algorithm solves the problem exactly but in a smart way so it doesn&#39;t evaluate many non-optimal combinations. The paper given in the references of the package gives more details. https://stats.stackexchange.com/questions/238134/what-is-change-point-analysis-for-at-most-one-changeamoc/238360?cid=456566#238360 Comment by adunaic on What is change point analysis for At Most One Change(AMOC)? adunaic https://stats.stackexchange.com/users/13409 2016-10-13T15:13:31Z 2016-10-13T15:13:31Z You can loop AMOC for multiple changes by dividing the data at the changepoint - that is precisely what BinSeg (Binary Segmentation) does. This is an approximate algorithm as the changepoints identified are conditional on the ones found at the previous steps. https://stats.stackexchange.com/questions/238134/what-is-change-point-analysis-for-at-most-one-changeamoc/238360?cid=453996#238360 Comment by adunaic on What is change point analysis for At Most One Change(AMOC)? adunaic https://stats.stackexchange.com/users/13409 2016-10-06T08:36:00Z 2016-10-06T08:36:00Z If you have or suspect multiple changes may be present then you should use a multiple changepoint algorithm such as PELT or BinSeg which you can specify using method=&quot;PELT&quot; for example. By definition of the test statistic the most significant would be selected but what significant means to a statistical algorithm and what significant means to the human eye are different. It may not correspond to the largest change (think about a large change followed quickly by smaller ones and you can create a counterexample). https://stats.stackexchange.com/questions/230772/switchpoint-calculation-in-linear-models-with-probabilistic-programming/231077?cid=453165#231077 Comment by adunaic on Switchpoint calculation in linear models with probabilistic programming adunaic https://stats.stackexchange.com/users/13409 2016-10-04T12:06:39Z 2016-10-04T12:06:39Z If you are still looking for regressions that are connected at the breakpoint then you might want to email Paul Fearnhead as a PhD student of his has just finished his thesis on precisely this problem. He will likely share the code: <a href="http://maths.lancs.ac.uk/~fearnhea/" rel="nofollow noreferrer">maths.lancs.ac.uk/~fearnhea</a>