Skip to contents

Fits a regression model in the \(L_q\) norm (also labeled as the \(L_p\) norm). In more detail, the optimization function \( \sum_i | y_i - x_i \beta | ^p\) is optimized. The nondifferentiable function is approximated by a differentiable approximation, i.e., we use \(|x| \approx \sqrt{x^2 + \varepsilon } \). The power \(p\) can also be estimated by using est_pow=TRUE, see Giacalone, Panarello and Mattera (2018). The algorithm iterates between estimating regression coefficients and the estimation of power values. The estimation of the power based on a vector of residuals e can be conducted using the function lq_fit_estimate_power.

Using the \(L_q\) norm in the regression is equivalent to assuming an expontial power function for residuals (Giacalone et al., 2018). The density function and a simulation function is provided by dexppow and rexppow, respectively. See also the normalp package.

Usage

lq_fit(y, X, w=NULL, pow=2, eps=0.001, beta_init=NULL, est_pow=FALSE, optimizer="optim",
    eps_vec=10^seq(0,-10, by=-.5), conv=1e-4, miter=20, lower_pow=.1, upper_pow=5)

lq_fit_estimate_power(e, pow_init=2, lower_pow=.1, upper_pow=10)

dexppow(x, mu=0, sigmap=1, pow=2, log=FALSE)

rexppow(n, mu=0, sigmap=1, pow=2, xbound=100, xdiff=.01)

Arguments

y

Dependent variable

X

Design matrix

w

Optional vector of weights

pow

Power \(p\) in \(L_q\) norm

est_pow

Logical indicating whether power should be estimated

eps

Parameter governing the differentiable approximation

e

Vector of resiuals

pow_init

Initial value of power

beta_init

Initial vector

optimizer

Can be "optim" or "nlminb".

eps_vec

Vector with decreasing \(\varepsilon\) values used in optimization

conv

Convergence criterion

miter

Maximum number of iterations

lower_pow

Lower bound for estimated power

upper_pow

Upper bound for estimated power

x

Vector

mu

Location parameter

sigmap

Scale parameter

log

Logical indicating whether the logarithm should be provided

n

Sample size

xbound

Lower and upper bound for density approximation

xdiff

Grid width for density approximation

Value

List with following several entries

coefficients

Vector of coefficients

res_optim

Results of optimization

...

More values

References

Giacalone, M., Panarello, D., & Mattera, R. (2018). Multicollinearity in regression: an efficiency comparison between $L_p$-norm and least squares estimators. Quality & Quantity, 52(4), 1831-1859. doi:10.1007/s11135-017-0571-y

Examples

#############################################################################
# EXAMPLE 1: Small simulated example with fixed power
#############################################################################

set.seed(98)
N <- 300
x1 <- stats::rnorm(N)
x2 <- stats::rnorm(N)
par1 <- c(1,.5,-.7)
y <- par1[1]+par1[2]*x1+par1[3]*x2 + stats::rnorm(N)
X <- cbind(1,x1,x2)

#- lm function in stats
mod1 <- stats::lm.fit(y=y, x=X)

#- use lq_fit function
mod2 <- sirt::lq_fit( y=y, X=X, pow=2, eps=1e-4)
mod1$coefficients
mod2$coefficients

if (FALSE) {
#############################################################################
# EXAMPLE 2: Example with estimated power values
#############################################################################

#*** simulate regression model with residuals from the exponential power distribution
#*** using a power of .30
set.seed(918)
N <- 2000
X <- cbind( 1, c(rep(1,N), rep(0,N)) )
e <- sirt::rexppow(n=2*N, pow=.3, xdiff=.01, xbound=200)
y <- X %*% c(1,.5) + e

#*** estimate model
mod <- sirt::lq_fit( y=y, X=X, est_pow=TRUE, lower_pow=.1)
mod1 <- stats::lm( y ~ 0 + X )
mod$coefficients
mod$pow
mod1$coefficients
}