Package 'clmstan'

Description

Combines multiple prior() objects into a single prior specification list.

Usage

## S3 method for class 'clm_prior_spec'
c(...)

Arguments

Value

An object of class "clm_prior_list" containing all prior specifications.

Examples

# Combine multiple priors
priors <- c(
  prior(normal(0, 2.5), class = "b"),
  prior(normal(0, 10), class = "Intercept"),
  prior(gamma(2, 0.1), class = "df")
)
print(priors)

Description

Creates a Cauchy distribution object for use with prior().

Usage

cauchy(mu = 0, sigma = 1)

Arguments

Value

An object of class "clm_dist" representing a Cauchy distribution.

See Also

prior(), normal(), gamma(), student_t()

Examples

# Create a Cauchy prior (weakly informative)
cauchy(0, 2.5)

# Use with prior()
prior(cauchy(0, 2.5), class = "b")

Description

Create prior specifications for cumulative link models in clmstan.

Default priors:

• Regression coefficients (beta): normal(0, 2.5)

• Cutpoints (c): normal(0, 10) for flexible, normal(0, 5) for symmetric

• Interval (d): gamma(2, 0.5) for equidistant threshold

Link parameter priors (when estimated):

Usage

clm_prior(
  beta_sd = NULL,
  c_sd = NULL,
  c1_mu = NULL,
  c1_sd = NULL,
  d_alpha = NULL,
  d_beta = NULL,
  cpos_sd = NULL,
  df_alpha = NULL,
  df_beta = NULL,
  lambda_ao_alpha = NULL,
  lambda_ao_beta = NULL,
  lambda_lg_mu = NULL,
  lambda_lg_sd = NULL,
  xi_mu = NULL,
  xi_sd = NULL,
  r_alpha = NULL,
  r_beta = NULL,
  theta1_alpha = NULL,
  theta1_beta = NULL,
  theta2_alpha = NULL,
  theta2_beta = NULL
)

Arguments

Value

An object of class "clm_prior" containing prior specifications.

Examples

# Create a prior object (does not require Stan)
my_prior <- clm_prior(beta_sd = 2, c_sd = 5)
print(my_prior)

## Not run: 
# Examples below require CmdStan and compiled Stan models
data(wine, package = "ordinal")

# Default priors (no customization needed)
fit <- clm_stan(rating ~ temp, data = wine,
                chains = 2, iter = 500, warmup = 250, refresh = 0)

# Custom prior for regression coefficients
fit2 <- clm_stan(rating ~ temp, data = wine,
                 prior = clm_prior(beta_sd = 1),
                 chains = 2, iter = 500, warmup = 250, refresh = 0)

## End(Not run)

Description

Fit a Cumulative Link Model using CmdStanR

Usage

clm_stan(
  formula,
  data,
  link = "logit",
  base = "logit",
  threshold = "flexible",
  link_param = NULL,
  prior = NULL,
  chains = 4,
  iter = 2000,
  warmup = NULL,
  ...
)

Arguments

Value

An object of class "clmstan"

Examples

## Not run: 
# Fit a proportional odds model
library(ordinal)
data(wine)
fit <- clm_stan(rating ~ temp + contact, data = wine, link = "logit")
print(fit)

# Fit with t-link (fixed df)
fit_t <- clm_stan(rating ~ temp, data = wine, link = "tlink",
                  link_param = list(df = 8))

# Fit with GEV link (estimate xi)
fit_gev <- clm_stan(rating ~ temp, data = wine, link = "gev",
                    link_param = list(xi = "estimate"))

## End(Not run)

Description

The clmstan class represents a fitted cumulative link model. It contains the CmdStanR fit object and additional metadata.

Slots

fit

The CmdStanMCMC object from cmdstanr

formula

The model formula

data

The original data frame

link

The link function used

base

The base distribution (for SP link)

threshold

The threshold structure

link_param

Link parameter settings (for flexible links)

full

TRUE if link parameters were estimated (Bayesian inference), FALSE if they were fixed at user-specified values

K

Number of response categories (cached from data)

N

Number of observations (extracted from data for efficiency)

P

Number of predictors (extracted from design matrix)

Description

Returns posterior point estimates (mean or median) for all model parameters.

Usage

## S3 method for class 'clmstan'
coef(object, type = c("mean", "median"), ...)

Arguments

Value

A named numeric vector with:

• Threshold coefficients (e.g., "1|2", "2|3", ...)

• Regression coefficients (variable names from formula)

Examples

## Not run: 
fit <- clm_stan(rating ~ temp, data = wine)
coef(fit)
coef(fit, type = "median")

## End(Not run)

Description

Provides a summary of MCMC convergence diagnostics including HMC-specific diagnostics (divergences, treedepth, E-BFMI) and general convergence measures (Rhat, ESS).

Usage

diagnostics(object, ...)

## S3 method for class 'clmstan'
diagnostics(
  object,
  detail = FALSE,
  rhat_threshold = 1.01,
  ess_threshold = 400,
  ...
)

Arguments

Details

The function checks for the following issues:

• Divergences: Number of divergent transitions (ideally 0)

• Treedepth: Transitions hitting max treedepth (efficiency issue)

• E-BFMI: Energy Bayesian Fraction of Missing Information (values < 0.3 indicate problems)

• Rhat: Potential scale reduction factor (values > 1.01 indicate lack of convergence)

• ESS: Effective sample size for bulk and tail (low values indicate high autocorrelation)

Value

Invisibly returns a list containing:

• hmc: HMC diagnostics from CmdStanMCMC$diagnostic_summary()

• convergence: Data frame of per-parameter Rhat and ESS values

• issues: Logical indicating whether any issues were detected

Examples

## Not run: 
fit <- clm_stan(rating ~ temp, data = wine)
diagnostics(fit)
diagnostics(fit, detail = TRUE)

## End(Not run)

Description

Computes autocorrelation function (ACF) values for MCMC chains and returns them in a tidy data frame format.

Usage

extract_acf(object, pars = NULL, lags = 20, ...)

Arguments

Details

The ACF measures how correlated each draw is with previous draws in the same chain. High autocorrelation at many lags indicates slow mixing and the need for more samples or reparameterization.

Ideally, ACF should drop to near zero within a few lags. Persistent high autocorrelation suggests the sampler is exploring the posterior slowly.

Value

A data frame with columns:

• parameter: Parameter name

• chain: Chain number

• lag: Lag value (0, 1, 2, ...)

• acf: Autocorrelation value

Examples

## Not run: 
fit <- clm_stan(rating ~ temp, data = wine)
acf_df <- extract_acf(fit)
head(acf_df)

# Plot ACF for specific parameters
library(ggplot2)
acf_df |>
  dplyr::filter(parameter == "beta[1]") |>
  ggplot(aes(x = lag, y = acf, color = factor(chain))) +
  geom_line() +
  geom_hline(yintercept = 0, linetype = "dashed")

## End(Not run)

Description

Returns expected category probabilities for each observation. This is equivalent to predict(object, type = "probs", summary = TRUE).

Usage

## S3 method for class 'clmstan'
fitted(
  object,
  newdata = NULL,
  summary = TRUE,
  robust = FALSE,
  probs = c(0.025, 0.975),
  ndraws = NULL,
  ...
)

Arguments

Value

If summary = TRUE (default): A data frame with N rows and columns for each category probability (P[Y=1], P[Y=2], etc.). If summary = FALSE: An S x N x K array of probability draws.

See Also

predict.clmstan(), posterior_predict.clmstan()

Description

Creates a flat (improper uniform) distribution object for use with prior(). A flat prior assigns equal probability density to all values, which is improper (does not integrate to 1) but can be used when the likelihood provides sufficient information for identification.

Usage

flat()

Value

An object of class "clm_dist" representing a flat distribution.

Note

Flat priors are supported for:

• Regression coefficients (class "b")

• Threshold classes ("Intercept", "c1", "cpos")

Using flat priors may lead to improper posteriors if the likelihood does not provide sufficient information. For thresholds with ordered constraints, Stan's internal transformation provides implicit regularization.

See Also

prior(), normal(), student_t(), cauchy()

Examples

# Create a flat prior for regression coefficients
prior(flat(), class = "b")

# Flat prior for thresholds (flexible)
prior(flat(), class = "Intercept")

Description

Creates a gamma distribution object for use with prior().

Usage

gamma(alpha, beta)

Arguments

Value

An object of class "clm_dist" representing a gamma distribution.

Note

This function masks base::gamma(). To use the base gamma function, use base::gamma() explicitly.

See Also

prior(), normal(), student_t(), cauchy()

Examples

# Create a gamma prior
gamma(2, 0.1)

# Use with prior() for degrees of freedom
prior(gamma(2, 0.1), class = "df")

Description

Check if object is clmstan

Usage

is.clmstan(x)

Arguments

Value

TRUE if x is a clmstan object

Description

clmstan supports the following link functions for cumulative link models:

Standard links (no additional parameters):

• "logit" - Logistic (proportional odds model)

• "probit" - Normal (latent variable interpretation)

• "cloglog" - Complementary log-log (proportional hazards)

• "loglog" - Log-log (Gumbel minimum)

• "cauchit" - Cauchy (heavy tails)

Flexible links (with additional parameters):

• "tlink" - Student-t (df > 0)

  • df = Inf: equals probit

  • df < 3: increasingly heavy tails; df > 30 is nearly normal

• "aranda_ordaz" - Aranda-Ordaz asymmetric (lambda > 0)

  • lambda = 1: equals logit

  • lambda -> 0: approaches cloglog

• "gev" - Generalized extreme value (shape parameter xi)

  • xi = 0: Gumbel (equals loglog)

  • xi < 0: Weibull (short tail)

  • xi > 0: Frechet (heavy tail)

• "sp" - Symmetric power (r > 0, base distribution)

  • r = 1: equals base distribution

  • 0 < r < 1: positively skewed

  • r > 1: negatively skewed

• "log_gamma" - Log-gamma (lambda)

  • lambda = 0: equals probit

  • lambda > 0 or < 0: asymmetric

• "aep" - Asymmetric exponential power (theta1 > 0, theta2 > 0)

  • alpha = 0.5 fixed for identifiability

  • theta1 = theta2: symmetric distribution

  • theta = 2: Gaussian kernel (but NOT equal to probit due to scaling)

  • theta < 2: heavy tails (leptokurtic)

  • theta > 2: light tails (platykurtic)

Link Parameter Specification

Flexible link parameters can be either fixed or estimated (inferred).

Fixed parameters: Specify a numeric value

clm_stan(y ~ x, link = "tlink", link_param = list(df = 8))
clm_stan(y ~ x, link = "gev", link_param = list(xi = 0))  # equals loglog
clm_stan(y ~ x, link = "aep", link_param = list(theta1 = 2, theta2 = 2))  # symmetric

Estimated parameters: Use "estimate" (with default prior)

clm_stan(y ~ x, link = "tlink", link_param = list(df = "estimate"))
clm_stan(y ~ x, link = "gev", link_param = list(xi = "estimate"))

Custom priors: Combine "estimate" with prior argument

clm_stan(y ~ x, link = "gev",
         link_param = list(xi = "estimate"),
         prior = prior(normal(0, 0.3), class = "xi"))

Default Priors for Link Parameters

When using "estimate", the following default priors are used:

SP Link Details (Li et al., 2019)

The Symmetric Power link uses a symmetric base distribution F_0, specified via the base argument. Supported bases:

• base = "logit": Logistic base (default)

• base = "probit": Normal base

• base = "cauchit": Cauchy base

• base = "tlink": Student-t base (requires df)

Note: Li et al. (2019) define F_0 as a CDF "whose corresponding PDF is symmetric about 0".

Description

Computes approximate leave-one-out cross-validation (LOO-CV) for a fitted cumulative link model using Pareto smoothed importance sampling (PSIS).

Usage

## S3 method for class 'clmstan'
loo(x, ..., r_eff = NULL, cores = getOption("mc.cores", 1), save_psis = FALSE)

Arguments

Details

The function extracts the log-likelihood matrix (log_lik) computed in the generated quantities block of the Stan model and passes it to loo.

Pareto k diagnostics: Observations with high Pareto k values (k > 0.7) indicate potential problems with the LOO approximation for those observations. Use plot() on the returned object to visualize the Pareto k values.

Model comparison: Use loo_compare to compare multiple models. Models with higher elpd_loo are preferred.

Value

An object of class c("psis_loo", "loo") containing:

• estimates: A matrix with columns Estimate and SE for elpd_loo, p_loo, and looic.

• pointwise: A matrix with pointwise contributions.

• diagnostics: A list with Pareto k values and effective sample sizes for each observation.

See Also

waic.clmstan for WAIC computation, loo for details on the LOO algorithm, loo_compare for model comparison.

Examples

## Not run: 
fit <- clm_stan(rating ~ temp, data = wine)
loo_result <- loo(fit)
print(loo_result)
plot(loo_result)

# Compare two models
fit1 <- clm_stan(rating ~ temp, data = wine, link = "logit")
fit2 <- clm_stan(rating ~ temp, data = wine, link = "probit")
loo::loo_compare(loo(fit1), loo(fit2))

## End(Not run)

Description

Creates a normal distribution object for use with prior().

Usage

normal(mu = 0, sigma = 1)

Arguments

Value

An object of class "clm_dist" representing a normal distribution.

See Also

prior(), gamma(), student_t(), cauchy()

Examples

# Create a normal prior
normal(0, 2.5)

# Use with prior()
prior(normal(0, 2.5), class = "b")

Description

Produces diagnostic plots using the bayesplot package.

Usage

## S3 method for class 'clmstan'
plot(
  x,
  type = c("trace", "dens", "hist", "areas", "intervals", "acf"),
  pars = NULL,
  ...
)

Arguments

Value

A ggplot object

Examples

## Not run: 
fit <- clm_stan(rating ~ temp, data = wine)
plot(fit)                    # trace plots
plot(fit, type = "dens")     # density plots
plot(fit, type = "intervals") # credible intervals
plot(fit, type = "acf")      # autocorrelation plots
plot(fit, pars = "beta")     # only beta parameters

## End(Not run)

Description

Draws from the posterior predictive distribution. For each posterior sample, a predicted category is sampled from the categorical distribution with the predicted probabilities.

Usage

posterior_predict.clmstan(object, newdata = NULL, ndraws = NULL, ...)

Arguments

Value

An integer matrix of dimension S x N containing predicted categories (1 to K), where S is the number of posterior draws and N is the number of observations.

See Also

predict.clmstan(), fitted.clmstan()

Description

Generates predictions from a fitted cumulative link model.

Usage

## S3 method for class 'clmstan'
predict(
  object,
  newdata = NULL,
  type = c("class", "probs"),
  summary = TRUE,
  robust = FALSE,
  probs = c(0.025, 0.975),
  ndraws = NULL,
  ...
)

Arguments

Value

Depending on type and summary:

• type = "class", summary = TRUE: A data frame with columns Estimate (mean/median predicted class), Est.Error (SD), quantile columns, and Class (modal predicted category).

• type = "class", summary = FALSE: An S x N integer matrix of predicted categories (1 to K), where S is the number of posterior draws and N is the number of observations.

• type = "probs", summary = TRUE: A data frame with columns for each category probability (P[Y=1], P[Y=2], etc.).

• type = "probs", summary = FALSE: An S x N x K array of predicted probabilities.

See Also

fitted.clmstan() for expected probabilities, posterior_predict.clmstan() for posterior predictive samples.

Description

Print method for clm_dist objects

Usage

## S3 method for class 'clm_dist'
print(x, ...)

Arguments

Value

Invisibly returns the input clm_dist object.

Description

Print method for clm_prior objects

Usage

## S3 method for class 'clm_prior'
print(x, ...)

Arguments

Value

Invisibly returns the input clm_prior object.

Description

Print method for clm_prior_list objects

Usage

## S3 method for class 'clm_prior_list'
print(x, ...)

Arguments

Value

Invisibly returns the input clm_prior_list object.

Description

Print method for clm_prior_spec objects

Usage

## S3 method for class 'clm_prior_spec'
print(x, ...)

Arguments

Value

Invisibly returns the input clm_prior_spec object.

Description

Print method for clmstan objects

Usage

## S3 method for class 'clmstan'
print(x, ...)

Arguments

Value

Invisibly returns x

Description

Print method for summary.clmstan objects

Usage

## S3 method for class 'summary.clmstan'
print(x, ...)

Arguments

Value

Invisibly returns x

Description

Specify prior distributions for model parameters using distribution functions.

Usage

prior(prior, class = "b", coef = "")

Arguments

Value

An object of class "clm_prior_spec" representing the prior specification.

See Also

normal(), gamma(), student_t(), cauchy(), clm_prior()

Examples

# Specify a normal prior for regression coefficients
prior(normal(0, 2.5), class = "b")

# Specify a gamma prior for degrees of freedom
prior(gamma(2, 0.1), class = "df")

# Combine multiple priors
c(
  prior(normal(0, 2.5), class = "b"),
  prior(normal(0, 10), class = "Intercept")
)

Description

Creates a Student-t distribution object for use with prior().

Usage

student_t(df = 3, mu = 0, sigma = 1)

Arguments

Value

An object of class "clm_dist" representing a Student-t distribution.

See Also

prior(), normal(), gamma(), cauchy()

Examples

# Create a Student-t prior with heavy tails
student_t(3, 0, 2.5)

# Use with prior()
prior(student_t(3, 0, 2.5), class = "b")

Description

Summary method for clmstan objects

Usage

## S3 method for class 'clmstan'
summary(object, probs = c(0.025, 0.5, 0.975), digits = 3, ...)

Arguments

Value

An object of class "summary.clmstan" containing:

• coefficients: Posterior summary for regression coefficients

• thresholds: Posterior summary for threshold parameters

• beta0: Posterior summary for intercept

• link_params: Posterior summary for link parameters (full model only)

• model_info: Model metadata (formula, link, threshold, K, N, P)

Description

Get supported link functions

Usage

supported_links(type = c("all", "standard", "flexible"))

Arguments

Value

A character vector of supported link function names

Examples

supported_links()
supported_links("standard")
supported_links("flexible")

Description

Get supported threshold structures

Usage

supported_thresholds()

Value

A character vector of supported threshold structure names

Examples

supported_thresholds()

Description

Computes the Widely Applicable Information Criterion (WAIC) for a fitted cumulative link model.

Usage

## S3 method for class 'clmstan'
waic(x, ...)

Arguments

Details

WAIC is an alternative to LOO-CV that is asymptotically equivalent to leave-one-out cross-validation. However, LOO-CV with PSIS is generally preferred because:

• It provides useful diagnostics (Pareto k values)

• It is more robust in finite samples

• It has been shown to be more reliable in practice

For most purposes, loo.clmstan is recommended over WAIC.

Value

An object of class c("waic", "loo") containing:

• estimates: A matrix with columns Estimate and SE for elpd_waic, p_waic, and waic.

• pointwise: A matrix with pointwise contributions.

See Also

loo.clmstan for LOO-CV (recommended), waic for details on WAIC computation, loo_compare for model comparison.

Examples

## Not run: 
fit <- clm_stan(rating ~ temp, data = wine)
waic_result <- waic(fit)
print(waic_result)

## End(Not run)