Type: Package
Title: Additive Model for Ordinal Data using Laplace P-Splines
Version: 0.9.1
Maintainer: Philippe Lambert <p.lambert@uliege.be>
BugReports: https://github.com/plambertULiege/ordgam/issues
Description: Additive proportional odds model for ordinal data using Laplace P-splines. The combination of Laplace approximations and P-splines enable fast and flexible inference in a Bayesian framework. Specific approximations are proposed to account for the asymmetry in the marginal posterior distributions of non-penalized parameters. For more details, see Lambert and Gressani (2023) <doi:10.1177/1471082X231181173> ; Preprint: <doi:10.48550/arXiv.2210.01668>).
License: GPL-3
URL: <https://github.com/plambertULiege/ordgam>
Encoding: UTF-8
LazyData: true
Imports: cubicBsplines, Matrix, mgcv, marqLevAlg, sn, MASS, numDeriv, ucminf
RoxygenNote: 7.2.1
Suggests: testthat (≥ 3.0.0)
Config/testthat/edition: 3
NeedsCompilation: no
Packaged: 2023-09-13 11:26:17 UTC; plambert
Author: Philippe Lambert ORCID iD [aut, cre]
Depends: R (≥ 3.5.0)
Repository: CRAN
Date/Publication: 2023-09-14 09:20:02 UTC

Internal function extracting design matrices from formulas in the DALSM function and computing penalty related matrices.

Description

Internal function extracting design matrices from formulas in the DALSM function and computing penalty related matrices.

Usage

DesignFormula(
  formula,
  data,
  K = 10,
  pen.order = 2,
  knots.x = NULL,
  n = NULL,
  nointercept = FALSE,
  verbose = TRUE
)

Compute the penalty matrix associated to a vector containing fixed (non-penalized) parameters and equal-size sub-vectors of penalized spline parameters

Description

Compute the penalty matrix associated to a vector containing fixed (non-penalized) parameters and equal-size sub-vectors of penalized spline parameters

Usage

Pcal.fun(nfixed, lambda, Pd.x)

Arguments

nfixed

the number of fixed (i.e. non-penalized) parameters

lambda

a vector of J penalty parameters where each component is associated to a sub-vector of spline parameters of length K

Pd.x

a penalty matrix of size J associated to a given sub-vector of spline parameters

Value

A block diagonal penalty matrix of size (nfixed+JK) given by Blockdiag(diag(0,nfixed), diag(lambda).kron.Pd.x)

References

Lambert, P. and Gressani, 0. (2023) Penalty parameter selection and asymmetry corrections to Laplace approximations in Bayesian P-splines models. Statistical Modelling. <doi:10.1177/1471082X231181173>. Preprint: <arXiv:2210.01668>.

Examples

Dd = diff(diag(1,5),diff=2) ## Difference penalty matrix for a vector of length 5
Pd = t(Dd) %*% Dd ## Penalty matrix of order 2
nfixed = 2 ## 2 unpenalized parameters
## Global penalty matrix when 2 unpenalized parameters and 2 additive terms with
##   2 vectors of 5 P-splines coefficients with lambda values 10 and 100 respectively.
Pcal.fun(nfixed=2,lambda=c(10,100),Pd)

Skew-Normal approximation to a density evaluated on a sparse grid

Description

Skew-Normal approximation to a density evaluated on a sparse grid

Usage

SNapprox(x, lfx)

Arguments

x

Vector containing a grid of values on the density support and covering the posterior mode.

lfx

Log density values on the grid x (possibly up to an additive constant)

Value

A list containing

Author(s)

Philippe Lambert p.lambert@uliege.be

References

Lambert, P. and Gressani, 0. (2023) Penalty parameter selection and asymmetry corrections to Laplace approximations in Bayesian P-splines models. Statistical Modelling. <doi:10.1177/1471082X231181173>. Preprint: <arXiv:2210.01668>.

See Also

STapprox.

Examples

library(ordgam)

## Density to be approximated by a Skew-Normal
dtarget = function(x) dgamma(x,10,2)
curve(dtarget(x),0,15,lwd=2,ylab="Density")

## Values of the target density on a sparse grid
ngrid = 6 ## Sparse grid size
xgrid = seq(2,8,length=ngrid) ## Grid
lfx = log(dtarget(xgrid)) ## Log values

## Skew-Normal approximation
dp = ordgam::SNapprox(xgrid,lfx)$dp
curve(sn::dsn(x,dp=dp),add=TRUE,lwd=2,lty=2,col=2)
points(xgrid,exp(lfx))
legend("topright",legend=c("Target density","Skew-Normal approx."),
       col=1:2,lty=1:2,lwd=2,bty="n")#' library(ordgam)

Skew-t approximation to a density evaluated on a sparse grid

Description

Skew-t approximation to a density evaluated on a sparse grid

Usage

STapprox(x, lfx)

Arguments

x

Vector containing a grid of values on the density support and covering the posterior mode.

lfx

Log density values on the grid x (possibly up to an additive constant)

Value

A list containing

Author(s)

Philippe Lambert p.lambert@uliege.be

References

Lambert, P. and Gressani, 0. (2023) Penalty parameter selection and asymmetry corrections to Laplace approximations in Bayesian P-splines models. Statistical Modelling. <doi:10.1177/1471082X231181173>. Preprint: <arXiv:2210.01668>.

See Also

SNapprox.

Examples

library(ordgam)

## Density to be approximated by a Skew-t
dtarget = function(x) dgamma(x,10,2)
curve(dtarget(x),0,15,lwd=2,ylab="Density")

## Values of the target density on a sparse grid
ngrid = 6 ## Sparse grid size
xgrid = seq(2,8,length=ngrid) ## Grid
lfx = log(dtarget(xgrid)) ## Log values

## Skew-t approximation
dp = ordgam::STapprox(xgrid,lfx)$dp
curve(sn::dst(x,dp=dp),add=TRUE,lwd=2,lty=2,col=2)
points(xgrid,exp(lfx))
legend("topright",legend=c("Target density","Skew-t approx."),
       col=1:2,lty=1:2,lwd=2,bty="n")

Generation of a recentered B-spline basis matrix in additive models

Description

Generation of a B-spline basis matrix with recentered columns to handle the identifiability constraint in additive models. See Wood (CRC Press 2017, pp. 175-176) for more details.

Usage

centeredBasis.gen(x, knots, cm = NULL, pen.order = 2, verbose = TRUE)

Arguments

x

vector of values where to compute the "recentered" B-spline basis

knots

vector of knots (that should cover the values in <x>)

cm

(Optional) values subtracted from each column of the original B-spline matrix

pen.order

penalty order for the B-spline parameters (Default: 2)

verbose

verbose indicator (Default: TRUE)

Value

List containing

Author(s)

Philippe Lambert p.lambert@uliege.be

References

Lambert, P. and Gressani, 0. (2023) Penalty parameter selection and asymmetry corrections to Laplace approximations in Bayesian P-splines models. Statistical Modelling. <doi:10.1177/1471082X231181173>. Preprint: <arXiv:2210.01668>.

Examples

x = seq(0,1,by=.01)
knots = seq(0,1,length=5)
obj = centeredBasis.gen(x,knots)
matplot(x,obj$B,type="l",ylab="Centered B-splines")
colMeans(obj$B)

Perception of gay men and lesbians in Wallonia, Belgium.

Description

Data extracted from the European Social Survey (2018) conducted in a series of European countries including Belgium and one of its three regions, Wallonia. Each of the 552 participants (aged at least 15) was asked to react to the following statement, "Gay men and lesbians should be free to live their own life as they wish", with a positioning on a Likert scale going from 1 (=Agree strongly) to 5 (=Disagree strongly), with 3 labelled as Neither agree nor disagree.

Usage

data(freehmsData)

Format

A data frame with 552 rows and 4 columns.

freehms

Ordinal response (1: Agree strongly to 5: Disagree strongly).

age

Age of the respondent.

gndr

Gender of the respondent.

eduyrs

Number of years of education completed.

References

Lambert, P. and Gressani, 0. (2023) Penalty parameter selection and asymmetry corrections to Laplace approximations in Bayesian P-splines models. Statistical Modelling. <doi:10.1177/1471082X231181173>. Preprint: <arXiv:2210.01668>.

European Social Survey Round 9 Data (2018) Data file edition 3.1. NSD - Norwegian Centre for Research Data, Norway.

Perception of gay men and lesbians in Belgium.

Description

Data extracted from the European Social Survey (2018) conducted in a series of European countries including Belgium. Each of the 1737 participants (aged at least 15) was asked to react to the following statement, "Gay men and lesbians should be free to live their own life as they wish", with a positioning on a Likert scale going from 1 (=Agree strongly) to 5 (=Disagree strongly), with 3 labelled as Neither agree nor disagree.

Usage

data(freehmsDataBE)

Format

A data frame with 1737 rows and 5 columns.

freehms

Ordinal response (1: Agree strongly to 5: Disagree strongly).

gndr

Gender of the respondent.

age

Age of the respondent.

eduyrs

Number of years of education completed.

region

Belgian region of residence: WAL (Wallonia), FL (Flanders) or BXL (Brussels).

References

Lambert, P. and Gressani, 0. (2023) Penalty parameter selection and asymmetry corrections to Laplace approximations in Bayesian P-splines models. Statistical Modelling. <doi:10.1177/1471082X231181173>. Preprint: <arXiv:2210.01668>.

European Social Survey Round 9 Data (2018) Data file edition 3.1. NSD - Norwegian Centre for Research Data, Norway.

Marginal posterior density function for a remapped non-penalized parameter in an ordgam model

Description

Marginal posterior density function for a remapped non-penalized parameter in an ordgam model

Usage

lmarg.gammaTilde(gamtk, k, model)

Arguments

gamtk

Remapped parameter value at which the marginal log posterior density for <gamma.tilde[k]> must be evaluated.

k

Targetted component in the vector of remapped non-penalized parameters <gamma.tilde>.

model

An ordgam.object.

Value

Log of p(gamma.tilde[k] | lambda,data)

Author(s)

Philippe Lambert p.lambert@uliege.be

References

Lambert, P. and Gressani, 0. (2023) Penalty parameter selection and asymmetry corrections to Laplace approximations in Bayesian P-splines models. Statistical Modelling. <doi:10.1177/1471082X231181173>. Preprint: <arXiv:2210.01668>.

See Also

ordgam.

Examples

library(ordgam)
data(freehmsData)
mod = ordgam(freehms ~ s(eduyrs) + s(age), data=freehmsData, descending=TRUE,
             lambda0=c(192,18),select.lambda=FALSE)
ngamma = with(mod, nalpha+nfixed) ## Number of non-penalized parms
k = 1 ## Focus on gamma.tilde[1]
x.grid = seq(-4,4,length=7) ## Grid of values for gamma.tilde[k]
lfy.grid = ordgam::lmarg.gammaTilde(x.grid,k=k,mod) ## log p(gamma.tilde[k] | D) on the grid
gamt.ST = ordgam::STapprox(x.grid,lfy.grid)$dp ## Approximate using a skew-t
## Plot the estimated marginal posterior for <gamma.tilde[k]>
xlab = bquote(tilde(gamma)[.(k)])
ylab = bquote(p(tilde(gamma)[.(k)]~ "|"~lambda~","~D))
xlim = sn::qst(c(.0001,.9999),dp=gamt.ST)
curve(sn::dst(x,dp=gamt.ST),xlim=xlim,
      xlab=xlab,ylab=ylab,col="blue",lwd=2,lty=1)

Posterior density function for the non-penalized parameters in an ordgam model

Description

Posterior density function for the non-penalized parameters in an ordgam model

Usage

lpost.gamma(model)

Arguments

model

An ordgam.object

Value

Log joint posterior density function for the non-penalized regression parameters.

Author(s)

Philippe Lambert p.lambert@uliege.be

References

Lambert, P. and Gressani, 0. (2023) Penalty parameter selection and asymmetry corrections to Laplace approximations in Bayesian P-splines models. Statistical Modelling. <doi:10.1177/1471082X231181173>. Preprint: <arXiv:2210.01668>.

See Also

ordgam.

Examples

library(ordgam)
data(freehmsData)
mod = ordgam(freehms ~ s(eduyrs) + s(age), data=freehmsData, descending=TRUE)
print(mod$theta) ## Model regression parameters
gam.hat = mod$theta[1:4] ## Non-penalized parameter estimates
ordgam::lpost.gamma(mod)(gam.hat)

Fit of an additive proportional odds model for ordinal data using Laplace approximations and P-splines

Description

Fit of an additive proportional odds model for ordinal data using Laplace approximations and P-splines

Usage

ordgam(
  formula,
  data,
  nc = NULL,
  K = 10,
  pen.order = 2,
  descending = TRUE,
  select.lambda = TRUE,
  lambda.family = "dgamma",
  lambda.optimizer = "nlminb",
  lprior.lambda = function(x) dgamma(x, 1, 1e-04, log = TRUE),
  theta0 = NULL,
  lambda0 = NULL,
  ci.level = 0.95,
  verbose = FALSE
)

Arguments

formula

A model formula

data

A data frame containing a column 'y' with the ordinal response (taking integer values) besides the covariates.

nc

(optional) Number of categories for the ordinal response.

K

Number of B-splines to model each additive term (Default: 10).

pen.order

Penalty order (Default: 2).

descending

Logical indicating if the odds of the response taking a value in the upper scale should be preferred over values in the lower scale (Default: TRUE).

select.lambda

Logical indicating if the penalty parameters should be tuned (Default: TRUE).

lambda.family

Prior for <lambda>. Possible choices are "none", "dgamma", "BetaPrime" or "myprior" for a user specified function for the prior of <lambda>.

lambda.optimizer

Algorithm used to maximize p(lambda|data). Possible choices are "nlminb","ucminf","nlm","LevMarq" (Default: "nlminb").

lprior.lambda

Log of the prior density for a <lambda> component if lambda.family set to "myprior".

theta0

(Optional) Vector containing starting values for the regression parameters.

lambda0

Vector of penalty parameters for the additive terms (Default: 10 for each additive term).

ci.level

Confidence levels of the computed credible intervals for the regression parameters.

verbose

Verbose mode (logical)

Value

an object of type ordgam.object.

Author(s)

Philippe Lambert p.lambert@uliege.be

References

Lambert, P. and Gressani, 0. (2023) Penalty parameter selection and asymmetry corrections to Laplace approximations in Bayesian P-splines models. Statistical Modelling. <doi:10.1177/1471082X231181173>. Preprint: <arXiv:2210.01668>.

See Also

ordregr, ordgam.object.

Examples

library(ordgam)
data(freehmsData)
mod = ordgam(freehms ~ gndr + s(eduyrs) + s(age),
             data=freehmsData, descending=TRUE)
print(mod)
plot(mod)

Object resulting from the fit of an additive proportional odds model using 'ordgam'

Description

An object returned by the ordgam function: this is a list with various components related to the fit of such a model.

Value

An ordgam object is a list with following elements:

Author(s)

Philippe Lambert p.lambert@uliege.be

References

Lambert, P. and Gressani, 0. (2023) Penalty parameter selection and asymmetry corrections to Laplace approximations in Bayesian P-splines models. Statistical Modelling. <doi:10.1177/1471082X231181173>. Preprint: <arXiv:2210.01668>.

See Also

ordgam, print.ordregr, plot.ordgam

Compute the additive terms estimated using an 'ordgam' model

Description

Compute the additive terms estimated using an 'ordgam' model

Usage

ordgam_additive(obj.ordgam, ngrid = 300, ci.level = 0.95)

Arguments

obj.ordgam

An object of class 'ordgam'.

ngrid

Number of grid points where the additive terms are computed.

ci.level

Credibility level for the pointwise credible region for the additive terms

Value

a list containing:

Author(s)

Philippe Lambert p.lambert@uliege.be

References

Lambert, P. and Gressani, 0. (2023) Penalty parameter selection and asymmetry corrections to Laplace approximations in Bayesian P-splines models. Statistical Modelling. <doi:10.1177/1471082X231181173>. Preprint: <arXiv:2210.01668>.

Examples

library(ordgam)
data(freehmsData)
mod = ordgam(freehms ~ gndr + s(eduyrs) + s(age),
             data=freehmsData, descending=TRUE)
obj = ordgam_additive(mod)
names(obj)
with(obj$f.grid$age,
      matplot(x, y.mat, lty=c(1,2,2),type="l",col=1,
              xlab="Age", ylab="f(Age)"))

Fit a proportional odds model for ordinal data

Description

Fit a proportional odds model for ordinal data

Usage

ordregr(
  y,
  nc = NULL,
  Xcal = Xcal,
  descending = FALSE,
  prior = list(mean = NULL, Prec = NULL),
  theta0 = NULL,
  ci.level = 0.95
)

Arguments

y

Vector containing the ordinal response (coded using integers in 1:nc).

nc

(optional) Maximum value of y.

Xcal

Design matrix (excluding intercept columns).

descending

Logical indicating if the odds of the response taking a value in the upper scale should be preferred over values in the lower scale.

prior

(optional) List giving the 'mean' and 'Prec'(ision) of the regression parameters.

theta0

(Optional) Vector containing starting values for the regression parameters.

ci.level

Confidence levels of the computed credible intervals for the regression parameters.

Value

An object of class ordregr.object.

Author(s)

Philippe Lambert p.lambert@uliege.be

References

Lambert, P. and Gressani, 0. (2023) Penalty parameter selection and asymmetry corrections to Laplace approximations in Bayesian P-splines models. Statistical Modelling. <doi:10.1177/1471082X231181173>. Preprint: <arXiv:2210.01668>.

See Also

ordgam, ordregr.object.

Examples

library(ordgam)
data(freehmsData)
Xcal = with(freehmsData, cbind(gndr,eduyrs,age))
mod = ordregr(y=freehmsData$freehms, Xcal=Xcal, descending=TRUE)
print(mod)

Object resulting from the fit of a proportional odds model using 'ordregr'

Description

An object returned by the ordregr function: this is a list with various components related to the fit of such a model.

Value

A ordregr object is a list with following elements:

Author(s)

Philippe Lambert p.lambert@uliege.be

References

Lambert, P. and Gressani, 0. (2023) Penalty parameter selection and asymmetry corrections to Laplace approximations in Bayesian P-splines models. Statistical Modelling. <doi:10.1177/1471082X231181173>. Preprint: <arXiv:2210.01668>.

See Also

ordregr, print.ordregr

Log-posterior function for a proportional odds model

Description

Log-posterior function for a proportional odds model

Usage

ordregr_lpost(
  y,
  nc,
  Xcal,
  theta,
  descending = FALSE,
  prior = list(mean = NULL, Prec = NULL),
  gradient = TRUE,
  Hessian = TRUE
)

Arguments

y

Vector containing the ordinal response (coded using integers in 1:nc).

nc

(optional) Maximum value of y.

Xcal

Design matrix.

theta

Vector c(alpha,beta) with intercepts <alpha> and regression parameters <beta>.

descending

Logical indicating if the odds of the response taking a value in the upper scale should be preferred over values in the lower scale.

prior

(optional) List given the mean and Prec(ision) of the regression parameters.

gradient

Logical indicating if the gradient of the log-posterior should be computed.

Hessian

Logical indicating if the Hessian of the log-posterior should be computed.

Value

The log-posterior with the following attributes:

References

Lambert, P. and Gressani, 0. (2023) Penalty parameter selection and asymmetry corrections to Laplace approximations in Bayesian P-splines models. Statistical Modelling. <doi:10.1177/1471082X231181173>. Preprint: <arXiv:2210.01668>.

Plot the the additive terms in an <ordgam> object with its credible regions

Description

Plot the the additive terms in an object generated by ordgam with its credible regions

Usage

## S3 method for class 'ordgam'
plot(x,ngrid=300,ci.level=.95,mfrow=NULL,...)

Arguments

x

An ordgam.object generated by ordgam.

ngrid

An integer indicating the number of gridpoints where the additive terms should be evaluated.

ci.level

Credibility level for the pointwise credible region of the additive terms.

mfrow

(Optional) A vector of the form c(nr, nc). Subsequent figures will be drawn in an nr-by-nc array on the device by rows.

...

Additional generic plotting arguments.

Value

In addition to the plots, an invisible object containing information on the estimated additive terms is returned, see the ordgam_additive function documentation for more details.

Author(s)

Philippe Lambert p.lambert@uliege.be

References

Lambert, P. and Gressani, 0. (2023) Penalty parameter selection and asymmetry corrections to Laplace approximations in Bayesian P-splines models. Statistical Modelling. <doi:10.1177/1471082X231181173>. Preprint: <arXiv:2210.01668>.

See Also

ordgam, ordgam.object, ordgam.additive, print.ordregr.

Examples

library(ordgam)
data(freehmsData)
mod = ordgam(freehms ~ gndr + s(eduyrs) + s(age),
             data=freehmsData, descending=TRUE)
print(mod)
plot(mod)

Print an 'ordregr' or an 'ordgam' object.

Description

Print a summary of the information contained in an ordregr.object or ordgam.object generated by ordregr or ordgam.

Usage

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

Arguments

x

An ordregr.object generated by ordregr or ordgam.

expEst

Logical indicating if the exponential of the regression coefficients should be printed (Default: TRUE)

...

Possible additional printing options.

Value

Print summary statistics.

Author(s)

Philippe Lambert p.lambert@uliege.be

References

Lambert, P. and Gressani, 0. (2023) Penalty parameter selection and asymmetry corrections to Laplace approximations in Bayesian P-splines models. Statistical Modelling. <doi:10.1177/1471082X231181173>. Preprint: <arXiv:2210.01668>.

See Also

ordregr, ordgam

Examples

library(ordgam)
data(freehmsData)
mod = ordgam(freehms ~ gndr + s(eduyrs) + s(age),
             data=freehmsData, descending=TRUE)
print(mod)
plot(mod)

Specification of smooth terms in formulas in the ordgam function.

Description

Specification of smooth terms in formulas in the ordgam function.

Usage

s(x)

Arguments

x

Name of the variable for which an additive term is requested

Value

The submitted variable for which an additive term is required

Significance test of an additive term

Description

Significance test of an additive term relying on the methodology in Wood (Biometrika 2013). It is extracted from a hidden function in the 'mgcv' package. The additive term is estimated using the product of a matrix <X> and a vector <p>.

Usage

testStat(p, X, V, rank = NULL, type = 0, res.df = -1)

Arguments

p

Vector of spline parameters used to estimate f(x)

X

(Centered) B-spline basis evaluated on a fine regular grid on the support of variable <x>

V

Posterior variance-covariance matrix of parameter <p>

rank

Effective dimension of <p>

type

0 value by default

res.df

-1 indicates that the scale is fixed (cf. ordinal response)

Value

Returns a list with following elements: