plot.envcpt: Plots optionally either ("fit") the data and fits from each...
In EnvCpt: Detection of Structural Changes in Climate and Environment Time Series
Description
Usage
Arguments
Details
Value
Author(s)
See Also
Examples
Description
Uses the output of the envcpt function and plots optionally ("fit") the original data and the fit from each of the 8 models or ("aic") the aic for each of the 8 models.
Usage
1
2
Arguments
x
A list produced as output from the envcpt function. In essence a list with the named elements: "mean","meancpt","meanar1","meanar2","meanar1cpt","meanar2cpt","trend","trendcpt","trendar1","trendar2","trendar1cpt","trendar2cpt" where each element is the output from the appropriate function call.
type
character vector. fit returns a plot of the data and fitted models. aic returns a bar chart of the aic values for each of the 8 models
lwd
Line width graphical parameter, see par for further details.
colors
colors for the individual models to be passed to the plotting functions. Note that this must be of length 12 regardless of how many models were evaluated. Color validity is checked using col2rgb.
...
Extra graphical parameters, passed to the original plot and the individual calls to lines or barchar.
data
This argument is only required when x doesn't contain any fits that provide the data which is required for the "fit" type option. Default is NA as most of the time this is likely not required.
Details
If type="fit", the function plots the data at the bottom and stacks the different fits for the 8 models from the envcpt function on top. No scale is given as all data and fits are scaled to be in (0,1). This is designed as an initial visualization tool for the fits only.
If type="aic" the function uses the AIC.envcpt function to calculate the AIC values for the envcpt output x. Then barcharts the AIC values in the same order as the type="fit" option. The minimum AIC is the preferred model and this is highlight by a solid block. This is designed as an initial visualization tool for the AIC values only.
If type="bic" the function uses the BIC.envcpt function to calculate the BIC values for the envcpt output x. Then barcharts the BIC values in the same order as the type="fit" option. The minimum BIC is the preferred model and this is highlight by a solid block. This is designed as an initial visualization tool for the BIC values only.
For ease of use, if "colours" is specified instead of (or in addition to) colors then this will take precedence.
Value
Returns the printed graphic to the active device.
Author(s)
Rebecca Killick & Claudie Beaulieu.
See Also
Examples
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
## Not run:
set.seed(1)
x=c(rnorm(100,0,1),rnorm(100,5,1))
out=envcpt(x) # run all models with default values
out[[1]] # first row is twice the negative log-likelihood for each model
# second row is the number of parameters
AIC(out) # returns AIC for each model.
which.min(AIC(out)) # gives meancpt (model 2) as the best model fit.
out$meancpt # gives the model fit for the meancpt model.
AICweights(out) # gives the AIC weights for each model
BIC(out) # returns the BIC for each model.
which.min(BIC(out)) # gives meancpt (model 2) as the best model fit too.
plot(out,type='fit') # plots the fits
plot(out,type="aic") # plots the aic values
plot(out,type="bic") # plots the bic values
set.seed(10)
x=c(0.01*(1:100),1.5-0.02*((101:250)-101))+rnorm(250,0,0.2)
out=envcpt(x,minseglen=10) # run all models with a minimum of 10 observations between changes
AIC(out) # returns the AIC for each model
which.min(AIC(out)) # gives trendcpt (model 8) as the best model fit.
out$trendcpt # gives the model fit for the trendcpt model.
AICweights(out) # gives the AIC weights for each model
BIC(out) # returns the BIC for each model.
which.min(BIC(out)) # gives trendcpt (model 8) as the best model fit too.
plot(out,type='fit') # plots the fits
plot(out,type="aic") # plots the aic values
plot(out,type="bic") # plots the bic values
set.seed(100)
x=arima.sim(model=list(ar=c(0.7,0.2)),n=500)+0.01*(1:500)
out=envcpt(x,models=c(3:6,9:12)) # runs a subset of models (those with AR components)
AIC(out) # returns the AIC for each model
which.min(AIC(out)) # gives trendar2 (model 10) as the best model fit.
out$trendar2 # gives the model fit for the trendar2 model. Notice that the trend is tiny but does
# produce a significantly better fit than the meanar2 model.
AICweights(out) # gives the AIC weights for each model
BIC(out) # returns the BIC for each model.
which.min(BIC(out)) # best fit is trendar2 (model 10) again.
plot(out,type='fit') # plots the fits
plot(out,type="aic") # plots the aic values
plot(out,type="bic") # plots the bic values
## End(Not run)
Example output
Loading required package: changepoint
Loading required package: zoo
Attaching package: 'zoo'
The following objects are masked from 'package:base':
as.Date, as.Date.numeric
Successfully loaded changepoint package version 2.2.2
NOTE: Predefined penalty values changed in version 2.2. Previous penalty values with a postfix 1 i.e. SIC1 are now without i.e. SIC and previous penalties without a postfix i.e. SIC are now with a postfix 0 i.e. SIC0. See NEWS and help files for further details.
Loading required package: MASS
Fitting 12 models
|
| | 0%
|
|====== | 8%
|
|============ | 17%
|
|================== | 25%
|
|======================= | 33%
|
|============================= | 42%
|
|=================================== | 50%
|
|========================================= | 58%
|
|=============================================== | 67%
|
|==================================================== | 75%
|
|========================================================== | 83%
|
|================================================================ | 92%
|
|======================================================================| 100%
mean meancpt meanar1 meanar2 meanar1cpt meanar2cpt trend
logLik 949.2517 535.486 683.8642 640.5247 532.9724 530.5236 735.4193
nparam 2.0000 5.000 3.0000 4.0000 7.0000 9.0000 3.0000
trendcpt trendar1 trendar2 trendar1cpt trendar2cpt
logLik 535.414 648.192 619.5107 532.7957 530.4183
nparam 7.000 4.000 5.0000 9.0000 11.0000
mean meancpt meanar1 meanar2 meanar1cpt meanar2cpt
953.2517 545.4860 689.8642 648.5247 546.9724 548.5236
trend trendcpt trendar1 trendar2 trendar1cpt trendar2cpt
741.4193 549.4140 656.1920 629.5107 550.7957 552.4183
meancpt
2
Class 'cpt' : Changepoint Object
~~ : S4 class containing 12 slots with names
cpttype date version data.set method test.stat pen.type pen.value minseglen cpts ncpts.max param.est
Created on : Thu May 4 11:08:19 2017
summary(.) :
----------
Created Using changepoint version 2.2.2
Changepoint type : Change in mean and variance
Method of analysis : PELT
Test Statistic : Normal
Type of penalty : MBIC with value, 21.19327
Minimum Segment Length : 5
Maximum no. of cpts : Inf
Changepoint Locations : 100
mean meancpt meanar1 meanar2 meanar1cpt meanar2cpt
1.471637e-89 5.164235e-01 2.299732e-32 2.179893e-23 2.456047e-01 1.130823e-01
trend trendcpt trendar1 trendar2 trendar1cpt trendar2cpt
1.467638e-43 7.245090e-02 4.715216e-25 2.932684e-19 3.630777e-02 1.613078e-02
mean meancpt meanar1 meanar2 meanar1cpt meanar2cpt
959.8483 561.9775 699.7591 661.7179 570.0606 578.2084
trend trendcpt trendar1 trendar2 trendar1cpt trendar2cpt
751.3143 572.5022 669.3852 646.0023 580.4806 588.6998
meancpt
2
Fitting 12 models
|
| | 0%
|
|====== | 8%
|
|============ | 17%
|
|================== | 25%
|
|======================= | 33%
|
|============================= | 42%
|
|=================================== | 50%
|
|========================================= | 58%
|
|=============================================== | 67%
|
|==================================================== | 75%
|
|========================================================== | 83%
|
|================================================================ | 92%
|
|======================================================================| 100%
mean meancpt meanar1 meanar2 meanar1cpt meanar2cpt
578.94707 -67.51376 27.29077 -29.55647 27.29077 29.29077
trend trendcpt trendar1 trendar2 trendar1cpt trendar2cpt
468.36581 -113.41437 20.92805 -33.70269 -109.94452 -105.29104
trendcpt
8
Class 'cpt.reg' : Changepoint Regression Object
~~ : S4 class containing 12 slots with names
cpttype date version data.set method test.stat pen.type pen.value minseglen cpts ncpts.max param.est
Created on : Thu May 4 11:08:19 2017
summary(.) :
----------
Created Using changepoint version 2.2.2
Changepoint type : Change in regression
Method of analysis : PELT
Test Statistic : Normal
Type of penalty : MBIC with value, 27.6073
Maximum no. of cpts : Inf
Changepoint Locations : 100
mean meancpt meanar1 meanar2 meanar1cpt
3.790989e-151 9.035235e-11 2.340974e-31 5.171508e-19 2.340974e-31
meanar2cpt trend trendcpt trendar1 trendar2
8.611961e-32 3.900936e-127 8.377773e-01 5.636952e-30 4.111095e-18
trendar1cpt trendar2cpt
1.477958e-01 1.442682e-02
mean meancpt meanar1 meanar2 meanar1cpt meanar2cpt
585.98999 24.04422 37.85515 -15.47062 37.85515 43.37661
trend trendcpt trendar1 trendar2 trendar1cpt trendar2cpt
478.93019 -88.76414 35.01389 -16.09538 -78.25138 -66.55497
trendcpt
8
Fitting 8 models
|
| | 0%
|
|========================== | 38%
|
|=================================== | 50%
|
|============================================ | 62%
|
|==================================================== | 75%
mean meancpt meanar1 meanar2 meanar1cpt meanar2cpt
NA NA 1478.915 1455.634 1478.915 1480.915
trend trendcpt trendar1 trendar2 trendar1cpt trendar2cpt
NA NA 1452.451 1430.101 1452.451 1430.101
trendar2
10
Call:
lm(formula = data[-c(1:2)] ~ c(1:(n - 2)) + data[2:(n - 1)] +
data[1:(n - 2)])
Coefficients:
(Intercept) c(1:(n - 2)) data[2:(n - 1)] data[1:(n - 2)]
-0.069419 0.001703 0.661337 0.186613
mean meancpt meanar1 meanar2 meanar1cpt meanar2cpt
NA NA 1.256735e-11 1.428051e-06 1.256735e-11 4.623269e-12
trend trendcpt trendar1 trendar2 trendar1cpt trendar2cpt
NA NA 7.011284e-06 4.999923e-01 7.011284e-06 4.999923e-01
mean meancpt meanar1 meanar2 meanar1cpt meanar2cpt
NA NA 1491.553 1472.484 1491.553 1497.765
trend trendcpt trendar1 trendar2 trendar1cpt trendar2cpt
NA NA 1469.302 1451.164 1469.302 1451.164
trendar2
10
EnvCpt documentation built on March 29, 2021, 5:11 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Description Usage Arguments Details Value Author(s) See Also Examples
Description
Uses the output of the envcpt function and plots optionally ("fit") the original data and the fit from each of the 8 models or ("aic") the aic for each of the 8 models.
Usage
1 2 |
Arguments
x |
A list produced as output from the |
type |
character vector. |
lwd |
Line width graphical parameter, see |
colors |
colors for the individual models to be passed to the plotting functions. Note that this must be of length 12 regardless of how many models were evaluated. Color validity is checked using col2rgb. |
... |
Extra graphical parameters, passed to the original plot and the individual calls to |
data |
This argument is only required when |
Details
If type="fit", the function plots the data at the bottom and stacks the different fits for the 8 models from the envcpt function on top. No scale is given as all data and fits are scaled to be in (0,1). This is designed as an initial visualization tool for the fits only.
If type="aic" the function uses the AIC.envcpt function to calculate the AIC values for the envcpt output x. Then barcharts the AIC values in the same order as the type="fit" option. The minimum AIC is the preferred model and this is highlight by a solid block. This is designed as an initial visualization tool for the AIC values only.
If type="bic" the function uses the BIC.envcpt function to calculate the BIC values for the envcpt output x. Then barcharts the BIC values in the same order as the type="fit" option. The minimum BIC is the preferred model and this is highlight by a solid block. This is designed as an initial visualization tool for the BIC values only.
For ease of use, if "colours" is specified instead of (or in addition to) colors then this will take precedence.
Value
Returns the printed graphic to the active device.
Author(s)
Rebecca Killick & Claudie Beaulieu.
See Also
Examples
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 | ## Not run:
set.seed(1)
x=c(rnorm(100,0,1),rnorm(100,5,1))
out=envcpt(x) # run all models with default values
out[[1]] # first row is twice the negative log-likelihood for each model
# second row is the number of parameters
AIC(out) # returns AIC for each model.
which.min(AIC(out)) # gives meancpt (model 2) as the best model fit.
out$meancpt # gives the model fit for the meancpt model.
AICweights(out) # gives the AIC weights for each model
BIC(out) # returns the BIC for each model.
which.min(BIC(out)) # gives meancpt (model 2) as the best model fit too.
plot(out,type='fit') # plots the fits
plot(out,type="aic") # plots the aic values
plot(out,type="bic") # plots the bic values
set.seed(10)
x=c(0.01*(1:100),1.5-0.02*((101:250)-101))+rnorm(250,0,0.2)
out=envcpt(x,minseglen=10) # run all models with a minimum of 10 observations between changes
AIC(out) # returns the AIC for each model
which.min(AIC(out)) # gives trendcpt (model 8) as the best model fit.
out$trendcpt # gives the model fit for the trendcpt model.
AICweights(out) # gives the AIC weights for each model
BIC(out) # returns the BIC for each model.
which.min(BIC(out)) # gives trendcpt (model 8) as the best model fit too.
plot(out,type='fit') # plots the fits
plot(out,type="aic") # plots the aic values
plot(out,type="bic") # plots the bic values
set.seed(100)
x=arima.sim(model=list(ar=c(0.7,0.2)),n=500)+0.01*(1:500)
out=envcpt(x,models=c(3:6,9:12)) # runs a subset of models (those with AR components)
AIC(out) # returns the AIC for each model
which.min(AIC(out)) # gives trendar2 (model 10) as the best model fit.
out$trendar2 # gives the model fit for the trendar2 model. Notice that the trend is tiny but does
# produce a significantly better fit than the meanar2 model.
AICweights(out) # gives the AIC weights for each model
BIC(out) # returns the BIC for each model.
which.min(BIC(out)) # best fit is trendar2 (model 10) again.
plot(out,type='fit') # plots the fits
plot(out,type="aic") # plots the aic values
plot(out,type="bic") # plots the bic values
## End(Not run)
|
Example output
Loading required package: changepoint
Loading required package: zoo
Attaching package: 'zoo'
The following objects are masked from 'package:base':
as.Date, as.Date.numeric
Successfully loaded changepoint package version 2.2.2
NOTE: Predefined penalty values changed in version 2.2. Previous penalty values with a postfix 1 i.e. SIC1 are now without i.e. SIC and previous penalties without a postfix i.e. SIC are now with a postfix 0 i.e. SIC0. See NEWS and help files for further details.
Loading required package: MASS
Fitting 12 models
|
| | 0%
|
|====== | 8%
|
|============ | 17%
|
|================== | 25%
|
|======================= | 33%
|
|============================= | 42%
|
|=================================== | 50%
|
|========================================= | 58%
|
|=============================================== | 67%
|
|==================================================== | 75%
|
|========================================================== | 83%
|
|================================================================ | 92%
|
|======================================================================| 100%
mean meancpt meanar1 meanar2 meanar1cpt meanar2cpt trend
logLik 949.2517 535.486 683.8642 640.5247 532.9724 530.5236 735.4193
nparam 2.0000 5.000 3.0000 4.0000 7.0000 9.0000 3.0000
trendcpt trendar1 trendar2 trendar1cpt trendar2cpt
logLik 535.414 648.192 619.5107 532.7957 530.4183
nparam 7.000 4.000 5.0000 9.0000 11.0000
mean meancpt meanar1 meanar2 meanar1cpt meanar2cpt
953.2517 545.4860 689.8642 648.5247 546.9724 548.5236
trend trendcpt trendar1 trendar2 trendar1cpt trendar2cpt
741.4193 549.4140 656.1920 629.5107 550.7957 552.4183
meancpt
2
Class 'cpt' : Changepoint Object
~~ : S4 class containing 12 slots with names
cpttype date version data.set method test.stat pen.type pen.value minseglen cpts ncpts.max param.est
Created on : Thu May 4 11:08:19 2017
summary(.) :
----------
Created Using changepoint version 2.2.2
Changepoint type : Change in mean and variance
Method of analysis : PELT
Test Statistic : Normal
Type of penalty : MBIC with value, 21.19327
Minimum Segment Length : 5
Maximum no. of cpts : Inf
Changepoint Locations : 100
mean meancpt meanar1 meanar2 meanar1cpt meanar2cpt
1.471637e-89 5.164235e-01 2.299732e-32 2.179893e-23 2.456047e-01 1.130823e-01
trend trendcpt trendar1 trendar2 trendar1cpt trendar2cpt
1.467638e-43 7.245090e-02 4.715216e-25 2.932684e-19 3.630777e-02 1.613078e-02
mean meancpt meanar1 meanar2 meanar1cpt meanar2cpt
959.8483 561.9775 699.7591 661.7179 570.0606 578.2084
trend trendcpt trendar1 trendar2 trendar1cpt trendar2cpt
751.3143 572.5022 669.3852 646.0023 580.4806 588.6998
meancpt
2
Fitting 12 models
|
| | 0%
|
|====== | 8%
|
|============ | 17%
|
|================== | 25%
|
|======================= | 33%
|
|============================= | 42%
|
|=================================== | 50%
|
|========================================= | 58%
|
|=============================================== | 67%
|
|==================================================== | 75%
|
|========================================================== | 83%
|
|================================================================ | 92%
|
|======================================================================| 100%
mean meancpt meanar1 meanar2 meanar1cpt meanar2cpt
578.94707 -67.51376 27.29077 -29.55647 27.29077 29.29077
trend trendcpt trendar1 trendar2 trendar1cpt trendar2cpt
468.36581 -113.41437 20.92805 -33.70269 -109.94452 -105.29104
trendcpt
8
Class 'cpt.reg' : Changepoint Regression Object
~~ : S4 class containing 12 slots with names
cpttype date version data.set method test.stat pen.type pen.value minseglen cpts ncpts.max param.est
Created on : Thu May 4 11:08:19 2017
summary(.) :
----------
Created Using changepoint version 2.2.2
Changepoint type : Change in regression
Method of analysis : PELT
Test Statistic : Normal
Type of penalty : MBIC with value, 27.6073
Maximum no. of cpts : Inf
Changepoint Locations : 100
mean meancpt meanar1 meanar2 meanar1cpt
3.790989e-151 9.035235e-11 2.340974e-31 5.171508e-19 2.340974e-31
meanar2cpt trend trendcpt trendar1 trendar2
8.611961e-32 3.900936e-127 8.377773e-01 5.636952e-30 4.111095e-18
trendar1cpt trendar2cpt
1.477958e-01 1.442682e-02
mean meancpt meanar1 meanar2 meanar1cpt meanar2cpt
585.98999 24.04422 37.85515 -15.47062 37.85515 43.37661
trend trendcpt trendar1 trendar2 trendar1cpt trendar2cpt
478.93019 -88.76414 35.01389 -16.09538 -78.25138 -66.55497
trendcpt
8
Fitting 8 models
|
| | 0%
|
|========================== | 38%
|
|=================================== | 50%
|
|============================================ | 62%
|
|==================================================== | 75%
mean meancpt meanar1 meanar2 meanar1cpt meanar2cpt
NA NA 1478.915 1455.634 1478.915 1480.915
trend trendcpt trendar1 trendar2 trendar1cpt trendar2cpt
NA NA 1452.451 1430.101 1452.451 1430.101
trendar2
10
Call:
lm(formula = data[-c(1:2)] ~ c(1:(n - 2)) + data[2:(n - 1)] +
data[1:(n - 2)])
Coefficients:
(Intercept) c(1:(n - 2)) data[2:(n - 1)] data[1:(n - 2)]
-0.069419 0.001703 0.661337 0.186613
mean meancpt meanar1 meanar2 meanar1cpt meanar2cpt
NA NA 1.256735e-11 1.428051e-06 1.256735e-11 4.623269e-12
trend trendcpt trendar1 trendar2 trendar1cpt trendar2cpt
NA NA 7.011284e-06 4.999923e-01 7.011284e-06 4.999923e-01
mean meancpt meanar1 meanar2 meanar1cpt meanar2cpt
NA NA 1491.553 1472.484 1491.553 1497.765
trend trendcpt trendar1 trendar2 trendar1cpt trendar2cpt
NA NA 1469.302 1451.164 1469.302 1451.164
trendar2
10
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.