A Base R, simple implementation of the No-Underrun Sampler. This package aims to mostly directly implement the algorithm as described by th
Setting up package
+234
-48
+8
-3
DESCRIPTION
+8
-3
DESCRIPTION
···
1
1
Package: R.NURS
2
-
Title: What the Package Does (One Line, Title Case)
2
+
Title: R Implementation of No-Underrun Sampler
3
3
Version: 0.0.0.9000
4
4
Authors@R:
5
-
person("First", "Last", , "first.last@example.com", role = c("aut", "cre"))
6
-
Description: What the package does (one paragraph).
5
+
person("Visruth", "Srimath Kandali", , "public@visruth.com", role = c("aut", "cre"))
6
+
Description: A bare-bones R implementation of the No-Underrun Sampler
7
+
(Nawaf Bou-Rabee, Bob Carpenter, Sifan Liu, Stefan Oberdörster. 2025.)
8
+
written in base R. This implementation attempts to stick very closely
9
+
to the algorithm described in Appendix A of the paper.
7
10
License: MIT + file LICENSE
11
+
Suggests:
12
+
matrixStats
8
13
Encoding: UTF-8
9
14
Roxygen: list(markdown = TRUE)
10
15
RoxygenNote: 7.3.2
+30
-3
R/NURS.R
+30
-3
R/NURS.R
···
1
1
log_sum_exp <- function(log_vals) {
2
+
# if matrixStats is available, rely on C version of logSumExp
2
3
if (requireNamespace("matrixStats", quietly = TRUE)) {
3
4
matrixStats::logSumExp(log_vals)
4
5
} else {
···
7
8
}
8
9
}
9
10
10
-
# No underrun condition
11
+
#' Check no underrun stopping condition
12
+
#'
13
+
#' @param log_vals orbit
14
+
#' @param log_eps_h log(epsilon) + log(h)
15
+
#' @returns True if Orbit satisfies stopping criterion.
11
16
NURS_stop <- function(log_vals, log_eps_h) {
12
17
max(log_vals[1], log_vals[length(log_vals)]) <=
13
18
log_eps_h + log_sum_exp(log_vals)
14
19
}
15
20
16
-
# Recursive sub stopping
21
+
#' Recursively check sub stopping
22
+
#'
23
+
#' @param log_vals sub orbit
24
+
#' @param log_eps_h log(epsilon) + log(h)
25
+
#' @returns True if sub orbit satisfies stopping criterion.
17
26
NURS_sub_stop <- function(log_vals, log_eps_h) {
18
27
n <- length(log_vals)
19
28
if (n < 2) return(FALSE)
···
23
32
NURS_sub_stop(log_vals[-left_indices], log_eps_h)
24
33
}
25
34
35
+
#' Single NURS step
36
+
#'
37
+
#' @param logpdf log (non-normalized) target density
38
+
#' @param theta current state
39
+
#' @param epsilon density threshold
40
+
#' @param h lattice size
41
+
#' @param M maximum number of doublings
42
+
#' @returns One NURS step from theta.
26
43
NURS_step <- function(logpdf, theta, epsilon, h, M) {
27
44
d <- length(theta)
28
45
# hit
···
81
98
)]]
82
99
}
83
100
101
+
#' NURS draws
102
+
#'
103
+
#' @param logpdf log (non-normalized) target density
104
+
#' @param theta_init initial state
105
+
#' @param n number of draws
106
+
#' @param epsilon density threshold
107
+
#' @param h lattice size
108
+
#' @param M maximum number of doublings
109
+
#' @export
110
+
#' @source <https://arxiv.org/abs/2501.18548v2>
84
111
NURS <- function(logpdf, theta_init, n, epsilon, h, M) {
85
112
d <- length(theta_init)
86
-
draws <- matrix(NA_real_, n, d)
113
+
draws <- matrix(NA, n, d)
87
114
draws[1, ] <- theta_init
88
115
for (i in 2:n) {
89
116
draws[i, ] <- NURS_step(logpdf, draws[i - 1, ], epsilon, h, M)
+56
README.Rmd
+56
README.Rmd
···
1
+
---
2
+
output: github_document
3
+
---
4
+
5
+
<!-- README.md is generated from README.Rmd. Please edit that file -->
6
+
7
+
```{r, include = FALSE}
8
+
knitr::opts_chunk$set(
9
+
collapse = TRUE,
10
+
comment = "#>",
11
+
fig.path = "man/figures/README-",
12
+
out.width = "100%"
13
+
)
14
+
```
15
+
16
+
# R.NURS
17
+
18
+
<!-- badges: start -->
19
+
<!-- badges: end -->
20
+
21
+
Base R implementation of the No-Underrun Sampler.
22
+
23
+
## Installation
24
+
25
+
You can install the development version of R.NURS from [GitHub](https://github.com/) with:
26
+
27
+
``` r
28
+
# install.packages("pak")
29
+
pak::pak("VisruthSK/R.NURS")
30
+
```
31
+
32
+
## Example: Neal's Funnel
33
+
34
+
```{r funnel}
35
+
library(R.NURS)
36
+
library(ggplot2)
37
+
38
+
set.seed(0)
39
+
logpdf_funnel <- function(theta) {
40
+
y <- theta[1]
41
+
dnorm(y, 0, 3, log = TRUE) + sum(dnorm(theta[-1], 0, exp(y / 2), log = TRUE))
42
+
}
43
+
samples <- NURS(
44
+
logpdf_funnel,
45
+
theta_init = rep(0, 15),
46
+
n = 5000,
47
+
epsilon = 0.1,
48
+
h = 0.5,
49
+
M = 5
50
+
)
51
+
52
+
data.frame(y = samples[, 1], x1 = samples[, 2]) |>
53
+
ggplot(aes(x = y, y = x1)) +
54
+
geom_point(alpha = 0.3) +
55
+
theme_minimal()
56
+
```
+48
-42
README.md
+48
-42
README.md
···
1
-
2
-
# R.NURS
3
-
4
-
Base R implementation of the No-Underrun Sampler.
5
-
6
-
## Installation
7
-
8
-
You can install the development version of R.NURS like so:
9
-
10
-
``` r
11
-
# install.packages("pak")
12
-
pak::pak("VisruthSK/R.NURS")
13
-
```
14
-
15
-
## Example
16
-
17
-
This is a basic example which shows you how to solve a common problem:
18
-
19
-
``` r
20
-
library(R.NURS)
21
-
library(ggplot2)
22
-
23
-
set.seed(0)
24
-
logpdf_funnel <- function(theta) {
25
-
y <- theta[1]
26
-
dnorm(y, 0, 3, log = TRUE) + sum(dnorm(theta[-1], 0, exp(y / 2), log = TRUE))
27
-
}
28
-
samples <- NURS(
29
-
logpdf_funnel,
30
-
theta_init = rep(0, 15),
31
-
n = 5000,
32
-
epsilon = 0.1,
33
-
h = 0.5,
34
-
M = 5
35
-
)
36
-
37
-
data.frame(y = samples[, 1], x1 = samples[, 2]) |>
38
-
ggplot(aes(x = y, y = x1)) +
39
-
geom_point(alpha = 0.3) +
40
-
theme_minimal()
41
-
```
42
-
1
+
2
+
<!-- README.md is generated from README.Rmd. Please edit that file -->
3
+
4
+
# R.NURS
5
+
6
+
<!-- badges: start -->
7
+
8
+
<!-- badges: end -->
9
+
10
+
Base R implementation of the No-Underrun Sampler.
11
+
12
+
## Installation
13
+
14
+
You can install the development version of R.NURS from
15
+
[GitHub](https://github.com/) with:
16
+
17
+
``` r
18
+
# install.packages("pak")
19
+
pak::pak("VisruthSK/R.NURS")
20
+
```
21
+
22
+
## Example: Neal’s Funnel
23
+
24
+
``` r
25
+
library(R.NURS)
26
+
library(ggplot2)
27
+
28
+
set.seed(0)
29
+
logpdf_funnel <- function(theta) {
30
+
y <- theta[1]
31
+
dnorm(y, 0, 3, log = TRUE) + sum(dnorm(theta[-1], 0, exp(y / 2), log = TRUE))
32
+
}
33
+
samples <- NURS(
34
+
logpdf_funnel,
35
+
theta_init = rep(0, 15),
36
+
n = 5000,
37
+
epsilon = 0.1,
38
+
h = 0.5,
39
+
M = 5
40
+
)
41
+
42
+
data.frame(y = samples[, 1], x1 = samples[, 2]) |>
43
+
ggplot(aes(x = y, y = x1)) +
44
+
geom_point(alpha = 0.3) +
45
+
theme_minimal()
46
+
```
47
+
48
+
<img src="man/figures/README-funnel-1.png" width="100%" />
+27
man/NURS.Rd
+27
man/NURS.Rd
···
1
+
% Generated by roxygen2: do not edit by hand
2
+
% Please edit documentation in R/NURS.R
3
+
\name{NURS}
4
+
\alias{NURS}
5
+
\title{NURS draws}
6
+
\source{
7
+
\url{https://arxiv.org/abs/2501.18548v2}
8
+
}
9
+
\usage{
10
+
NURS(logpdf, theta_init, n, epsilon, h, M)
11
+
}
12
+
\arguments{
13
+
\item{logpdf}{log (non-normalized) target density}
14
+
15
+
\item{theta_init}{initial state}
16
+
17
+
\item{n}{number of draws}
18
+
19
+
\item{epsilon}{density threshold}
20
+
21
+
\item{h}{lattice size}
22
+
23
+
\item{M}{maximum number of doublings}
24
+
}
25
+
\description{
26
+
NURS draws
27
+
}
+25
man/NURS_step.Rd
+25
man/NURS_step.Rd
···
1
+
% Generated by roxygen2: do not edit by hand
2
+
% Please edit documentation in R/NURS.R
3
+
\name{NURS_step}
4
+
\alias{NURS_step}
5
+
\title{Single NURS step}
6
+
\usage{
7
+
NURS_step(logpdf, theta, epsilon, h, M)
8
+
}
9
+
\arguments{
10
+
\item{logpdf}{log (non-normalized) target density}
11
+
12
+
\item{theta}{current state}
13
+
14
+
\item{epsilon}{density threshold}
15
+
16
+
\item{h}{lattice size}
17
+
18
+
\item{M}{maximum number of doublings}
19
+
}
20
+
\value{
21
+
One NURS step from theta.
22
+
}
23
+
\description{
24
+
Single NURS step
25
+
}
+19
man/NURS_stop.Rd
+19
man/NURS_stop.Rd
···
1
+
% Generated by roxygen2: do not edit by hand
2
+
% Please edit documentation in R/NURS.R
3
+
\name{NURS_stop}
4
+
\alias{NURS_stop}
5
+
\title{Check no underrun stopping condition}
6
+
\usage{
7
+
NURS_stop(log_vals, log_eps_h)
8
+
}
9
+
\arguments{
10
+
\item{log_vals}{orbit}
11
+
12
+
\item{log_eps_h}{log(epsilon) + log(h)}
13
+
}
14
+
\value{
15
+
True if Orbit satisfies stopping criterion.
16
+
}
17
+
\description{
18
+
Check no underrun stopping condition
19
+
}
+19
man/NURS_sub_stop.Rd
+19
man/NURS_sub_stop.Rd
···
1
+
% Generated by roxygen2: do not edit by hand
2
+
% Please edit documentation in R/NURS.R
3
+
\name{NURS_sub_stop}
4
+
\alias{NURS_sub_stop}
5
+
\title{Recursively check sub stopping}
6
+
\usage{
7
+
NURS_sub_stop(log_vals, log_eps_h)
8
+
}
9
+
\arguments{
10
+
\item{log_vals}{sub orbit}
11
+
12
+
\item{log_eps_h}{log(epsilon) + log(h)}
13
+
}
14
+
\value{
15
+
True if sub orbit satisfies stopping criterion.
16
+
}
17
+
\description{
18
+
Recursively check sub stopping
19
+
}
man/figures/README-funnel-1.png
man/figures/README-funnel-1.png
This is a binary file and will not be displayed.