-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathREADME.Rmd
92 lines (61 loc) · 3 KB
/
README.Rmd
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
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
---
title: "SCPME"
output: github_document
---
```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE, tidy = TRUE)
```
[](https://travis-ci.org/MGallow/SCPME)
[](https://cran.r-project.org/package=SCPME)
## Overview
`SCPME` is an implementation of the methods described in "Shrinking Characteristics of Precision Matrix Estimators" ([link](https://doi.org/10.1093/biomet/asy023)). It estimates a penalized precision matrix via a modified alternating direction method of multipliers (ADMM) algorithm.
<p align="center">
<img src = "https://github.com/MGallow/SCPME/raw/master/vignettes/images/gif.gif"/>
</p>
A (possibly incomplete) list of functions contained in the package can be found below:
* `shrink()` computes the estimated precision matrix
* `data_gen()` data generation function (for convenience)
* `plot.shrink()` produces a heat map or line graph for cross validation errors
See package [website](https://mgallow.github.io/SCPME/) or [manual](https://github.com/MGallow/ADMMsigma/blob/master/SCPME.pdf).
## Installation
```{r, eval = FALSE}
# The easiest way to install is from GitHub:
# install.packages("devtools")
devtools::install_github("MGallow/shrink")
```
If there are any issues/bugs, please let me know: [github](https://github.com/MGallow/SCPME/issues). You can also contact me via my [website](https://mgallow.github.io/). Pull requests are welcome!
## Usage
```{r, message = FALSE}
library(SCPME)
set.seed(123)
# generate data from a sparse oracle precision matrix
# we can use the built-in `data_gen` function
# generate 100 x 5 X data matrix and 100 x 1 Y data matrix
data = data_gen(p = 5, n = 100, r = 1)
# the default regression coefficients are sparse
data$betas
# default oracle precision matrix is also sparse
round(qr.solve(data$SigmaX), 5)
# now suppose we are interested in estimating the precision
# print marginal sample precision matrix for X
# this is perhaps a bad estimate (not sparse)
sample = (nrow(data$X) - 1)/nrow(data$X)*cov(data$X)
round(qr.solve(sample), 5)
# now use SCPME to estimate preicison matrix (omega) assuming sparsity
# note that this is simply a lasso penalized preicision matrix
shrink(data$X, lam = 0.5, crit.cv = "loglik")
# what if we instead assumed sparsity in beta? (print estimated omega)
# recall that beta is a product of marginal precision of X and cov(X, Y)
lam_max = max(abs(t(data$X) %*% data$Y))
(shrink = shrink(data$X, data$Y, B = cov(data$X, data$Y), nlam = 20, lam.max = lam_max))
# print estimated beta
shrink$Z
# we could also assume sparsity in beta AND omega (print estimated omega)
(shrink2 = shrink(data$X, data$Y, B = cbind(cov(data$X, data$Y), diag(ncol(data$X))), nlam = 20, lam.max = 10, lam.min.ratio = 1e-4))
# print estimated beta
shrink2$Z[,1, drop = FALSE]
# produce CV heat map for shrink
plot(shrink, type = "heatmap")
# produce line graph for CV errors for shrink
plot(shrink, type = "line")
```