EnvelopeSize: Envelope Sizing and Optimization
In glmbayes: Bayesian Generalized Linear Models (IID Samples)

EnvelopeSize

R Documentation

Envelope Sizing and Optimization

Description

EnvelopeSize() is the high-level entry point that constructs per-dimension grids and expected draw counts, while EnvelopeOpt() performs the adaptive optimization used when Gridtype = 2.

Usage

EnvelopeSize(a, G1, Gridtype = 2L, n = 1000L, n_envopt = -1,
                    use_opencl = FALSE, verbose = FALSE)

EnvelopeOpt(a1,n,core_cnt=1L)

Arguments

`a`	Numeric vector of diagonal precisions for the log-likelihood (posterior precision is `1 + a_i`).
`G1`	Numeric matrix of candidate grid points (3 * l1).
`Gridtype`	Integer code controlling grid sizing logic: 1 = static threshold test 2 = adaptive optimization via `EnvelopeOpt()` 3 = always three-point grid 4 = always single-point grid
`n`	Integer; number of posterior draws to generate (used for grid sizing).
`n_envopt`	Integer; effective sample size passed to `EnvelopeOpt`. Defaults to -1, which means "use n".
`use_opencl`	Logical; if `TRUE`, attempt GPU acceleration.
`verbose`	Logical; if `TRUE`, print progress messages.
`a1`	Numeric vector of diagonal elements of the data precision matrix (used by `EnvelopeOpt`).
`core_cnt`	Integer; number of OpenCL cores or parallel workers available (default 1). When >1, envelope build cost is scaled down to reflect parallel construction.

Details

These functions implement the grid sizing logic used in envelope construction for rejection sampling. They make use of the theory described in \insertCiteNygren2006glmbayes and the general implementation outlined in \insertCiteglmbayesSimmethodsglmbayes.

EnvelopeSize() returns the constructed grid (G2), index vectors (GIndex1), expected draw count (E_draws), and the per-dimension grid index.

EnvelopeOpt() implements the adaptive optimization used in Gridtype = 2, ranking dimensions by posterior variance and promoting them to three-point tangents when the tradeoff is favorable.

Value

EnvelopeSize(): A list with components G2, GIndex1, E_draws, and gridindex.
EnvelopeOpt(): An integer vector of length l1 with entries 1 (single-point) or 3 (three-point).

Gridtype Logic and Candidates per Draw

The envelope sizing logic follows the analysis of \insertCiteNygren2006glmbayes.

Gridtype 1: Static Threshold

For each dimension i, if \sqrt{1 + a_i} \leq 2/\sqrt{\pi} \approx 1.128379, then a single tangent at the posterior mode suffices. Expected candidates per draw in that dimension: \sqrt{1 + a_i}. Otherwise, a symmetric three-point envelope is used at (\theta^\star_i - \omega_i, \theta^\star_i, \theta^\star_i + \omega_i), with expected candidates per draw bounded above by 2/\sqrt{\pi}.

Gridtype 2: Adaptive Optimization

Each dimension is assigned either a single-point or three-point envelope by minimizing

T_\mathrm{total}(g_i) = T_\mathrm{build}(g_i) + T_\mathrm{sample}(n, acc_i(g_i)).

The optimizer balances build cost (grows with number of tangents) against sampling cost (decreases as acceptance improves). Expected candidates per draw: \prod_j \mathrm{scaleest}_{i,j}, where each factor is either \sqrt{1+a_j} (single-point) or 2/\sqrt{\pi} (three-point), depending on the optimization outcome.

Gridtype 3: Always Three-Point

Every dimension uses a symmetric three-point envelope. Expected candidates per draw:

\left(\tfrac{2}{\sqrt{\pi}}\right)^k

for k dimensions, as shown in Theorem 3 of \insertCiteNygren2006glmbayes.

Gridtype 4: Always Single-Point

Every dimension uses a single tangent at the posterior mode. Expected candidates per draw:

\prod_{i=1}^k \sqrt{1 + a_i}

(Example 1 in \insertCiteNygren2006glmbayes).

References

\insertAllCited

Examples

data(menarche,package="MASS")

summary(menarche)
plot(Menarche/Total ~ Age, data=menarche)

Age2=menarche$Age-13

x<-matrix(as.numeric(1.0),nrow=length(Age2),ncol=2)
x[,2]=Age2

y=menarche$Menarche/menarche$Total
wt=menarche$Total

mu<-matrix(as.numeric(0.0),nrow=2,ncol=1)
mu[2,1]=(log(0.9/0.1)-log(0.5/0.5))/3

V1<-1*diag(as.numeric(2.0))

# 2 standard deviations for prior estimate at age 13 between 0.1 and 0.9
## Specifies uncertainty around the point estimates

V1[1,1]<-((log(0.9/0.1)-log(0.5/0.5))/2)^2 
V1[2,2]=(3*mu[2,1]/2)^2  # Allows slope to be up to 1 times as large as point estimate 

famfunc<-glmbfamfunc(binomial(logit))

f1<-famfunc$f1
f2<-famfunc$f2
f3<-famfunc$f3
f5<-famfunc$f5
f6<-famfunc$f6

dispersion2<-as.numeric(1.0)
start <- mu
offset2=rep(as.numeric(0.0),length(y))
P=solve(V1)
n=1000



###### Adjust weight for dispersion

wt2=wt/dispersion2

######################### Shift mean vector to offset so that adjusted model has 0 mean

alpha=x%*%as.vector(mu)+offset2
mu2=0*as.vector(mu)
P2=P
x2=x


#####  Optimization step to find posterior mode and associated Precision

parin=start-mu

opt_out=optim(parin,f2,f3,y=as.vector(y),x=as.matrix(x),mu=as.vector(mu2),
              P=as.matrix(P),alpha=as.vector(alpha),wt=as.vector(wt2),
              method="BFGS",hessian=TRUE
)

bstar=opt_out$par  ## Posterior mode for adjusted model
bstar
bstar+as.vector(mu)  # mode for actual model
A1=opt_out$hessian # Approximate Precision at mode

## Standardize Model

Standard_Mod=glmb_Standardize_Model(y=as.vector(y), x=as.matrix(x),
P=as.matrix(P),bstar=as.matrix(bstar,ncol=1), A1=as.matrix(A1))

bstar2=Standard_Mod$bstar2
A=Standard_Mod$A
x2=Standard_Mod$x2
mu2=Standard_Mod$mu2
P2=Standard_Mod$P2
L2Inv=Standard_Mod$L2Inv
L3Inv=Standard_Mod$L3Inv

## Derive a and G1 (as EnvelopeBuild does internally)
a <- diag(A)
omega <- (sqrt(2) - exp(-1.20491 - 0.7321*sqrt(0.5 + a))) / sqrt(1 + a)
b2 <- as.vector(bstar2)
G1 <- rbind(b2 - omega, b2, b2 + omega)

## EnvelopeOpt: standalone call (used by EnvelopeSize when Gridtype=2)
grid_opt <- EnvelopeOpt(a, n)
grid_opt

## EnvelopeSize for each Gridtype
size_1 <- EnvelopeSize(a, G1, Gridtype=1L, n=n)   # static threshold
size_2 <- EnvelopeSize(a, G1, Gridtype=2L, n=n)   # uses EnvelopeOpt
size_3 <- EnvelopeSize(a, G1, Gridtype=3L, n=n)   # always 3-point
size_4 <- EnvelopeSize(a, G1, Gridtype=4L, n=n)   # always single-point

## EnvelopeBuild for each Gridtype
Env_1 <- EnvelopeBuild(as.vector(bstar2), as.matrix(A), y, as.matrix(x2),
  as.matrix(mu2,ncol=1), as.matrix(P2), as.vector(alpha), as.vector(wt2),
  family="binomial", link="logit", Gridtype=1L, n=as.integer(n), sortgrid=FALSE)
Env_2 <- EnvelopeBuild(as.vector(bstar2), as.matrix(A), y, as.matrix(x2),
  as.matrix(mu2,ncol=1), as.matrix(P2), as.vector(alpha), as.vector(wt2),
  family="binomial", link="logit", Gridtype=2L, n=as.integer(n), sortgrid=FALSE)
Env_3 <- EnvelopeBuild(as.vector(bstar2), as.matrix(A), y, as.matrix(x2),
  as.matrix(mu2,ncol=1), as.matrix(P2), as.vector(alpha), as.vector(wt2),
  family="binomial", link="logit", Gridtype=3L, n=as.integer(n), sortgrid=FALSE)
Env_4 <- EnvelopeBuild(as.vector(bstar2), as.matrix(A), y, as.matrix(x2),
  as.matrix(mu2,ncol=1), as.matrix(P2), as.vector(alpha), as.vector(wt2),
  family="binomial", link="logit", Gridtype=4L, n=as.integer(n), sortgrid=FALSE)

Env_3

glmbayes documentation built on Aug. 5, 2026, 1:07 a.m.