Type: Package
Title: Copula Based Stochastic Frontier Quantile Model
Version: 0.1.0
Maintainer: Woraphon Yamaka <woraphon.econ@gmail.com>
Description: Provides estimation procedures for copula-based stochastic frontier quantile models for cross-sectional data. The package implements maximum likelihood estimation of quantile regression models allowing flexible dependence structures between error components through various copula families (e.g., Gaussian and Student-t). It enables estimation of conditional quantile effects, dependence parameters, log-likelihood values, and information criteria (AIC and BIC). The framework combines quantile regression methodology introduced by Koenker and Bassett (1978) <doi:10.2307/1913643> with copula theory described in Joe (2014, ISBN:9781466583221). This approach allows modeling heterogeneous effects across quantiles while capturing nonlinear dependence structures between variables.
License: GPL-3
Encoding: UTF-8
Imports: ald, VineCopula, stats, graphics, MASS
RoxygenNote: 7.3.3
NeedsCompilation: no
Packaged: 2026-02-27 07:39:33 UTC; Acer
Author: Woraphon Yamaka [aut, cre], Paravee Maneejuk [aut], Nuttaphong Kaewtathip [aut]
Repository: CRAN
Date/Publication: 2026-03-04 09:50:08 UTC

Technical Efficiency Measure

Description

Computes technical efficiency for a copula-based stochastic frontier quantile model using simulation-based conditional expectations.

Usage

TE(theta, Y, X, family, tau,
   nSim = 200,
   rho2 = NULL,
   seed = NULL,
   plot = FALSE)

Arguments

theta

Numeric vector of estimated model parameters. The expected ordering is: regression coefficients (including intercept), \sigma_v, \sigma_u, copula parameter rho, and optionally rho2 if a two-parameter copula is used.

Y

Numeric vector of dependent variable observations.

X

Numeric matrix (or object coercible to a matrix) of independent variables (without intercept column).

family

Integer specifying the copula family (see VineCopula::BiCopPDF). Some families require a second parameter rho2.

tau

Quantile level in (0,1) for the asymmetric Laplace distribution.

nSim

Number of Monte Carlo draws used to approximate the conditional expectation. Larger values reduce simulation noise but increase computation time.

rho2

Optional second copula parameter. If NULL and theta contains an additional element beyond rho, the function will use that element as rho2.

seed

Optional integer seed for reproducibility in simulation-based computation. If NULL, no seed is set.

plot

Logical. If TRUE, produces a plot of sorted technical efficiency values.

Details

Technical efficiency is computed as a simulated conditional expectation:

TE_i = E[\exp(-U_i) \mid w_i],

where w_i = Y_i - x_i^\top\beta and the weights are constructed using the asymmetric Laplace density for the noise term and the copula density capturing dependence between the inefficiency component and the noise component.

The inefficiency term U is generated from a truncated asymmetric Laplace distribution on (0,\infty).

If plot = TRUE, the function produces a line plot of sorted technical efficiency values.

The following copula families are supported, together with their parameter bounds: 

1  = Gaussian copula (par: (LB = -0.99, UB = 0.99))
2  = Student t copula (par: (LB = -0.99, UB = 0.99); par2: (LB = 0, UB = Inf))
3  = Clayton copula (par: (LB = 0.1, UB = Inf))
4  = Gumbel copula (par: [LB = 0.99, UB = Inf))
5  = Frank copula (par: (LB = -Inf, UB =  0) U (LB = 0, UB = Inf))
6  = Joe copula (par: (LB = 0.99, UB = Inf))
7  = BB1 (Clayton-Gumbel) copula (par: (LB = 0, UB = Inf); par2: [LB = 0.99, UB = Inf))
8  = BB6 copula (par: (LB = 0.99, UB = Inf); par2: (LB = 0, UB = Inf))
9  = BB7 (Joe-Clayton) copula (par: [LB = 0.99, UB = Inf); par2: (LB = 0, UB = Inf))
10 = BB8 copula (par: (LB = 0, UB = 1); par2: (LB = 0, UB =  Inf))
13 = Survival Clayton (180 degrees rotation; par: (LB = 0, UB = Inf))
14 = Survival Gumbel (180 degrees rotation; par: [LB = 0.99, UB = Inf))
16 = Survival Joe (180 degrees rotation; par: (LB = 0.99, UB = Inf))
17 = Survival BB1 (par: (LB = 0, UB = Inf); par2: [LB = 0.99, UB = Inf))
18 = Survival BB6 (par: (LB = 0.99, UB = Inf); par2: (0, UB = Inf))
19 = Survival BB7 (par: [LB = 0.99, UB = Inf); par2: (0, UB = Inf))
20 = Survival BB8 (par: (LB = 0, 0.99); par2: (LB = 0, UB = Inf))
23 = Rotated Clayton (90 degrees; par: (LB = -Inf, UB = 0))
24 = Rotated Gumbel (90 degrees; par: (LB = -Inf, UB = -0.99])
26 = Rotated Joe (90 degrees; par: (LB = -Inf, UB = -0.99))
27 = Rotated BB1 (90 degrees; par: (LB = -Inf, UB = 0); par2: [LB = 0.99, UB = Inf))
28 = Rotated BB6 (90 degrees; par: (LB = -Inf, UB = -0.99]; par2: (LB = 0, UB = Inf))
29 = Rotated BB7 (90 degrees; par: (LB = -Inf, UB = -0.99]; par2: (LB = 0, UB = Inf))
30 = Rotated BB8 (90 degrees; par: (LB = -0.99, UB = 0); par2: (LB = 0, UB = Inf))
33 = Rotated Clayton (270 degrees; par: (LB = -Inf, UB = 0))
34 = Rotated Gumbel (270 degrees; par: (LB = -Inf, UB = -0.99])
36 = Rotated Joe (270 degrees; par: (LB = -Inf, UB = -0.99))
37 = Rotated BB1 (270 degrees; par: (LB = -Inf, UB = 0); par2: [LB = 0.99, UB = Inf))
38 = Rotated BB6 (270 degrees; par: (LB = -Inf, UB = -0.99]; par2: (LB = 0, UB = Inf))
39 = Rotated BB7 (270 degrees; par: (LB = -Inf, UB = -0.99]; par2: (LB = 0, UB = Inf))
40 = Rotated BB8 (270 degrees; par: (LB = -0.99, UB = 0); par2: (LB = 0, UB = Inf))

See the VineCopula package (Nagler et al., 2025) for additional details on copula parameterization.

Value

A numeric vector of technical efficiency values (one per observation).

Author(s)

Woraphon Yamaka, Paravee Maneejuk, and Nuttaphong Kaewtathip

References

Pipitpojanakarn, V., Maneejuk, P., Yamaka, W., & Sriboonchitta, S. (2016). Analysis of agricultural production in Asia and measurement of technical efficiency using copula-based stochastic frontier quantile model. In International Symposium on Integrated Uncertainty in Knowledge Modelling and Decision Making (pp. 701–714). Springer.

Nagler, T., Schepsmeier, U., Stoeber, J., Brechmann, E. C., Graeler, B., Erhardt, T., ... & Killiches, M. (2025, July). Package ‘VineCopula’.

Examples



sim_data <- sfa_simu_quantile(
  n = 100,
  beta = c(1, 0.5, -0.3),
  sigV = 1,
  sigU = 1,
  tau = 0.5,
  family = 1,
  rho = 0.3,
  seed = 123
)

model <- copSQM(
  Y = sim_data$Y,
  X = sim_data$X,
  family = 1,
  tau = 0.5,
  RHO = 0.3,
  LB = -0.99,
  UB = 0.99,
  nSim = 50,
  seed = 123
)

te_values <- TE(
  theta = model$result[, "Estimate"],
  Y = sim_data$Y,
  X = sim_data$X,
  family = 1,
  tau = 0.5,
  nSim = 200,
  seed = 123,
  plot = TRUE
)

head(te_values)

Copula-Based Stochastic Frontier Quantile Model

Description

Estimates a copula-based stochastic frontier quantile model using simulated maximum likelihood. Dependence between the two-sided noise component and the one-sided inefficiency component is captured via a bivariate copula from the VineCopula package (Nagler et al., 2025).

Usage

copSQM(Y, X, family, tau,
       RHO, LB, UB,
       RHO2 = NULL, LB2 = NULL, UB2 = NULL,
       nSim = 50,
       seed = NULL,
       maxit = 10000)

Arguments

Y

A numeric vector of dependent variable observations.

X

A numeric matrix (or object coercible to a matrix) of explanatory variables. The intercept is added internally.

family

Integer specifying the bivariate copula family (see VineCopula::BiCopPDF). Examples include 1 = Gaussian, 2 = Student-t. Some copula families require a second parameter.

tau

Quantile level in (0,1) for the asymmetric Laplace distribution (ALD).

RHO

Initial value of the first copula parameter (passed to par in VineCopula). Note that this parameter is not always a correlation; its interpretation depends on the chosen copula family.

LB

Lower bound for the first copula parameter RHO.

UB

Upper bound for the first copula parameter RHO.

RHO2

Optional initial value of the second copula parameter (passed to par2 in VineCopula). This is required for copula families that have two parameters (e.g. Student-t where RHO2 represents degrees of freedom, and BB families where RHO2 is a shape parameter). If NULL, a one-parameter copula is assumed.

LB2

Optional lower bound for the second copula parameter RHO2. Must be supplied when RHO2 is not NULL.

UB2

Optional upper bound for the second copula parameter RHO2. Must be supplied when RHO2 is not NULL.

nSim

Number of Monte Carlo draws used to approximate the likelihood integral for each observation. Larger values improve accuracy but increase computation time.

seed

Optional integer seed for reproducibility of the simulation draws used in the likelihood. If NULL, no seed is set.

maxit

Maximum number of iterations for the optimizer (stats::optim with "L-BFGS-B").

Details

The model follows the stochastic frontier decomposition

Y_i = x_i^\top \beta + v_i - u_i,

where v_i is a two-sided noise term and u_i \ge 0 is the inefficiency term. Both components are modeled using the asymmetric Laplace distribution (ALD) at quantile level tau, with u_i truncated to (0,\infty).

Dependence between u_i and v_i is introduced through a bivariate copula density c(\cdot,\cdot) from the VineCopula package. The log-likelihood is evaluated by simulating draws of u_i from the truncated ALD and approximating the likelihood integral by Monte Carlo averaging (nSim draws per observation).

When a two-parameter copula family is used, RHO2 must be provided along with bounds LB2 and UB2.

The following copula families are supported, together with their parameter bounds: 

1  = Gaussian copula (par: (LB = -0.99, UB = 0.99))
2  = Student t copula (par: (LB = -0.99, UB = 0.99); par2: (LB = 0, UB = Inf))
3  = Clayton copula (par: (LB = 0.1, UB = Inf))
4  = Gumbel copula (par: [LB = 0.99, UB = Inf))
5  = Frank copula (par: (LB = -Inf, UB =  0) U (LB = 0, UB = Inf))
6  = Joe copula (par: (LB = 0.99, UB = Inf))
7  = BB1 (Clayton-Gumbel) copula (par: (LB = 0, UB = Inf); par2: [LB = 0.99, UB = Inf))
8  = BB6 copula (par: (LB = 0.99, UB = Inf); par2: (LB = 0, UB = Inf))
9  = BB7 (Joe-Clayton) copula (par: [LB = 0.99, UB = Inf); par2: (LB = 0, UB = Inf))
10 = BB8 copula (par: (LB = 0, UB = 1); par2: (LB = 0, UB =  Inf))
13 = Survival Clayton (180 degrees rotation; par: (LB = 0, UB = Inf))
14 = Survival Gumbel (180 degrees rotation; par: [LB = 0.99, UB = Inf))
16 = Survival Joe (180 degrees rotation; par: (LB = 0.99, UB = Inf))
17 = Survival BB1 (par: (LB = 0, UB = Inf); par2: [LB = 0.99, UB = Inf))
18 = Survival BB6 (par: (LB = 0.99, UB = Inf); par2: (0, UB = Inf))
19 = Survival BB7 (par: [LB = 0.99, UB = Inf); par2: (0, UB = Inf))
20 = Survival BB8 (par: (LB = 0, 0.99); par2: (LB = 0, UB = Inf))
23 = Rotated Clayton (90 degrees; par: (LB = -Inf, UB = 0))
24 = Rotated Gumbel (90 degrees; par: (LB = -Inf, UB = -0.99])
26 = Rotated Joe (90 degrees; par: (LB = -Inf, UB = -0.99))
27 = Rotated BB1 (90 degrees; par: (LB = -Inf, UB = 0); par2: [LB = 0.99, UB = Inf))
28 = Rotated BB6 (90 degrees; par: (LB = -Inf, UB = -0.99]; par2: (LB = 0, UB = Inf))
29 = Rotated BB7 (90 degrees; par: (LB = -Inf, UB = -0.99]; par2: (LB = 0, UB = Inf))
30 = Rotated BB8 (90 degrees; par: (LB = -0.99, UB = 0); par2: (LB = 0, UB = Inf))
33 = Rotated Clayton (270 degrees; par: (LB = -Inf, UB = 0))
34 = Rotated Gumbel (270 degrees; par: (LB = -Inf, UB = -0.99])
36 = Rotated Joe (270 degrees; par: (LB = -Inf, UB = -0.99))
37 = Rotated BB1 (270 degrees; par: (LB = -Inf, UB = 0); par2: [LB = 0.99, UB = Inf))
38 = Rotated BB6 (270 degrees; par: (LB = -Inf, UB = -0.99]; par2: (LB = 0, UB = Inf))
39 = Rotated BB7 (270 degrees; par: (LB = -Inf, UB = -0.99]; par2: (LB = 0, UB = Inf))
40 = Rotated BB8 (270 degrees; par: (LB = -0.99, UB = 0); par2: (LB = 0, UB = Inf))

See the VineCopula package (Nagler et al., 2025) for additional details on copula parameterization.

Value

A list with components:

result

A matrix of parameter estimates, standard errors, z-values, and p-values.

AIC

Akaike Information Criterion.

BIC

Bayesian Information Criterion.

Loglikelihood

Maximized log-likelihood value.

convergence

Convergence code returned by stats::optim. A value of 0 indicates successful convergence.

Author(s)

Woraphon Yamaka, Paravee Maneejuk and Nuttaphong Kaewtathip

References

Pipitpojanakarn, V., Maneejuk, P., Yamaka, W., & Sriboonchitta, S. (2016). Analysis of agricultural production in Asia and measurement of technical efficiency using copula-based stochastic frontier quantile model. In International Symposium on Integrated Uncertainty in Knowledge Modelling and Decision Making (pp. 701–714). Springer.

Pipitpojanakarn, V., Yamaka, W., Sriboonchitta, S., & Maneejuk, P. (2017). Frontier Quantile Model Using a Generalized Class of Skewed Distributions. Advanced Science Letters, 23(11), 10737–10742.

Nagler, T., Schepsmeier, U., Stoeber, J., Brechmann, E. C., Graeler, B., Erhardt, T., ... & Killiches, M. (2025, July). Package ‘VineCopula’.

Examples


set.seed(123)

sim_data <- sfa_simu_quantile(
  n = 50,
  beta = c(1, 0.5, -0.3),
  sigV = 1,
  sigU = 1,
  tau = 0.5,
  family = 1,
  rho = 0.3,
  seed = 123
)

model <- copSQM(
  Y = sim_data$Y,
  X = sim_data$X,
  family = 1,
  tau = 0.5,
  RHO = 0.3,
  LB = -0.99,
  UB = 0.99,
  nSim = 50,
  seed = 123
)
model

## Example with a two-parameter copula (e.g., Student-t: family = 2)
## Here RHO2 typically represents degrees of freedom and must be bounded.
model_t <- copSQM(
  Y = sim_data$Y,
  X = sim_data$X,
  family = 2,
  tau = 0.5,
  RHO = 0.3,
  LB = -0.99,
  UB = 0.99,
  RHO2 = 4,
  LB2 = 2.01,
  UB2 = 50,
  nSim = 50,
  seed = 123
)
model_t

Simulation of Copula-Based Stochastic Frontier Quantile Model

Description

Generates simulated data from a copula-based stochastic frontier quantile model with asymmetric Laplace noise and a truncated asymmetric Laplace inefficiency term.

Usage

sfa_simu_quantile(n,
                  beta,
                  sigV,
                  sigU,
                  tau = 0.5,
                  family = 1,
                  rho = 0.5,
                  rho2 = NULL,
                  seed = NULL)

Arguments

n

Integer. Number of observations to simulate.

beta

Numeric vector of regression coefficients including an intercept. Its length must match the number of regressors plus one (intercept). This simulator generates two regressors, so length(beta) must be 3.

sigV

Positive numeric value. Scale parameter of the two-sided noise term V.

sigU

Positive numeric value. Scale parameter of the one-sided inefficiency term U.

tau

Quantile level in (0,1). Default is 0.5.

family

Integer specifying the bivariate copula family (see VineCopula::BiCopSim). Examples include 1 (Gaussian) and 2 (Student-t). Some families require a second parameter rho2.

rho

First copula parameter (passed to par in VineCopula). Interpretation depends on the copula family and is not always a correlation.

rho2

Optional second copula parameter (passed to par2 in VineCopula). This must be provided for copula families that require a second parameter (e.g., degrees of freedom for the Student-t copula, or a shape parameter for BB copulas).

seed

Optional integer seed for reproducibility. If NULL, no seed is set.

Details

This function simulates data from the stochastic frontier model:

Y = X\beta + V - U,

where:

The simulator constructs two regressors internally (denoted x_1 and x_2).

The following copula families are supported, together with their parameter bounds: 

1  = Gaussian copula (par: (LB = -0.99, UB = 0.99))
2  = Student t copula (par: (LB = -0.99, UB = 0.99); par2: (LB = 0, UB = Inf))
3  = Clayton copula (par: (LB = 0.1, UB = Inf))
4  = Gumbel copula (par: [LB = 0.99, UB = Inf))
5  = Frank copula (par: (LB = -Inf, UB =  0) U (LB = 0, UB = Inf))
6  = Joe copula (par: (LB = 0.99, UB = Inf))
7  = BB1 (Clayton-Gumbel) copula (par: (LB = 0, UB = Inf); par2: [LB = 0.99, UB = Inf))
8  = BB6 copula (par: (LB = 0.99, UB = Inf); par2: (LB = 0, UB = Inf))
9  = BB7 (Joe-Clayton) copula (par: [LB = 0.99, UB = Inf); par2: (LB = 0, UB = Inf))
10 = BB8 copula (par: (LB = 0, UB = 1); par2: (LB = 0, UB =  Inf))
13 = Survival Clayton (180 degrees rotation; par: (LB = 0, UB = Inf))
14 = Survival Gumbel (180 degrees rotation; par: [LB = 0.99, UB = Inf))
16 = Survival Joe (180 degrees rotation; par: (LB = 0.99, UB = Inf))
17 = Survival BB1 (par: (LB = 0, UB = Inf); par2: [LB = 0.99, UB = Inf))
18 = Survival BB6 (par: (LB = 0.99, UB = Inf); par2: (0, UB = Inf))
19 = Survival BB7 (par: [LB = 0.99, UB = Inf); par2: (0, UB = Inf))
20 = Survival BB8 (par: (LB = 0, 0.99); par2: (LB = 0, UB = Inf))
23 = Rotated Clayton (90 degrees; par: (LB = -Inf, UB = 0))
24 = Rotated Gumbel (90 degrees; par: (LB = -Inf, UB = -0.99])
26 = Rotated Joe (90 degrees; par: (LB = -Inf, UB = -0.99))
27 = Rotated BB1 (90 degrees; par: (LB = -Inf, UB = 0); par2: [LB = 0.99, UB = Inf))
28 = Rotated BB6 (90 degrees; par: (LB = -Inf, UB = -0.99]; par2: (LB = 0, UB = Inf))
29 = Rotated BB7 (90 degrees; par: (LB = -Inf, UB = -0.99]; par2: (LB = 0, UB = Inf))
30 = Rotated BB8 (90 degrees; par: (LB = -0.99, UB = 0); par2: (LB = 0, UB = Inf))
33 = Rotated Clayton (270 degrees; par: (LB = -Inf, UB = 0))
34 = Rotated Gumbel (270 degrees; par: (LB = -Inf, UB = -0.99])
36 = Rotated Joe (270 degrees; par: (LB = -Inf, UB = -0.99))
37 = Rotated BB1 (270 degrees; par: (LB = -Inf, UB = 0); par2: [LB = 0.99, UB = Inf))
38 = Rotated BB6 (270 degrees; par: (LB = -Inf, UB = -0.99]; par2: (LB = 0, UB = Inf))
39 = Rotated BB7 (270 degrees; par: (LB = -Inf, UB = -0.99]; par2: (LB = 0, UB = Inf))
40 = Rotated BB8 (270 degrees; par: (LB = -0.99, UB = 0); par2: (LB = 0, UB = Inf))

See the VineCopula package (Nagler et al., 2025) for additional details on copula parameterization.

Value

A list containing:

Y

Numeric vector of simulated dependent variable values.

X

Matrix of simulated independent variables (without intercept column).

V

Simulated two-sided noise term.

U

Simulated non-negative inefficiency term.

copula_uniforms

The simulated dependent uniforms from the copula (two columns).

Author(s)

Woraphon Yamaka, Paravee Maneejuk and Nuttaphong Kaewtathip

References

Koenker, R. (2005). Quantile Regression. Cambridge University Press.

Joe, H. (2014). Dependence Modeling with Copulas. Chapman & Hall/CRC.

Nagler, T., Schepsmeier, U., Stoeber, J., Brechmann, E. C., Graeler, B., Erhardt, T., ... & Killiches, M. (2025, July). Package ‘VineCopula’.

Examples

sim_data <- sfa_simu_quantile(
  n = 100,
  beta = c(1, 0.5, -0.3),
  sigV = 1,
  sigU = 1,
  tau = 0.5,
  family = 1,
  rho = 0.3,
  seed = 123
)
str(sim_data)