ADMM on l1 regularized quadratic program

Description

Applies Alternating Direction Method of Multipliers to the l1-regularized quadratic program

f(\mathbf{x}) + g(\mathbf{x}) = \frac{1}{2}\mathbf{x}^TA\mathbf{x} - d^T\mathbf{x} + \lambda |\mathbf{x}|_1

Usage

ADMM_EN2(R, d, x0, lam, mu, maxits, tol, quiet, selector = rep(1, dim(x)[1]))

Arguments

Details

This function is used by other functions and should only be called explicitly for debugging purposes.

Value

ADMM_EN2 returns an object of class "ADMM_EN2" including a list with the following named components

call

The matched call.

x

Found solution.

y

Dual solution.

z

Slack variables.

k

Number of iterations used.

See Also

Used by: SDAD and the SDADcv cross-validation version.

ADMM on l1 regularized quadratic program

Description

Applies Alternating Direction Method of Multipliers to the l1-regularized quadratic program

f(\mathbf{x}) + g(\mathbf{x}) = \frac{1}{2}\mathbf{x}^TA\mathbf{x} - d^T\mathbf{x} + \lambda |\mathbf{x}|_1

Usage

ADMM_EN_SMW(Ainv, V, R, d, x0, lam, mu, maxits, tol, quiet, selector)

Arguments

Details

This function is used by other functions and should only be called explicitly for debugging purposes.

Value

ADMM_EN_SMW returns an object of class "ADMM_EN_SMW" including a list with the following named components

call

The matched call.

x

Found solution.

y

Dual solution.

z

Slack variables.

k

Number of iterations used.

See Also

Used by: SDAD and the SDADcv cross-validation version.

Accelerated Proximal Gradient on l1 regularized quadratic program

Description

Applies accelerated proximal gradient algorithm to the l1-regularized quadratic program

f(\mathbf{x}) + g(\mathbf{x}) = \frac{1}{2}\mathbf{x}^TA\mathbf{x} - d^T\mathbf{x} + \lambda |\mathbf{x}|_1

Usage

APG_EN2(A, d, x0, lam, alpha, maxits, tol, selector = rep(1, dim(x0)[1]))

Arguments

Details

This function is used by other functions and should only be called explicitly for debugging purposes.

Value

APG_EN2 returns an object of class "APG_EN2" including a list with the following named components

call

The matched call.

x

Found solution.

k

Number of iterations used.

See Also

Used by: SDAAP and the SDAAPcv cross-validation version.

Accelerated Proximal Gradient (with backtracking) on l1 regularized quadratic program

Description

Applies accelerated proximal gradient algorithm (with backtracking) to the l1-regularized quadratic program

f(\mathbf{x}) + g(\mathbf{x}) = \frac{1}{2}\mathbf{x}^TA\mathbf{x} - d^T\mathbf{x} + \lambda |\mathbf{x}|_1

Usage

APG_EN2bt(
  A,
  Xt,
  Om,
  gamma,
  d,
  x0,
  lam,
  L,
  eta,
  maxits,
  tol,
  selector = rep(1, dim(x0)[1])
)

Arguments

Details

This function is used by other functions and should only be called explicitly for debugging purposes.

Value

APG_EN2bt returns an object of class "APG_EN2bt" including a list with the following named components

call

The matched call.

x

Found solution.

k

Number of iterations used.

See Also

Used by: SDAAP and the SDAAPcv cross-validation version.

Accelerated Proximal Gradient on l1 regularized quadratic program

Description

Applies accelerated proximal gradient algorithm to the l1-regularized quadratic program (with rank reduced Omega inside A)

f(\mathbf{x}) + g(\mathbf{x}) = \frac{1}{2}\mathbf{x}^TA\mathbf{x} - d^T\mathbf{x} + \lambda |\mathbf{x}|_1

Usage

APG_EN2rr(A, d, x0, lam, alpha, maxits, tol, selector = rep(1, dim(x0)[1]))

Arguments

Details

This function is used by other functions and should only be called explicitly for debugging purposes.

Value

APG_EN2 returns an object of class "APG_EN2" including a list with the following named components

call

The matched call.

x

Found solution.

k

Number of iterations used.

See Also

Used by: SDAAP and the SDAAPcv cross-validation version.

Accelerated Sparse Discriminant Analysis

Description

Applies accelerated proximal gradient algorithm, proximal gradient algorithm or alternating direction methods of multipliers algorithm to the optimal scoring formulation of sparse discriminant analysis proposed by Clemmensen et al. 2011.

argmin{|(Y_t\theta-X_t\beta)|_2^2 + t|\beta|_1 + \lambda|\beta|_2^2}

Usage

ASDA(Xt, ...)

## Default S3 method:
ASDA(
  Xt,
  Yt,
  Om = diag(p),
  gam = 0.001,
  lam = 1e-06,
  q = K - 1,
  method = "SDAAP",
  control = list(),
  ...
)

Arguments

Details

The control list contains the following entries to further tune the algorithms.

PGsteps

Maximum number if inner proximal gradient/ADMM algorithm for finding beta. Default value is 1000.

PGtol

Stopping tolerance for inner method. If the method is SDAD, then this must be a vector of two values, absolute (first element) and relative tolerance (second element). Default value is 1e-5 for both absolute and relative tolerances.

maxits

Number of iterations to run. Default value is 250.

tol

Stopping tolerance. Default value is 1e-3.

mu

Penalty parameter for augmented Lagrangian term, must be greater than zero and only needs to be specified when using method SDAD. Default value is 1.

CV

Logical value which is TRUE if cross validation is supposed to be performed. If cross-validation is performed, then lam should be specified as a vector containing the regularization values to be tested. Default value is FALSE.

folds

Integer determining the number of folds in cross-validation. Not needed if CV is not specified. Default value is 5.

feat

Maximum fraction of nonzero features desired in validation scheme. Not needed if CV is not specified. Default value is 0.15.

quiet

Set to FALSE if status updates are supposed to be printed to the R console. Default value is TRUE. Note that this triggers a lot of printing to the console.

ordinal

Set to TRUE if the labels are ordinal. Only available for methods SDAAP and SDAD.

initTheta

Option to set the initial theta vector, by default it is a vector of all ones for the first theta.

bt

Logical indicating whether backtracking should be used, only applies to the Proximal Gradient based methods. By default, backtracking is not used.

L

Initial estimate for Lipshitz constant used for backtracking. Default value is 0.25.

eta

Scalar for Lipshitz constant. Default value is 1.25.

rankRed

Boolean indicating whether Om is factorized, such that R^t*R=Om, currently only applicable for accelerated proximal gradient.

Value

ASDA returns an object of class "ASDA" including a list with the following named components:

call

The matched call.

B

p by q matrix of discriminant vectors, i.e. sparse loadings.

Q

K by q matrix of scoring vectors, i.e. optimal scores.

varNames

Names of the predictors used, i.e. column names of Xt.

origP

Number of variables in Xt.

fit

Output from function lda on projected data. This is NULL the trivial solution is found, i.e. B is all zeroes. Use lower values of lam if that is the case.

classes

The classes in Yt.

lambda

The lambda/lam used, best value found by cross- validation if CV is TRUE.

NULL

Note

The input matrix Xt should be normalized, i.e. each column corresponding to a variable should have its mean subtracted and scaled to unit length. The functions normalize and normalizetest are supplied for this purpose in the package.

See Also

SDAAP, SDAP and SDAD

Examples

    set.seed(123)
    # Prepare training and test set
    train <- c(1:40,51:90,101:140)
    Xtrain <- iris[train,1:4]
    nX <- normalize(Xtrain)
    Xtrain <- nX$Xc
    Ytrain <- iris[train,5]
    Xtest <- iris[-train,1:4]
    Xtest <- normalizetest(Xtest,nX)
    Ytest <- iris[-train,5]

    # Define parameters for Alternating Direction Method of Multipliers (SDAD)
    Om <- diag(4)+0.1*matrix(1,4,4) #elNet coef mat
    gam <- 0.0001
    lam <- 0.0001
    method <- "SDAD"
    q <- 2
    control <- list(PGsteps = 100,
                    PGtol = c(1e-5,1e-5),
                    mu = 1,
                    maxits = 100,
                    tol = 1e-3,
                    quiet = FALSE)

    # Run the algorithm
    res <- ASDA(Xt = Xtrain,
                Yt = Ytrain,
                Om = Om,
                gam = gam ,
                lam = lam,
                q = q,
                method = method,
                control = control)

    # Can also just use the defaults, which is Accelerated Proximal Gradient (SDAAP):
    resDef <- ASDA(Xtrain,Ytrain)

    # Some example on simulated data
    # Generate Gaussian data on three classes with plenty of redundant variables

    # This example shows the basic steps on how to apply this to data, i.e.:
    #  1) Setup training data
    #  2) Normalize
    #  3) Train
    #  4) Predict
    #  5) Plot projected data
    #  6) Accuracy on test set

    P <- 300 # Number of variables
    N <- 50 # Number of samples per class

    # Mean for classes, they are zero everywhere except the first 3 coordinates
    m1 <- rep(0,P)
    m1[1] <- 3

    m2 <- rep(0,P)
    m2[2] <- 3

    m3 <- rep(0,P)
    m3[3] <- 3

    # Sample dummy data
    Xtrain <- rbind(MASS::mvrnorm(n=N,mu = m1, Sigma = diag(P)),
                   MASS::mvrnorm(n=N,mu = m2, Sigma = diag(P)),
                   MASS::mvrnorm(n=N,mu = m3, Sigma = diag(P)))

    Xtest <- rbind(MASS::mvrnorm(n=N,mu = m1, Sigma = diag(P)),
                   MASS::mvrnorm(n=N,mu = m2, Sigma = diag(P)),
                   MASS::mvrnorm(n=N,mu = m3, Sigma = diag(P)))

    # Generate the labels
    Ytrain <- factor(rep(1:3,each=N))
    Ytest <- Ytrain

    # Normalize the data
    Xt <- accSDA::normalize(Xtrain)
    Xtrain <- Xt$Xc # Use the centered and scaled data
    Xtest <- accSDA::normalizetest(Xtest,Xt)

    # Train the classifier and increase the sparsity parameter from the default
    # so we penalize more for non-sparse solutions.
    res <- accSDA::ASDA(Xtrain,Ytrain,lam=0.01)

    # Plot the projected training data, it is projected to
    # 2-dimension because we have 3 classes. The number of discriminant
    # vectors is maximum number of classes minus 1.
    XtrainProjected <- Xtrain%*%res$beta

    plot(XtrainProjected[,1],XtrainProjected[,2],col=Ytrain)

    # Predict on the test data
    preds <- predict(res, newdata = Xtest)

    # Plot projected test data with predicted and correct labels
    XtestProjected <- Xtest%*%res$beta

    plot(XtestProjected[,1],XtestProjected[,2],col=Ytest,
         main="Projected test data with original labels")
    plot(XtestProjected[,1],XtestProjected[,2],col=preds$class,
         main="Projected test data with predicted labels")

    # Calculate accuracy
    sum(preds$class == Ytest)/(3*N) # We have N samples per class, so total 3*N

barplot for ASDA objects

Description

This is a function to visualize the discriminant vector from the ASDA method. The plot is constructed as a ggplot barplot and the main purpose of it is to visually inspect the sparsity of the discriminant vectors. The main things to look for are how many parameters are non-zero and if there is any structure in the ones that are non-zero, but the structure is dependent on the order you specify your variables. For time-series data, this could mean that a chunk of variables are non-zero that are close in time, meaning that there is some particular event that is best for discriminating between the classes that you have.

Usage

ASDABarPlot(asdaObj, numDVs = 1, xlabel, ylabel, getList = FALSE, main, ...)

Arguments

Value

barplot.ASDA returns either a single combined plot or a list of individual ggplot objects.

Note

This function is used as a quick diagnostics tool for the output from the ASDA function. Feel free to look at the code to customize the plots in any way you like.

See Also

ASDA

Examples

    # Generate and ASDA object with your data, e.g.
    # Prepare training and test set
    # This is a very small data set, I advise you to try it on something with more
    # variables, e.g. something from this source: http://www.cs.ucr.edu/~eamonn/time_series_data/
    # or possibly run this on the Gaussian data example from the ASDA function
    train <- c(1:40,51:90,101:140)
    Xtrain <- iris[train,1:4]
    nX <- normalize(Xtrain)
    Xtrain <- nX$Xc
    Ytrain <- iris[train,5]
    Xtest <- iris[-train,1:4]
    Xtest <- normalizetest(Xtest,nX)
    Ytest <- iris[-train,5]
    # Run the method
    resIris <- ASDA(Xtrain,Ytrain)

    # Look at the barplots of the DVs
    ASDABarPlot(resIris)

Sparse Discriminant Analysis solved via Accelerated Proximal Gradient

Description

Applies accelerated proximal gradient algorithm to the optimal scoring formulation of sparse discriminant analysis proposed by Clemmensen et al. 2011.

Usage

SDAAP(Xt, ...)

## Default S3 method:
SDAAP(
  Xt,
  Yt,
  Om,
  gam,
  lam,
  q,
  PGsteps,
  PGtol,
  maxits,
  tol,
  selector = rep(1, dim(Xt)[2]),
  initTheta,
  bt = FALSE,
  L,
  eta,
  rankRed = FALSE,
  ...
)

Arguments

Value

SDAAP returns an object of class "SDAAP" including a list with the following named components:

call

The matched call.

B

p by q matrix of discriminant vectors.

Q

K by q matrix of scoring vectors.

subits

Total number of iterations in proximal gradient subroutine.

totalits

Number coordinate descent iterations for all discriminant vectors

NULL

See Also

SDAAPcv, SDAP and SDAD

Sparse Discriminant Analysis solved via ADMM

Description

Applies alternating direction methods of multipliers algorithm to the optimal scoring formulation of sparse discriminant analysis proposed by Clemmensen et al. 2011.

Usage

SDAD(Xt, ...)

## Default S3 method:
SDAD(
  Xt,
  Yt,
  Om,
  gam,
  lam,
  mu,
  q,
  PGsteps,
  PGtol,
  maxits,
  tol,
  selector = rep(1, dim(Xt)[2]),
  initTheta,
  ...
)

Arguments

Value

SDAD returns an object of class "SDAD" including a list with the following named components: (More will be added later to handle the predict function)

call

The matched call.

B

p by q matrix of discriminant vectors.

Q

K by q matrix of scoring vectors.

subits

Total number of iterations in proximal gradient subroutine.

totalits

Number coordinate descent iterations for all discriminant vectors

NULL

See Also

SDADcv, SDAAP and SDAP

Sparse Discriminant Analysis solved via Proximal Gradient

Description

Applies proximal gradient algorithm to the optimal scoring formulation of sparse discriminant analysis proposed by Clemmensen et al. 2011.

Usage

SDAP(Xt, ...)

## Default S3 method:
SDAP(
  Xt,
  Yt,
  Om,
  gam,
  lam,
  q,
  PGsteps,
  PGtol,
  maxits,
  tol,
  initTheta,
  bt = FALSE,
  L,
  eta,
  ...
)

Arguments

Value

SDAP returns an object of class "SDAP" including a list with the following named components: (More will be added later to handle the predict function)

call

The matched call.

B

p by q matrix of discriminant vectors.

Q

K by q matrix of scoring vectors.

subits

Total number of iterations in proximal gradient subroutine.

totalits

Number coordinate descent iterations for all discriminant vectors

NULL

See Also

SDAPcv, SDAAP and SDAD

Sparse Zero Variance Discriminant Analysis

Description

Applies SZVD heuristic for sparse zero-variance discriminant analysis to given training set.

Usage

SZVD(train, ...)

## Default S3 method:
SZVD(
  train,
  gamma,
  D,
  penalty = TRUE,
  scaling = TRUE,
  tol = list(abs = 1e-04, rel = 1e-04),
  maxits = 2000,
  beta = 1,
  quiet = TRUE,
  ...
)

Arguments

Details

This function will currently solve as a standalone function in accSDA for time comparison. A wrapper function like ASDA will be created to use the functionality of plots and such. Maybe call it ASZDA. For that purpose the individual ZVD function will need to be implemented.

Value

SZVD returns an object of class "SZVD" including a list with the following named components:

DVs

Discriminant vectors.

its

Number of iterations required to find DVs.

pen_scal

Weights used in reweighted l1-penalty.

N

Basis for the null-space of the sample within-class covariance.

means

Training class-means.

mus

Training meand and variance scaling/centering terms.

w0

unpenalized zero-variance discriminants (initial solutions) plus B and W, etc.

NULL

See Also

Used by: SZVDcv.

Examples

  set.seed(123)
  P <- 300 # Number of variables
  N <- 50 # Number of samples per class

  # Mean for classes, they are zero everywhere except the first 3 coordinates
  m1 <- rep(0,P)
  m1[1] <- 3

  m2 <- rep(0,P)
  m2[2] <- 3

  m3 <- rep(0,P)
  m3[3] <- 3

  # Sample dummy data
  Xtrain <- rbind(MASS::mvrnorm(n=N,mu = m1, Sigma = diag(P)),
                 MASS::mvrnorm(n=N,mu = m2, Sigma = diag(P)),
                MASS::mvrnorm(n=N,mu = m3, Sigma = diag(P)))


  # Generate the labels
  Ytrain <- rep(1:3,each=N)

  # Normalize the data
  Xt <- accSDA::normalize(Xtrain)
  Xtrain <- Xt$Xc

  # Train the classifier and increase the sparsity parameter from the default
  # so we penalize more for non-sparse solutions.
  res <- accSDA::SZVD(cbind(Ytrain,Xtrain),beta=2.5,
                     maxits=1000,tol = list(abs = 1e-04, rel = 1e-04))

Alternating Direction Method of Multipliers for SZVD

Description

Iteratively solves the problem

\min(-1/2*x^TB^Tx + \gamma p(y): ||x||_2 \leq 1, DNx = y)

Usage

SZVD_ADMM(B, N, D, sols0, pen_scal, gamma, beta, tol, maxits, quiet = TRUE)

Arguments

Details

This function is used by other functions and should only be called explicitly for debugging purposes.

Value

SZVD_ADMM returns an object of class "SZVD_ADMM" including a list with the following named components

x,y,z

Iterates at termination.

its

Number of iterations required to converge.

errtol

Stopping error bound at termination

See Also

Used by: SZVDcv.

Cross-validation of sparse zero variance discriminant analysis

Description

Applies alternating direction methods of multipliers to solve sparse zero variance discriminant analysis.

Usage

SZVD_kFold_cv(X, ...)

## Default S3 method:
SZVD_kFold_cv(
  X,
  Y,
  folds,
  gams,
  beta,
  D,
  q,
  maxits,
  tol,
  ztol,
  feat,
  penalty,
  quiet,
  ...
)

Arguments

Details

Add how max_gamma is calculated from the ZVD solution. This function might require a wrapper similar to ASDA.

Value

SZVDcv returns an object of class "SZVDcv" including a list with the named components DVs and gambest. Where DVs are the discriminant vectors for the best l1 regularization parameter and gambest is the best regularization parameter found in the cross-validation.

NULL

See Also

Used by: SZVDcv.

Examples

  P <- 150 # Number of variables
  N <- 20 # Number of samples per class

  # Mean for classes, they are zero everywhere except the first 3 coordinates
  m1 <- rep(0,P)
  m1[1] <- 3

  m2 <- rep(0,P)
  m2[2] <- 3

  m3 <- rep(0,P)
  m3[3] <- 3

  # Sample dummy data
  Xtrain <- rbind(MASS::mvrnorm(n=N,mu = m1, Sigma = diag(P)),
                 MASS::mvrnorm(n=N,mu = m2, Sigma = diag(P)),
                MASS::mvrnorm(n=N,mu = m3, Sigma = diag(P)))


  # Generate the labels
  Ytrain <- cbind(c(rep(1,N),rep(0,2*N)),
                  c(rep(0,N),rep(1,N),rep(0,N)),
                  c(rep(0,2*N),rep(1,N)))

  # Normalize the data
  Xt <- accSDA::normalize(Xtrain)
  Xtrain <- Xt$Xc

  # Train the classifier and increase the sparsity parameter from the default
  # so we penalize more for non-sparse solutions.
  res <- accSDA::SZVD_kFold_cv(Xtrain,Ytrain,folds=2,gams=2,beta=2.5,q=1, D=diag(P),
                              maxits=50,tol=list(abs=1e-2,rel=1e-2),
                              ztol=1e-3,feat=0.3,quiet=FALSE,penalty=TRUE)

Cross-validation of sparse zero variance discriminant analysis

Description

Applies alternating direction methods of multipliers to solve sparse zero variance discriminant analysis.

Usage

SZVDcv(Atrain, ...)

## Default S3 method:
SZVDcv(
  Atrain,
  Aval,
  k,
  num_gammas,
  g_mults,
  D,
  sparsity_pen,
  scaling,
  penalty,
  beta,
  tol,
  ztol,
  maxits,
  quiet,
  ...
)

Arguments

Details

This function might require a wrapper similar to ASDA.

Value

SZVDcv returns an object of class "SZVDcv" including a list with the following named components:

DVs

Discriminant vectors for the best choice of gamma.

all_DVs

Discriminant vectors for all choices of gamma.

l0_DVs

Discriminant vectors for gamma minimizing cardinality.

mc_DVs

Discriminant vector minimizing misclassification.

gamma

Choice of gamma minimizing validation criterion.

gammas

Set of all gammas trained on.

max_g

Maximum value of gamma guaranteed to yield a nontrivial solution.

ind

Index of best gamma.

w0

unpenalized zero-variance discriminants (initial solutions) plus B and W, etc. from ZVD

NULL

See Also

Non CV version: SZVD.

Examples

  P <- 300 # Number of variables
  N <- 50 # Number of samples per class

  # Mean for classes, they are zero everywhere except the first 3 coordinates
  m1 <- rep(0,P)
  m1[1] <- 3

  m2 <- rep(0,P)
  m2[2] <- 3

  m3 <- rep(0,P)
  m3[3] <- 3

  # Sample dummy data
  Xtrain <- rbind(MASS::mvrnorm(n=N,mu = m1, Sigma = diag(P)),
                 MASS::mvrnorm(n=N,mu = m2, Sigma = diag(P)),
                MASS::mvrnorm(n=N,mu = m3, Sigma = diag(P)))
 Xval <- rbind(MASS::mvrnorm(n=N,mu = m1, Sigma = diag(P)),
                 MASS::mvrnorm(n=N,mu = m2, Sigma = diag(P)),
                MASS::mvrnorm(n=N,mu = m3, Sigma = diag(P)))

  # Generate the labels
  Ytrain <- rep(1:3,each=N)
  Yval <- rep(1:3,each=N)


  # Train the classifier and increase the sparsity parameter from the default
  # so we penalize more for non-sparse solutions.

  res <- accSDA::SZVDcv(cbind(Ytrain,Xtrain),cbind(Yval,Xval),num_gammas=4,
                        g_mults = c(0,1),beta=2.5,
                        D=diag(P), maxits=100,tol=list(abs=1e-3,rel=1e-3), k = 3,
                        ztol=1e-4,sparsity_pen=0.3,quiet=FALSE,penalty=TRUE,scaling=TRUE)

Zero Variance Discriminant Analysis

Description

Implements the ZVD algorithm to solve dicriminant vectors.

Usage

ZVD(A, ...)

## Default S3 method:
ZVD(A, scaling = FALSE, get_DVs = FALSE, ...)

Arguments

Details

This function should potentially be made internal for the release.

Value

SZVDcv returns an object of class "ZVD" including a list with the following named components:

dvs

discriminant vectors (optional).

B

sample between-class covariance.

W

sample within-class covariance.

N

basis for the null space of the sample within-class covariance.

mu

training mean and variance scaling/centering terms

means

vectors of sample class-means.

k

number of classes in given data set.

labels

list of classes.

obs

matrix of data observations.

class_obs

Matrices of observations of each class.

NULL

See Also

Used by: SZVDcv.

Examples

  # Generate Gaussian data on three classes with bunch of redundant variables

  P <- 300 # Number of variables
  N <- 50 # Number of samples per class

  # Mean for classes, they are zero everywhere except the first 3 coordinates
  m1 <- rep(0,P)
  m1[1] <- 3

  m2 <- rep(0,P)
  m2[2] <- 3

  m3 <- rep(0,P)
  m3[3] <- 3

  # Sample dummy data
  Xtrain <- rbind(MASS::mvrnorm(n=N,mu = m1, Sigma = diag(P)),
              MASS::mvrnorm(n=N,mu = m2, Sigma = diag(P)),
              MASS::mvrnorm(n=N,mu = m3, Sigma = diag(P)))


  # Generate the labels
  Ytrain <- rep(1:3,each=N)

  # Normalize the data
  Xt <- accSDA::normalize(Xtrain)
  Xtrain <- Xt$Xc

  # Train the classifier and increase the sparsity parameter from the default
  # so we penalize more for non-sparse solutions.
  res <- accSDA::ZVD(cbind(Ytrain,Xtrain))

accSDA: A package for performing sparse discriminant analysis in various ways.

Description

The accSDA package provides functions to perform sparse discriminant analysis using a selection of three optimization methods, proximal gradient (PG), accelerated proximal gradient (APG) and alternating direction method of multipliers (ADMM). The package is intended to extend the available tools to perform sparse discriminant analysis in R. The three methods can be called from the function ASDA. Cross validation is also implemented for the L1 regularization parameter. Functions for doing predictions, summary, printing and simple plotting are also provided. The sparse discriminant functions perform lda on the projected data by default, using the lda function in the MASS package. The functions return an object of the same class as the name of the function and provide the lda solution, along with the projected data, thus other kinds of classification algorithms can be employed on the projected data.

Generate data for ordinal examples in the package

Description

Given the parameters, the function creates a dataset for testing the ordinal functionality of the package. The data is samples from multivariate Gaussians with different means, where the mean varies along a sinusoidal curve w.r.t. the class label.

Usage

genDat(numClasses, numObsPerClass, mu, sigma)

Arguments

Details

This function is used to demonstrate the usage of the ordinal classifier.

Value

genDat Returns a list with the following attributes:

X

A matrix with two columns and numObsPerClass*numClasses rows.

Y

Labels for the rows of X.

Author(s)

Gudmundur Einarsson

See Also

ordASDA

Examples

set.seed(123)

    # You can play around with these values to generate some 2D data to test one
    numClasses <- 15
    sigma <- matrix(c(1,-0.2,-0.2,1),2,2)
    mu <- c(0,0)
    numObsPerClass <- 5

    # Generate the data, can access with train$X and train$Y
    train <- accSDA::genDat(numClasses,numObsPerClass,mu,sigma)
    test <- accSDA::genDat(numClasses,numObsPerClass*2,mu,sigma)

    # Visualize it, only using the first variable gives very good separation
    plot(train$X[,1],train$X[,2],col = factor(train$Y),asp=1,main="Training Data")

Normalize training data

Description

Normalize a vector or matrix to zero mean and unit length columns.

Usage

normalize(X)

Arguments

Details

This function can e.g. be used for the training data in the ASDA function.

Value

normalize Returns a list with the following attributes:

Xc

The normalized data

mx

Mean of columns of X.

vx

Length of columns of X.

Id

Logical vector indicating which variables are included in X. If some of the columns have zero length they are omitted

Author(s)

Line Clemmensen

References

Clemmensen, L., Hastie, T. and Ersboell, K. (2008) "Sparse discriminant analysis", Technical report, IMM, Technical University of Denmark

See Also

normalizetest, predict.ASDA, ASDA

Examples

## Data
X<-matrix(sample(seq(3),12,replace=TRUE),nrow=3)

## Normalize data
Nm<-normalize(X)
print(Nm$Xc)

## See if any variables have been removed
which(!Nm$Id)

Normalize training data

Description

Normalize test data using output from the normalize() of the training data

Usage

normalizetest(Xtst, Xn)

Arguments

Details

This function can e.g. be used for the test data in the predict.ASDA function.

Value

normalizetest returns the normalized test data Xtst

Author(s)

Line Clemmensen

References

Clemmensen, L., Hastie, T. and Ersboell, K. (2008) "Sparse discriminant analysis", Technical report, IMM, Technical University of Denmark

See Also

normalize, predict.ASDA, ASDA

Examples

## Data
Xtr<-matrix(sample(seq(3),12,replace=TRUE),nrow=3)
Xtst<-matrix(sample(seq(3),12,replace=TRUE),nrow=3)

## Normalize training data
Nm<-normalize(Xtr)

## Normalize test data
Xtst<-normalizetest(Xtst,Nm)

Finding null space of linear operator

Description

Finds the null space of a linear operator A in R^{n \times m}. The null space is given as a matrix, where the columns form an orthonormal basis for the nullspace. This function emulates the null function in matlab, it works exactly the same, but the basis vectors may be different, i.e. rotated.

Usage

nullSp(A)

Arguments

Details

This function is used by other functions and should only be called explicitly for debugging purposes.

Value

nullSp returns a matrix whose columns span the nullspace of A.

See Also

Alternative Null function in MASS package.

Ordinal Accelerated Sparse Discriminant Analysis

Description

Applies accelerated proximal gradient algorithm to the optimal scoring formulation of sparse discriminant analysis proposed by Clemmensen et al. 2011. The problem is further casted to a binary classification problem as described in "Learning to Classify Ordinal Data: The Data Replication Method" by Cardoso and da Costa to handle the ordinal labels. This function serves as a wrapper for the ASDA function, where the appropriate data augmentation is performed. Since the problem is casted into a binary classication problem, only a single discriminant vector comes from the result. The first *p* entries correspond to the variables/coefficients for the predictors, while the following K-1 entries correspond to biases for the found hyperplane, to separate the classes. The resulting object is of class ordASDA and has an accompanying predict function. The paper by Cardoso and dat Costa can be found here: (http://www.jmlr.org/papers/volume8/cardoso07a/cardoso07a.pdf).

Usage

ordASDA(Xt, ...)

## Default S3 method:
ordASDA(
  Xt,
  Yt,
  s = 1,
  Om,
  gam = 0.001,
  lam = 1e-06,
  method = "SDAAP",
  control,
  ...
)

Arguments

Value

ordASDA returns an object of class "ordASDA" including a list with the same components as an ASDA objects and:

h

Scalar value for biases.

K

Number of classes.

NULL

Note

Remember to normalize the data.

See Also

ASDA.

Examples

    set.seed(123)

    # You can play around with these values to generate some 2D data to test one
    numClasses <- 5
    sigma <- matrix(c(1,-0.2,-0.2,1),2,2)
    mu <- c(0,0)
    numObsPerClass <- 5

    # Generate the data, can access with train$X and train$Y
    train <- accSDA::genDat(numClasses,numObsPerClass,mu,sigma)
    test <- accSDA::genDat(numClasses,numObsPerClass*2,mu,sigma)

    # Visualize it, only using the first variable gives very good separation
    plot(train$X[,1],train$X[,2],col = factor(train$Y),asp=1,main="Training Data")

    # Train the ordinal based model
    res <- accSDA::ordASDA(train$X,train$Y,s=2,h=1, gam=1e-6, lam=1e-3)
    vals <- predict(object = res,newdata = test$X) # Takes a while to run ~ 10 seconds
    sum(vals==test$Y)/length(vals) # Get accuracy on test set
    #plot(test$X[,1],test$X[,2],col = factor(test$Y),asp=1,
    #      main="Test Data with correct labels")
    #plot(test$X[,1],test$X[,2],col = factor(vals),asp=1,
    #    main="Test Data with predictions from ordinal classifier")

Predict method for sparse discriminant analysis

Description

Predicted values based on fit from the function ASDA. This function is used to classify new observations based on their explanatory variables/features.

Usage

## S3 method for class 'ASDA'
predict(object, newdata = NULL, ...)

Arguments

Value

A list with components:

class

The classification (a factor)

posterior

posterior probabilities for the classes

x

the scores

Note

The input matrix newdata should be normalized w.r.t. the normalization of the training data

See Also

SDAAP, SDAP and SDAD

Examples

    # Prepare training and test set
    train <- c(1:40,51:90,101:140)
    Xtrain <- iris[train,1:4]
    nX <- normalize(Xtrain)
    Xtrain <- nX$Xc
    Ytrain <- iris[train,5]
    Xtest <- iris[-train,1:4]
    Xtest <- normalizetest(Xtest,nX)
    Ytest <- iris[-train,5]

    # Define parameters for SDAD
    Om <- diag(4)+0.1*matrix(1,4,4) #elNet coef mat
    gam <- 0.01
    lam <- 0.01
    method <- "SDAD"
    q <- 2
    control <- list(PGsteps = 100,
                    PGtol = c(1e-5,1e-5),
                    mu = 1,
                    maxits = 100,
                    tol = 1e-3,
                    quiet = FALSE)

    # Run the algorithm
    res <- ASDA(Xt = Xtrain,
                Yt = Ytrain,
                Om = Om,
                gam = gam ,
                lam = lam,
                q = q,
                method = method,
                control = control)

    # Do the predictions on the test set
    preds <- predict(object = res, newdata = Xtest)

Predict method for ordinal sparse discriminant analysis

Description

Predicted values based on fit from the function ordASDA. This function is used to classify new observations based on their explanatory variables/features. There is no need to normalize the data, the data is normalized based on the normalization data from the ordASDA object.

Usage

## S3 method for class 'ordASDA'
predict(object, newdata = NULL, ...)

Arguments

Value

A vector of predictions.

See Also

ordASDA

Examples

    set.seed(123)

    # You can play around with these values to generate some 2D data to test one
    numClasses <- 5
    sigma <- matrix(c(1,-0.2,-0.2,1),2,2)
    mu <- c(0,0)
    numObsPerClass <- 5

    # Generate the data, can access with train$X and train$Y
    train <- accSDA::genDat(numClasses,numObsPerClass,mu,sigma)
    test <- accSDA::genDat(numClasses,numObsPerClass*2,mu,sigma)

    # Visualize it, only using the first variable gives very good separation
    plot(train$X[,1],train$X[,2],col = factor(train$Y),asp=1,main="Training Data")

    # Train the ordinal based model
    res <- accSDA::ordASDA(train$X,train$Y,s=2,h=1, gam=1e-6, lam=1e-3)
    vals <- predict(object = res,newdata = test$X) # Takes a while to run ~ 10 seconds
    sum(vals==test$Y)/length(vals) # Get accuracy on test set
    #plot(test$X[,1],test$X[,2],col = factor(test$Y),asp=1,
    #      main="Test Data with correct labels")
    #plot(test$X[,1],test$X[,2],col = factor(vals),asp=1,
    #     main="Test Data with predictions from ordinal classifier")

Print method for ASDA object

Description

Prints a summary of the output from the ASDA function. The output summarizes the discriminant analysis in human readable format.

Usage

## S3 method for class 'ASDA'
print(x, digits = max(3, getOption("digits") - 3), numshow = 5, ...)

Arguments

Value

An invisible copy of x.

See Also

ASDA, predict.ASDA and SDAD

Examples

    # Prepare training and test set
    train <- c(1:40,51:90,101:140)
    Xtrain <- iris[train,1:4]
    nX <- normalize(Xtrain)
    Xtrain <- nX$Xc
    Ytrain <- iris[train,5]
    Xtest <- iris[-train,1:4]
    Xtest <- normalizetest(Xtest,nX)
    Ytest <- iris[-train,5]

    # Run the algorithm
    resDef <- ASDA(Xtrain,Ytrain)

    # Print
    print(resDef)

Accelerated Proximal Gradient on l1 regularized quadratic program

Description

Applies accelerated proximal gradient algorithm to the l1-regularized quadratic program

f(\mathbf{x}) + g(\mathbf{x}) = \frac{1}{2}\mathbf{x}^TA\mathbf{x} - d^T\mathbf{x} + \lambda |\mathbf{x}|_1

Usage

prox_EN(A, d, x0, lam, alpha, maxits, tol)

Arguments

Details

This function is used by other functions and should only be called explicitly for debugging purposes.

Value

prox_EN returns an object of class "prox_EN" including a list with the following named components

call

The matched call.

x

Found solution.

k

Number of iterations used.

See Also

Used by: SDAP and the SDAPcv cross-validation version.

Accelerated Proximal Gradient on l1 regularized quadratic program with backtracking

Description

Applies accelerated proximal gradient (with backtracking) algorithm to the l1-regularized quadratic program

f(\mathbf{x}) + g(\mathbf{x}) = \frac{1}{2}\mathbf{x}^TA\mathbf{x} - d^T\mathbf{x} + \lambda |\mathbf{x}|_1

Usage

prox_ENbt(A, Xt, Om, gamma, d, x0, lam, L, eta, maxits, tol)

Arguments

Details

This function is used by other functions and should only be called explicitly for debugging purposes.

Value

prox_ENbt returns an object of class "prox_ENbt" including a list with the following named components

call

The matched call.

x

Found solution.

k

Number of iterations used.

See Also

Used by: SDAP and the SDAPcv cross-validation version.

Classify test data using nearest centroid classification and discriminant vectors learned from the training set.

Description

This function is used in SZVDcv and is only meant for internal use at this stage. Will potentially be released in future versions.

Usage

test_ZVD(w, test, classMeans, mus, scaling, ztol)

Arguments

Details

This function is used by other functions and should only be called explicitly for debugging purposes. Potential release in the future. This function should potentially be made internal for the release.

Value

test_ZVD returns an object of class "test_ZVD" including a list with the following named components

stats

list containing number of misclassified observations, l0 and l1 norms of discriminants.

pred_labs

predicted class labels according to nearest centroid and the discriminants.

See Also

Used by: SZVDcv.

Softmax for SZVD ADMM iterations

Description

Applies softmax the soft thresholding shrinkage operator to v with tolerance a. That is, output is the vector with entries with absolute value v_i - a if |v_i| > a and zero otherwise, with sign pattern matching that of v.

Usage

vec_shrink(v, a)

Arguments

Details

This function is used by other functions and should only be called explicitly for debugging purposes.

Value

thresholded v vector.

x,y,z

Iterates at termination.

its

Number of iterations required to converge.

errtol

Stopping error bound at termination

See Also

Used by: SZVD_ADMM.