Last updated: 2018-10-18
workflowr checks: (Click a bullet for more information) ✔ R Markdown file: up-to-date
Great! Since the R Markdown file has been committed to the Git repository, you know the exact version of the code that produced these results.
✔ Environment: empty
Great job! The global environment was empty. Objects defined in the global environment can affect the analysis in your R Markdown file in unknown ways. For reproduciblity it’s best to always run the code in an empty environment.
✔ Seed:
set.seed(20180414)
The command set.seed(20180414)
was run prior to running the code in the R Markdown file. Setting a seed ensures that any results that rely on randomness, e.g. subsampling or permutations, are reproducible.
✔ Session information: recorded
Great job! Recording the operating system, R version, and package versions is critical for reproducibility.
✔ Repository version: 04dec25
wflow_publish
or wflow_git_commit
). workflowr only checks the R Markdown file, but you know if there are other scripts or data files that it depends on. Below is the status of the Git repository when the results were generated:
Ignored files:
Ignored: .DS_Store
Ignored: .Rhistory
Ignored: .Rproj.user/
Ignored: analysis/.Rhistory
Untracked files:
Untracked: analysis/null.Rmd
Untracked: analysis/test.Rmd
Untracked: data/geneMatrix.tsv
Untracked: data/liter_data_4_summarize_ld_1_lm_less_3.rds
Untracked: data/meta.tsv
Untracked: docs/figure/test.Rmd/
Note that any generated files, e.g. HTML, png, CSS, etc., are not included in this status report because it is ok for generated content to have uncommitted changes.
File | Version | Author | Date | Message |
---|---|---|---|---|
Rmd | 04dec25 | stephens999 | 2018-10-18 | workflowr::wflow_publish(“analysis/changepoint.Rmd”) |
Compare susier on some example change point problems and the changepoint package.
This example comes from Killick and Eckley
library(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.
set.seed(10)
m.data=c(rnorm(100,0,1),rnorm(100,1,1),rnorm(100,0,1),rnorm(100,0.2,1))
ts.plot(m.data,xlab="Index")
Here we apply susier to this example. It finds 2 (out of the three) changepoints.
library("susieR")
susie_cp = function(y,...){
n=length(y)
X = matrix(0,nrow=n,ncol=n-1)
for(j in 1:(n-1)){
for(i in (j+1):n){
X[i,j] = 1
}
}
s = susie(X,y,min_abs_corr=0.9,...)
return(s)
}
s = susie_cp(m.data)
ts.plot(m.data,xlab="Index")
lines(predict(s),col=2,lwd=2)
Note that susieR can provide credible sets for the changepoints. Here it finds two, shown in this plot:
#plot a time series y with confidence sets from susie fit s overlaid
# does +- 0.5 so that singletons show up
plot_cp = function(s,y){
library("ggplot2")
df<-data.frame(x = 1:length(y),y = y)
CS = s$sets$cs
p= ggplot(df) + geom_line(mapping=aes_string(x="x", y="y"))
for(i in 1:length(CS)){
p = p+ annotate("rect", fill = "red", alpha = 0.5,
xmin = min(CS[[i]])-0.5, xmax = max(CS[[i]])+0.5,
ymin = -Inf, ymax = Inf)
}
p
}
plot_cp(s,m.data)
data(Lai2005fig4)
Lai.default=cpt.mean(Lai2005fig4[,5],method="PELT")
plot(Lai.default,pch=20,col="grey",cpt.col="black",type="p",xlab="Index")
cpts(Lai.default)
[1] 81 85 89 96 123 133
coef(Lai.default)
$mean
[1] 0.2468910 4.6699210 0.4495538 4.5902489 0.2079891 4.2913844 0.2291286
We see in this case susie seems to “miss” one of the changepoints.
fitted = susie_cp(Lai2005fig4[,5])
plot(Lai2005fig4[,5])
lines(predict(fitted),col=2)
See if this is maybe a convergence issue solved by susie_auto
:
susie_cp_auto = function(y,...){
n=length(y)
X = matrix(0,nrow=n,ncol=n-1)
for(j in 1:(n-1)){
for(i in (j+1):n){
X[i,j] = 1
}
}
s = susie_auto(X,y,min_abs_corr=0.9,...)
return(s)
}
fitted_a = susie_cp_auto(Lai2005fig4[,5])
plot(Lai2005fig4[,5])
lines(predict(fitted_a),col=2)
Plot confidence sets:
plot_cp(fitted_a,Lai2005fig4[,5])
sessionInfo()
R version 3.5.1 (2018-07-02)
Platform: x86_64-apple-darwin15.6.0 (64-bit)
Running under: OS X El Capitan 10.11.6
Matrix products: default
BLAS: /Library/Frameworks/R.framework/Versions/3.5/Resources/lib/libRblas.0.dylib
LAPACK: /Library/Frameworks/R.framework/Versions/3.5/Resources/lib/libRlapack.dylib
locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] ggplot2_3.0.0 susieR_0.5.0.0347 changepoint_2.2.2 zoo_1.8-4
loaded via a namespace (and not attached):
[1] Rcpp_0.12.19 bindr_0.1.1 compiler_3.5.1
[4] pillar_1.3.0 git2r_0.23.0 plyr_1.8.4
[7] workflowr_1.1.1 R.methodsS3_1.7.1 R.utils_2.7.0
[10] tools_3.5.1 digest_0.6.18 evaluate_0.12
[13] tibble_1.4.2 gtable_0.2.0 lattice_0.20-35
[16] pkgconfig_2.0.2 rlang_0.2.2 Matrix_1.2-14
[19] yaml_2.2.0 bindrcpp_0.2.2 withr_2.1.2
[22] stringr_1.3.1 dplyr_0.7.7 knitr_1.20
[25] tidyselect_0.2.5 rprojroot_1.3-2 grid_3.5.1
[28] glue_1.3.0 R6_2.3.0 rmarkdown_1.10
[31] purrr_0.2.5 magrittr_1.5 whisker_0.3-2
[34] backports_1.1.2 scales_1.0.0 matrixStats_0.54.0
[37] htmltools_0.3.6 assertthat_0.2.0 colorspace_1.3-2
[40] labeling_0.3 stringi_1.2.4 lazyeval_0.2.1
[43] munsell_0.5.0 crayon_1.3.4 R.oo_1.22.0
This reproducible R Markdown analysis was created with workflowr 1.1.1